diff --git a/IO/Image/CMakeLists.txt b/IO/Image/CMakeLists.txt index ac89187d69a12b0b43c916aef0c24b93ce2d5539..085a4d76ad4473485962c0376aafe2fd15b6fe73 100644 --- a/IO/Image/CMakeLists.txt +++ b/IO/Image/CMakeLists.txt @@ -19,6 +19,9 @@ set(Module_SRCS vtkMedicalImageReader2.cxx vtkMetaImageReader.cxx vtkMetaImageWriter.cxx + vtkNIFTIImageHeader.cxx + vtkNIFTIImageReader.cxx + vtkNIFTIImageWriter.cxx vtkNrrdReader.cxx vtkPNGReader.cxx vtkPNGWriter.cxx @@ -37,4 +40,9 @@ set_source_files_properties( ABSTRACT ) +set_source_files_properties( + vtkNIFTIPrivate.h + WRAP_EXCLUDE + ) + vtk_module_library(vtkIOImage ${Module_SRCS}) diff --git a/IO/Image/Testing/Cxx/CMakeLists.txt b/IO/Image/Testing/Cxx/CMakeLists.txt index 79d510711457b01a7f055a43ba74ee2a51de069d..db65a224986de1aa372328c0de44c2b190485719 100644 --- a/IO/Image/Testing/Cxx/CMakeLists.txt +++ b/IO/Image/Testing/Cxx/CMakeLists.txt @@ -1,6 +1,9 @@ vtk_add_test_cxx(${vtk-module}CxxTests data_tests # TestImageReader2Factory.cxx # fixme (deps not satisfied) TestNrrdReader.cxx + TestNIFTIReaderWriter.cxx + TestNIFTIReaderAnalyze.cxx + TestNIFTI2.cxx ) set(TestMetaIO_ARGS "DATA{${VTK_TEST_INPUT_DIR}/HeadMRVolume.mhd,HeadMRVolume.raw}") vtk_add_test_cxx(${vtk-module}CxxTests tests diff --git a/IO/Image/Testing/Cxx/TestNIFTI2.cxx b/IO/Image/Testing/Cxx/TestNIFTI2.cxx new file mode 100644 index 0000000000000000000000000000000000000000..7e6da488330551e2e20cbcef3805b85983ef6dc2 --- /dev/null +++ b/IO/Image/Testing/Cxx/TestNIFTI2.cxx @@ -0,0 +1,183 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: TestNIFTI2.cxx + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ +/* +Test NIFTI-2 support in VTK. +*/ + +#include "vtkNew.h" +#include "vtkTestUtilities.h" +#include "vtkRegressionTestImage.h" + +#include "vtkCamera.h" +#include "vtkImageData.h" +#include "vtkImageMathematics.h" +#include "vtkImageProperty.h" +#include "vtkImageSlice.h" +#include "vtkImageSliceMapper.h" +#include "vtkRenderWindow.h" +#include "vtkRenderWindowInteractor.h" +#include "vtkRenderer.h" + +#include "vtkNIFTIImageHeader.h" +#include "vtkNIFTIImageReader.h" +#include "vtkNIFTIImageWriter.h" + +#include <string> + +static const char *dispfile = "Data/avg152T1_RL_nifti2.nii.gz"; + +static void TestDisplay( + vtkRenderWindow *renwin, const char *infile, const char *tempDir) +{ + std::string outpath = tempDir; + outpath += "/"; + outpath += "avg152T1_RL_nifti2.nii.gz"; + + vtkNew<vtkNIFTIImageReader> reader1; + if (!reader1->CanReadFile(infile)) + { + cerr << "CanReadFile failed for " << infile << "\n"; + exit(1); + } + + reader1->SetFileName(infile); + reader1->Update(); + + vtkNew<vtkNIFTIImageWriter> writer; + writer->SetInputConnection(reader1->GetOutputPort()); + writer->SetFileName(outpath.c_str()); + writer->SetNIFTIHeader(reader1->GetNIFTIHeader()); + writer->SetSFormMatrix(reader1->GetSFormMatrix()); + writer->SetNIFTIVersion(2); + writer->Update(); + + vtkNew<vtkNIFTIImageReader> reader; + reader->SetFileName(outpath.c_str()); + reader->Update(); + + vtkNIFTIImageHeader *header = + reader->GetNIFTIHeader(); + std::string magic = header->GetMagic(); + if (magic != "n+2") + { + cerr << "File is not a NIFTIv2 file\n"; + exit(1); + } + + int size[3]; + double center[3], spacing[3]; + reader->GetOutput()->GetDimensions(size); + reader->GetOutput()->GetCenter(center); + reader->GetOutput()->GetSpacing(spacing); + double center1[3] = { center[0], center[1], center[2] }; + double center2[3] = { center[0], center[1], center[2] }; + if (size[2] % 2 == 1) + { + center1[2] += 0.5*spacing[2]; + } + if (size[0] % 2 == 1) + { + center2[0] += 0.5*spacing[0]; + } + double vrange[2]; + reader->GetOutput()->GetScalarRange(vrange); + + vtkNew<vtkImageSliceMapper> map1; + map1->BorderOn(); + map1->SliceAtFocalPointOn(); + map1->SliceFacesCameraOn(); + map1->SetInputConnection(reader->GetOutputPort()); + vtkNew<vtkImageSliceMapper> map2; + map2->BorderOn(); + map2->SliceAtFocalPointOn(); + map2->SliceFacesCameraOn(); + map2->SetInputConnection(reader->GetOutputPort()); + + vtkNew<vtkImageSlice> slice1; + slice1->SetMapper(map1.GetPointer()); + slice1->GetProperty()->SetColorWindow(vrange[1]-vrange[0]); + slice1->GetProperty()->SetColorLevel(0.5*(vrange[0]+vrange[1])); + + vtkNew<vtkImageSlice> slice2; + slice2->SetMapper(map2.GetPointer()); + slice2->GetProperty()->SetColorWindow(vrange[1]-vrange[0]); + slice2->GetProperty()->SetColorLevel(0.5*(vrange[0]+vrange[1])); + + double ratio = size[0]*1.0/(size[0]+size[2]); + + vtkNew<vtkRenderer> ren1; + ren1->SetViewport(0,0,ratio,1.0); + + vtkNew<vtkRenderer> ren2; + ren2->SetViewport(ratio,0.0,1.0,1.0); + ren1->AddViewProp(slice1.GetPointer()); + ren2->AddViewProp(slice2.GetPointer()); + + vtkCamera *cam1 = ren1->GetActiveCamera(); + cam1->ParallelProjectionOn(); + cam1->SetParallelScale(0.5*spacing[1]*size[1]); + cam1->SetFocalPoint(center1[0], center1[1], center1[2]); + cam1->SetPosition(center1[0], center1[1], center1[2] - 100.0); + + vtkCamera *cam2 = ren2->GetActiveCamera(); + cam2->ParallelProjectionOn(); + cam2->SetParallelScale(0.5*spacing[1]*size[1]); + cam2->SetFocalPoint(center2[0], center2[1], center2[2]); + cam2->SetPosition(center2[0] + 100.0, center2[1], center2[2]); + + renwin->SetSize(size[0] + size[2], size[1]); + renwin->AddRenderer(ren1.GetPointer()); + renwin->AddRenderer(ren2.GetPointer()); +}; + +int TestNIFTI2(int argc, char *argv[]) +{ + // perform the display test + char *infile = + vtkTestUtilities::ExpandDataFileName(argc, argv, dispfile); + if (!infile) + { + cerr << "Could not locate input file " << dispfile << "\n"; + return 1; + } + std::string inpath = infile; + delete [] infile; + + char *tempDir = vtkTestUtilities::GetArgOrEnvOrDefault( + "-T", argc, argv, "VTK_TEMP_DIR", "Testing/Temporary"); + if (!tempDir) + { + cerr << "Could not determine temporary directory.\n"; + return 1; + } + std::string tmppath = tempDir; + delete [] tempDir; + + vtkNew<vtkRenderWindow> renwin; + vtkNew<vtkRenderWindowInteractor> iren; + iren->SetRenderWindow(renwin.GetPointer()); + + TestDisplay(renwin.GetPointer(), inpath.c_str(), tmppath.c_str()); + + int retVal = vtkRegressionTestImage(renwin.GetPointer()); + if (retVal == vtkRegressionTester::DO_INTERACTOR) + { + renwin->Render(); + iren->Start(); + retVal = vtkRegressionTester::PASSED; + } + + return (retVal != vtkRegressionTester::PASSED); +} diff --git a/IO/Image/Testing/Cxx/TestNIFTIReaderAnalyze.cxx b/IO/Image/Testing/Cxx/TestNIFTIReaderAnalyze.cxx new file mode 100644 index 0000000000000000000000000000000000000000..e8ebc8e8e7edfb9e8218053aa554542db6d1a22f --- /dev/null +++ b/IO/Image/Testing/Cxx/TestNIFTIReaderAnalyze.cxx @@ -0,0 +1,146 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: TestNIFTIReaderAnalyze.cxx + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ +/* +Test compatibility of the vtkNIFTIImageReader with Analyze 7.5 files. +*/ + +#include "vtkNew.h" +#include "vtkTestUtilities.h" +#include "vtkRegressionTestImage.h" + +#include "vtkCamera.h" +#include "vtkImageData.h" +#include "vtkImageMathematics.h" +#include "vtkImageProperty.h" +#include "vtkImageSlice.h" +#include "vtkImageSliceMapper.h" +#include "vtkRenderWindow.h" +#include "vtkRenderWindowInteractor.h" +#include "vtkRenderer.h" + +#include "vtkNIFTIImageHeader.h" +#include "vtkNIFTIImageReader.h" +#include "vtkNIFTIImageWriter.h" + +#include <string> + +static const char *dispfile = "Data/ANALYZE.HDR"; + +static void TestDisplay(vtkRenderWindow *renwin, const char *infile) +{ + vtkNew<vtkNIFTIImageReader> reader; + if (!reader->CanReadFile(infile)) + { + cerr << "CanReadFile failed for " << infile << "\n"; + exit(1); + } + reader->SetFileName(infile); + reader->Update(); + + int size[3]; + double center[3], spacing[3]; + reader->GetOutput()->GetDimensions(size); + reader->GetOutput()->GetCenter(center); + reader->GetOutput()->GetSpacing(spacing); + double center1[3] = { center[0], center[1], center[2] }; + double center2[3] = { center[0], center[1], center[2] }; + if (size[2] % 2 == 1) + { + center1[2] += 0.5*spacing[2]; + } + if (size[0] % 2 == 1) + { + center2[0] += 0.5*spacing[0]; + } + double vrange[2]; + reader->GetOutput()->GetScalarRange(vrange); + + vtkNew<vtkImageSliceMapper> map1; + map1->BorderOn(); + map1->SliceAtFocalPointOn(); + map1->SliceFacesCameraOn(); + map1->SetInputConnection(reader->GetOutputPort()); + vtkNew<vtkImageSliceMapper> map2; + map2->BorderOn(); + map2->SliceAtFocalPointOn(); + map2->SliceFacesCameraOn(); + map2->SetInputConnection(reader->GetOutputPort()); + + vtkNew<vtkImageSlice> slice1; + slice1->SetMapper(map1.GetPointer()); + slice1->GetProperty()->SetColorWindow(vrange[1]-vrange[0]); + slice1->GetProperty()->SetColorLevel(0.5*(vrange[0]+vrange[1])); + + vtkNew<vtkImageSlice> slice2; + slice2->SetMapper(map2.GetPointer()); + slice2->GetProperty()->SetColorWindow(vrange[1]-vrange[0]); + slice2->GetProperty()->SetColorLevel(0.5*(vrange[0]+vrange[1])); + + double ratio = size[0]*1.0/(size[0]+size[2]); + + vtkNew<vtkRenderer> ren1; + ren1->SetViewport(0,0,ratio,1.0); + + vtkNew<vtkRenderer> ren2; + ren2->SetViewport(ratio,0.0,1.0,1.0); + ren1->AddViewProp(slice1.GetPointer()); + ren2->AddViewProp(slice2.GetPointer()); + + vtkCamera *cam1 = ren1->GetActiveCamera(); + cam1->ParallelProjectionOn(); + cam1->SetParallelScale(0.5*spacing[1]*size[1]); + cam1->SetFocalPoint(center1[0], center1[1], center1[2]); + cam1->SetPosition(center1[0], center1[1], center1[2] - 100.0); + + vtkCamera *cam2 = ren2->GetActiveCamera(); + cam2->ParallelProjectionOn(); + cam2->SetParallelScale(0.5*spacing[1]*size[1]); + cam2->SetFocalPoint(center2[0], center2[1], center2[2]); + cam2->SetPosition(center2[0] + 100.0, center2[1], center2[2]); + + renwin->SetSize(size[0] + size[2], size[1]); + renwin->AddRenderer(ren1.GetPointer()); + renwin->AddRenderer(ren2.GetPointer()); +}; + +int TestNIFTIReaderAnalyze(int argc, char *argv[]) +{ + // perform the display test + char *infile = + vtkTestUtilities::ExpandDataFileName(argc, argv, dispfile); + if (!infile) + { + cerr << "Could not locate input file " << dispfile << "\n"; + return 1; + } + std::string inpath = infile; + delete [] infile; + + vtkNew<vtkRenderWindow> renwin; + vtkNew<vtkRenderWindowInteractor> iren; + iren->SetRenderWindow(renwin.GetPointer()); + + TestDisplay(renwin.GetPointer(), inpath.c_str()); + + int retVal = vtkRegressionTestImage(renwin.GetPointer()); + if (retVal == vtkRegressionTester::DO_INTERACTOR) + { + renwin->Render(); + iren->Start(); + retVal = vtkRegressionTester::PASSED; + } + + return (retVal != vtkRegressionTester::PASSED); +} diff --git a/IO/Image/Testing/Cxx/TestNIFTIReaderWriter.cxx b/IO/Image/Testing/Cxx/TestNIFTIReaderWriter.cxx new file mode 100644 index 0000000000000000000000000000000000000000..2184e6ab8297637985c7111f625b30eb43984db8 --- /dev/null +++ b/IO/Image/Testing/Cxx/TestNIFTIReaderWriter.cxx @@ -0,0 +1,345 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: TestNIFTIReaderWriter.cxx + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ +/* +Test NIFTI support in VTK by reading a file, writing it, and +then re-reading it to ensure that the contents are identical. +*/ + +#include "vtkNew.h" +#include "vtkTestUtilities.h" +#include "vtkRegressionTestImage.h" + +#include "vtkCamera.h" +#include "vtkImageData.h" +#include "vtkImageMathematics.h" +#include "vtkImageProperty.h" +#include "vtkImageSlice.h" +#include "vtkImageSliceMapper.h" +#include "vtkMatrix4x4.h" +#include "vtkRenderWindow.h" +#include "vtkRenderWindowInteractor.h" +#include "vtkRenderer.h" +#include "vtkStringArray.h" + +#include "vtkNIFTIImageHeader.h" +#include "vtkNIFTIImageReader.h" +#include "vtkNIFTIImageWriter.h" + +#include <string> + +static const char *testfiles[6][2] = { + { "Data/minimal.nii.gz", "out_minimal.nii.gz" }, + { "Data/minimal.img.gz", "out_minimal.hdr" }, + { "Data/nifti_rgb.nii.gz", "out_nifti_rgb.nii" }, + { "Data/filtered_func_data.nii.gz", "out_filtered_func_data.nii.gz" }, + { "Data/minimal.hdr.gz", "out_minimal_2.nii" }, + { "Data/minimal.img.gz", "" } +}; + +static const char *dispfile = "Data/avg152T1_RL_nifti.nii.gz"; + +static void TestDisplay(vtkRenderWindow *renwin, const char *infile) +{ + vtkNew<vtkNIFTIImageReader> reader; + reader->SetFileName(infile); + reader->Update(); + + int size[3]; + double center[3], spacing[3]; + reader->GetOutput()->GetDimensions(size); + reader->GetOutput()->GetCenter(center); + reader->GetOutput()->GetSpacing(spacing); + double center1[3] = { center[0], center[1], center[2] }; + double center2[3] = { center[0], center[1], center[2] }; + if (size[2] % 2 == 1) + { + center1[2] += 0.5*spacing[2]; + } + if (size[0] % 2 == 1) + { + center2[0] += 0.5*spacing[0]; + } + double vrange[2]; + reader->GetOutput()->GetScalarRange(vrange); + + vtkNew<vtkImageSliceMapper> map1; + map1->BorderOn(); + map1->SliceAtFocalPointOn(); + map1->SliceFacesCameraOn(); + map1->SetInputConnection(reader->GetOutputPort()); + vtkNew<vtkImageSliceMapper> map2; + map2->BorderOn(); + map2->SliceAtFocalPointOn(); + map2->SliceFacesCameraOn(); + map2->SetInputConnection(reader->GetOutputPort()); + + vtkNew<vtkImageSlice> slice1; + slice1->SetMapper(map1.GetPointer()); + slice1->GetProperty()->SetColorWindow(vrange[1]-vrange[0]); + slice1->GetProperty()->SetColorLevel(0.5*(vrange[0]+vrange[1])); + + vtkNew<vtkImageSlice> slice2; + slice2->SetMapper(map2.GetPointer()); + slice2->GetProperty()->SetColorWindow(vrange[1]-vrange[0]); + slice2->GetProperty()->SetColorLevel(0.5*(vrange[0]+vrange[1])); + + double ratio = size[0]*1.0/(size[0]+size[2]); + + vtkNew<vtkRenderer> ren1; + ren1->SetViewport(0,0,ratio,1.0); + + vtkNew<vtkRenderer> ren2; + ren2->SetViewport(ratio,0.0,1.0,1.0); + ren1->AddViewProp(slice1.GetPointer()); + ren2->AddViewProp(slice2.GetPointer()); + + vtkCamera *cam1 = ren1->GetActiveCamera(); + cam1->ParallelProjectionOn(); + cam1->SetParallelScale(0.5*spacing[1]*size[1]); + cam1->SetFocalPoint(center1[0], center1[1], center1[2]); + cam1->SetPosition(center1[0], center1[1], center1[2] - 100.0); + + vtkCamera *cam2 = ren2->GetActiveCamera(); + cam2->ParallelProjectionOn(); + cam2->SetParallelScale(0.5*spacing[1]*size[1]); + cam2->SetFocalPoint(center2[0], center2[1], center2[2]); + cam2->SetPosition(center2[0] + 100.0, center2[1], center2[2]); + + renwin->SetSize(size[0] + size[2], size[1]); + renwin->AddRenderer(ren1.GetPointer()); + renwin->AddRenderer(ren2.GetPointer()); +}; + +static double TestReadWriteRead( + const char *infile, const char *infile2, const char *outfile) +{ + // read a NIFTI file + vtkNew<vtkNIFTIImageReader> reader; + if (infile2 == 0) + { + reader->SetFileName(infile); + } + else + { + vtkNew<vtkStringArray> filenames; + filenames->InsertNextValue(infile); + filenames->InsertNextValue(infile2); + reader->SetFileNames(filenames.GetPointer()); + } + reader->TimeAsVectorOn(); + reader->Update(); + + vtkNew<vtkNIFTIImageWriter> writer; + writer->SetInputConnection(reader->GetOutputPort()); + writer->SetFileName(outfile); + // copy most information directly from the header + vtkNIFTIImageHeader *header = writer->GetNIFTIHeader(); + header->DeepCopy(reader->GetNIFTIHeader()); + header->SetDescrip("VTK Test Data"); + // this information will override the reader's header + writer->SetQFac(reader->GetQFac()); + writer->SetTimeDimension(reader->GetTimeDimension()); + writer->SetTimeSpacing(reader->GetTimeSpacing()); + writer->SetRescaleSlope(reader->GetRescaleSlope()); + writer->SetRescaleIntercept(reader->GetRescaleIntercept()); + writer->SetQFormMatrix(reader->GetQFormMatrix()); + if (reader->GetSFormMatrix()) + { + writer->SetSFormMatrix(reader->GetSFormMatrix()); + } + else + { + writer->SetSFormMatrix(reader->GetQFormMatrix()); + } + writer->Write(); + + // to exercise PrintSelf + vtkOStrStreamWrapper s; + reader->Print(s); + header->Print(s); + writer->Print(s); + + vtkNew<vtkNIFTIImageReader> reader2; + reader2->SetFileName(outfile); + reader2->TimeAsVectorOn(); + reader2->Update(); + + vtkNew<vtkImageMathematics> diff; + diff->SetOperationToSubtract(); + diff->SetInputConnection(0,reader->GetOutputPort()); + diff->SetInputConnection(1,reader2->GetOutputPort()); + diff->Update(); + double diffrange[2]; + diff->GetOutput()->GetScalarRange(diffrange); + double differr = diffrange[0]*diffrange[0] + diffrange[1]*diffrange[1]; + + return differr; +} + +static int TestNIFTIHeader() +{ + vtkNew<vtkNIFTIImageHeader> header1; + vtkNew<vtkNIFTIImageHeader> header2; + + header1->SetIntentCode(vtkNIFTIImageHeader::IntentZScore); + header1->SetIntentName("ZScore"); + header1->SetIntentP1(1.0); + header1->SetIntentP2(2.0); + header1->SetIntentP3(3.0); + header1->SetSclSlope(2.0); + header1->SetSclInter(1024.0); + header1->SetCalMin(-1024.0); + header1->SetCalMax(3072.0); + header1->SetSliceDuration(1.0); + header1->SetSliceStart(2); + header1->SetSliceEnd(14); + header1->SetXYZTUnits( + vtkNIFTIImageHeader::UnitsMM | vtkNIFTIImageHeader::UnitsSec); + header1->SetDimInfo(0); + header1->SetDescrip("Test header"); + header1->SetAuxFile("none"); + header1->SetQFormCode(vtkNIFTIImageHeader::XFormScannerAnat); + header1->SetQuaternB(0.0); + header1->SetQuaternC(1.0); + header1->SetQuaternD(0.0); + header1->SetQOffsetX(10.0); + header1->SetQOffsetY(30.0); + header1->SetQOffsetZ(20.0); + header1->SetSFormCode(vtkNIFTIImageHeader::XFormAlignedAnat); + double matrix[16]; + vtkMatrix4x4::Identity(matrix); + header1->SetSRowX(matrix); + header1->SetSRowY(matrix+4); + header1->SetSRowZ(matrix+8); + + header2->DeepCopy(header1.GetPointer()); + bool success = true; + success &= (header2->GetIntentCode() == vtkNIFTIImageHeader::IntentZScore); + success &= (strcmp(header2->GetIntentName(), "ZScore") == 0); + success &= (header2->GetIntentP1() == 1.0); + success &= (header2->GetIntentP2() == 2.0); + success &= (header2->GetIntentP3() == 3.0); + success &= (header2->GetSclSlope() == 2.0); + success &= (header2->GetSclInter() == 1024.0); + success &= (header2->GetCalMin() == -1024.0); + success &= (header2->GetCalMax() == 3072.0); + success &= (header2->GetSliceDuration() == 1.0); + success &= (header2->GetSliceStart() == 2); + success &= (header2->GetSliceEnd() == 14); + success &= (header2->GetXYZTUnits() == + (vtkNIFTIImageHeader::UnitsMM | vtkNIFTIImageHeader::UnitsSec)); + success &= (header2->GetDimInfo() == 0); + success &= (strcmp(header2->GetDescrip(), "Test header") == 0); + success &= (strcmp(header2->GetAuxFile(), "none") == 0); + success &= (header2->GetQFormCode() == vtkNIFTIImageHeader::XFormScannerAnat); + success &= (header2->GetQuaternB() == 0.0); + success &= (header2->GetQuaternC() == 1.0); + success &= (header2->GetQuaternD() == 0.0); + success &= (header2->GetQOffsetX() == 10.0); + success &= (header2->GetQOffsetY() == 30.0); + success &= (header2->GetQOffsetZ() == 20.0); + header2->GetSRowX(matrix); + header2->GetSRowY(matrix + 4); + header2->GetSRowZ(matrix + 8); + success &= (matrix[0] == 1.0 && matrix[5] == 1.0 && matrix[10] == 1.0); + + return success; +} + +int TestNIFTIReaderWriter(int argc, char *argv[]) +{ + // perform the header test + if (!TestNIFTIHeader()) + { + cerr << "Failed TestNIFTIHeader\n"; + return 1; + } + + // perform the read/write test + for (int i = 0; i < 5; i++) + { + char *infile2 = 0; + char *infile = + vtkTestUtilities::ExpandDataFileName(argc, argv, testfiles[i][0]); + if (i == 4) + { + infile2 = + vtkTestUtilities::ExpandDataFileName(argc, argv, testfiles[5][0]); + } + if (!infile) + { + cerr << "Could not locate input file " << testfiles[i][0] << "\n"; + return 1; + } + + char *tempDir = vtkTestUtilities::GetArgOrEnvOrDefault( + "-T", argc, argv, "VTK_TEMP_DIR", "Testing/Temporary"); + if (!tempDir) + { + cerr << "Could not determine temporary directory.\n"; + return 1; + } + + std::string outpath = tempDir; + outpath += "/"; + outpath += testfiles[i][1]; + delete [] tempDir; + + vtkNew<vtkNIFTIImageReader> testReader; + testReader->GetFileExtensions(); + testReader->GetDescriptiveName(); + if (!testReader->CanReadFile(infile)) + { + cerr << "CanReadFile() failed for " << infile << "\n"; + return 1; + } + + double err = TestReadWriteRead(infile, infile2, outpath.c_str()); + if (err != 0.0) + { + cerr << "Input " << infile << " differs from output " << outpath.c_str() + << "\n"; + return 1; + } + delete [] infile; + delete [] infile2; + } + + // perform the display test + char *infile = + vtkTestUtilities::ExpandDataFileName(argc, argv, dispfile); + if (!infile) + { + cerr << "Could not locate input file " << dispfile << "\n"; + return 1; + } + + vtkNew<vtkRenderWindow> renwin; + vtkNew<vtkRenderWindowInteractor> iren; + iren->SetRenderWindow(renwin.GetPointer()); + + TestDisplay(renwin.GetPointer(), infile); + delete [] infile; + + int retVal = vtkRegressionTestImage(renwin.GetPointer()); + if (retVal == vtkRegressionTester::DO_INTERACTOR) + { + renwin->Render(); + iren->Start(); + retVal = vtkRegressionTester::PASSED; + } + + return (retVal != vtkRegressionTester::PASSED); +} diff --git a/IO/Image/Testing/Data/Baseline/TestNIFTI2.png.md5 b/IO/Image/Testing/Data/Baseline/TestNIFTI2.png.md5 new file mode 100644 index 0000000000000000000000000000000000000000..1f0fb3416dad1c69bae1bd0c286d683221af81f0 --- /dev/null +++ b/IO/Image/Testing/Data/Baseline/TestNIFTI2.png.md5 @@ -0,0 +1 @@ +6d7083f60dfc333062876ad36836159e diff --git a/IO/Image/Testing/Data/Baseline/TestNIFTIReaderAnalyze.png.md5 b/IO/Image/Testing/Data/Baseline/TestNIFTIReaderAnalyze.png.md5 new file mode 100644 index 0000000000000000000000000000000000000000..4b49491d17cf08fce8fb510da9ed63cfe95e91b7 --- /dev/null +++ b/IO/Image/Testing/Data/Baseline/TestNIFTIReaderAnalyze.png.md5 @@ -0,0 +1 @@ +a21ab5999bf7df1c0d2c013f097fff18 diff --git a/IO/Image/Testing/Data/Baseline/TestNIFTIReaderWriter.png.md5 b/IO/Image/Testing/Data/Baseline/TestNIFTIReaderWriter.png.md5 new file mode 100644 index 0000000000000000000000000000000000000000..7e0a61edf707c69313dda0ffb63b29d2c5aa6859 --- /dev/null +++ b/IO/Image/Testing/Data/Baseline/TestNIFTIReaderWriter.png.md5 @@ -0,0 +1 @@ +49e8210a47e8e2eb657db6c4b7973d86 diff --git a/IO/Image/Testing/Python/CMakeLists.txt b/IO/Image/Testing/Python/CMakeLists.txt index 6c1773e6df954d702a8bf34304110e92e6b90a6c..75e49868362d61c1a78d461bd6d0fdd05d639cbd 100644 --- a/IO/Image/Testing/Python/CMakeLists.txt +++ b/IO/Image/Testing/Python/CMakeLists.txt @@ -6,6 +6,7 @@ vtk_add_test_python( TestTIFFReader.py dem.py TestMetaImage2D.py + TestNIFTIReaderWriter.py TestSetFileNames.py TestImageJSONWriter.py,NO_VALID ) diff --git a/IO/Image/Testing/Python/TestNIFTIReaderWriter.py b/IO/Image/Testing/Python/TestNIFTIReaderWriter.py new file mode 100644 index 0000000000000000000000000000000000000000..42400986a59cf69303eb669330895c572352c6e8 --- /dev/null +++ b/IO/Image/Testing/Python/TestNIFTIReaderWriter.py @@ -0,0 +1,154 @@ +#! /usr/bin/env python + +""" +Test NIFTI support in VTK by reading a file, writing it, and +then re-reading it to ensure that the contents are identical. +""" + +import vtk +from vtk.test import Testing +from vtk.util.misc import vtkGetDataRoot +from vtk.util.misc import vtkGetTempDir + +VTK_DATA_ROOT = vtkGetDataRoot() +VTK_TEMP_DIR = vtkGetTempDir() + +import sys +import os + +testfiles = [ + ["minimal.nii.gz", "out_minimal.nii.gz"], + ["minimal.img.gz", "out_minimal.hdr"] +] + +dispfile = "avg152T1_RL_nifti.nii.gz" + +def TestDisplay(file1): + """Display the output""" + + inpath = os.path.join(str(VTK_DATA_ROOT), "Data", file1) + + reader = vtk.vtkNIFTIImageReader() + reader.SetFileName(inpath) + reader.Update() + + size = reader.GetOutput().GetDimensions() + center = reader.GetOutput().GetCenter() + spacing = reader.GetOutput().GetSpacing() + center1 = (center[0], center[1], center[2]) + center2 = (center[0], center[1], center[2]) + if size[2] % 2 == 1: + center1 = (center[0], center[1], center[2] + 0.5*spacing[2]) + if size[0] % 2 == 1: + center2 = (center[0] + 0.5*spacing[0], center[1], center[2]) + vrange = reader.GetOutput().GetScalarRange() + + map1 = vtk.vtkImageSliceMapper() + map1.BorderOn() + map1.SliceAtFocalPointOn() + map1.SliceFacesCameraOn() + map1.SetInputConnection(reader.GetOutputPort()) + map2 = vtk.vtkImageSliceMapper() + map2.BorderOn() + map2.SliceAtFocalPointOn() + map2.SliceFacesCameraOn() + map2.SetInputConnection(reader.GetOutputPort()) + + slice1 = vtk.vtkImageSlice() + slice1.SetMapper(map1) + slice1.GetProperty().SetColorWindow(vrange[1]-vrange[0]) + slice1.GetProperty().SetColorLevel(0.5*(vrange[0]+vrange[1])) + slice2 = vtk.vtkImageSlice() + slice2.SetMapper(map2) + slice2.GetProperty().SetColorWindow(vrange[1]-vrange[0]) + slice2.GetProperty().SetColorLevel(0.5*(vrange[0]+vrange[1])) + + ratio = size[0]*1.0/(size[0]+size[2]) + + ren1 = vtk.vtkRenderer() + ren1.SetViewport(0,0,ratio,1.0) + ren2 = vtk.vtkRenderer() + ren2.SetViewport(ratio,0.0,1.0,1.0) + ren1.AddViewProp(slice1) + ren2.AddViewProp(slice2) + + cam1 = ren1.GetActiveCamera() + cam1.ParallelProjectionOn() + cam1.SetParallelScale(0.5*spacing[1]*size[1]) + cam1.SetFocalPoint(center1[0], center1[1], center1[2]) + cam1.SetPosition(center1[0], center1[1], center1[2] - 100.0) + + cam2 = ren2.GetActiveCamera() + cam2.ParallelProjectionOn() + cam2.SetParallelScale(0.5*spacing[1]*size[1]) + cam2.SetFocalPoint(center2[0], center2[1], center2[2]) + cam2.SetPosition(center2[0] + 100.0, center2[1], center2[2]) + + if "-I" in sys.argv: + style = vtk.vtkInteractorStyleImage() + style.SetInteractionModeToImageSlicing() + + iren = vtk.vtkRenderWindowInteractor() + iren.SetInteractorStyle(style) + + renwin = vtk.vtkRenderWindow() + renwin.SetSize(size[0] + size[2], size[1]) + renwin.AddRenderer(ren1) + renwin.AddRenderer(ren2) + + renwin.Render() + + if "-I" in sys.argv: + renwin.SetInteractor(iren) + iren.Initialize() + iren.Start() + + return renwin + +def TestReadWriteRead(infile, outfile): + """Read, write, and re-read a file, return difference.""" + + inpath = os.path.join(str(VTK_DATA_ROOT), "Data", infile) + outpath = os.path.join(str(VTK_TEMP_DIR), outfile) + + # read a NIFTI file + reader = vtk.vtkNIFTIImageReader() + reader.SetFileName(inpath) + reader.TimeAsVectorOn() + reader.Update() + + writer = vtk.vtkNIFTIImageWriter() + writer.SetInputConnection(reader.GetOutputPort()) + writer.SetFileName(outpath) + # copy most information directoy from the header + writer.SetNIFTIHeader(reader.GetNIFTIHeader()) + # this information will override the reader's header + writer.SetQFac(reader.GetQFac()) + writer.SetTimeDimension(reader.GetTimeDimension()) + writer.SetQFormMatrix(reader.GetQFormMatrix()) + writer.SetSFormMatrix(reader.GetSFormMatrix()) + writer.Write() + + reader2 = vtk.vtkNIFTIImageReader() + reader2.SetFileName(outpath) + reader2.TimeAsVectorOn() + reader2.Update() + + diff = vtk.vtkImageMathematics() + diff.SetOperationToSubtract() + diff.SetInputConnection(0,reader.GetOutputPort()) + diff.SetInputConnection(1,reader2.GetOutputPort()) + diff.Update() + diffrange = diff.GetOutput().GetScalarRange() + differr = diffrange[0]**2 + diffrange[1]**2 + + return differr + +for infile, outfile in testfiles: + err = TestReadWriteRead(infile, outfile) + if err: + sys.stderr.write( + "Input " + infile + " differs from outfile " + outfile) + sys.exit(1) + +renWin = TestDisplay(dispfile) diff --git a/IO/Image/module.cmake b/IO/Image/module.cmake index 751f053c7581608c470f9f8addddb9cd71275af3..086e41e5cfc17cba78748d04c01deaa90390dc25 100644 --- a/IO/Image/module.cmake +++ b/IO/Image/module.cmake @@ -20,6 +20,7 @@ vtk_module(vtkIOImage vtkTestingCore vtkImagingSources vtkImagingMath + vtkInteractionStyle vtkRenderingContext2D vtkRenderingCore vtkTestingCore diff --git a/IO/Image/vtkNIFTIImageHeader.cxx b/IO/Image/vtkNIFTIImageHeader.cxx new file mode 100644 index 0000000000000000000000000000000000000000..19581d8d2374894ff7a189ecc59965fda1667b32 --- /dev/null +++ b/IO/Image/vtkNIFTIImageHeader.cxx @@ -0,0 +1,453 @@ +/*========================================================================= + +Program: Visualization Toolkit +Module: vtkNIFTIImageHeader.cxx + +Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen +All rights reserved. +See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ + +#include "vtkNIFTIImageHeader.h" +#include "vtkNIFTIImagePrivate.h" + +#include <vtkObjectFactory.h> + +#include <string.h> +#include <float.h> +#include <math.h> +#include <ctype.h> + +vtkStandardNewMacro(vtkNIFTIImageHeader); + +//---------------------------------------------------------------------------- +namespace { + +// utility function to normalize floats that are close to zero +double vtkNIFTINormalizeFloat(double d) +{ + return (fabs(d) < FLT_MIN ? 0.0 : d); +} + +double vtkNIFTINormalizeDouble(double d) +{ + return (fabs(d) < DBL_MIN ? 0.0 : d); +} +} // end anonymous namespace + +//---------------------------------------------------------------------------- +vtkNIFTIImageHeader::vtkNIFTIImageHeader() +{ + this->Initialize(); +} + +//---------------------------------------------------------------------------- +vtkNIFTIImageHeader::~vtkNIFTIImageHeader() +{ +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::Initialize() +{ + memset(this->Magic, '\0', sizeof(this->Magic)); + this->VoxOffset = 0; + this->DataType = 0; + this->BitPix = 0; + for (int i = 0; i < 8; i++) + { + this->Dim[i] = 0; + this->PixDim[i] = 0.0; + } + this->IntentCode = 0; + memset(this->IntentName, '\0', sizeof(this->IntentName)); + this->IntentP1 = 0.0; + this->IntentP2 = 0.0; + this->IntentP3 = 0.0; + this->SclSlope = 0.0; + this->SclInter = 0.0; + this->CalMin = 0.0; + this->CalMax = 0.0; + this->SliceDuration = 0.0; + this->TOffset = 0.0; + this->SliceStart = 0; + this->SliceEnd = 0; + this->SliceCode = 0; + this->XYZTUnits = 0; + this->DimInfo = 0; + memset(this->Descrip, '\0', sizeof(this->Descrip)); + memset(this->AuxFile, '\0', sizeof(this->AuxFile)); + this->QFormCode = 0; + this->SFormCode = 0; + this->QuaternB = 0.0; + this->QuaternC = 0.0; + this->QuaternD = 0.0; + this->QOffsetX = 0.0; + this->QOffsetY = 0.0; + this->QOffsetZ = 0.0; + for (int i = 0; i < 4; i++) + { + this->SRowX[i] = 0.0; + this->SRowY[i] = 0.0; + this->SRowZ[i] = 0.0; + } +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::SetHeader(const nifti_1_header *hdr) +{ + // clear all fields (in case supplied header is Analyze 7.5) + this->Initialize(); + + // check if header is NIfTI (vs. Analyze 7.5) + bool isnifti = (hdr->magic[0] == 'n' && + (hdr->magic[1] == '+' || hdr->magic[1] == 'i') && + hdr->magic[2] == '1' && + hdr->magic[3] == '\0'); + + if (isnifti) + { + memcpy(this->Magic, hdr->magic, sizeof(hdr->magic)); + } + this->VoxOffset = static_cast<vtkTypeInt64>(hdr->vox_offset); + this->DataType = hdr->datatype; + this->BitPix = hdr->bitpix; + for (int i = 0; i < 8; i++) + { + this->Dim[i] = hdr->dim[i]; + this->PixDim[i] = hdr->pixdim[i]; + } + if (isnifti) + { + this->IntentCode = hdr->intent_code; + strncpy(this->IntentName, hdr->intent_name, sizeof(hdr->intent_name)); + this->IntentP1 = hdr->intent_p1; + this->IntentP2 = hdr->intent_p2; + this->IntentP3 = hdr->intent_p3; + this->SclSlope = hdr->scl_slope; + this->SclInter = hdr->scl_inter; + } + this->CalMin = hdr->cal_min; + this->CalMax = hdr->cal_max; + if (isnifti) + { + this->SliceDuration = hdr->slice_duration; + this->TOffset = hdr->toffset; + this->SliceStart = hdr->slice_start; + this->SliceEnd = hdr->slice_end; + this->SliceCode = hdr->slice_code; + } + this->XYZTUnits = hdr->xyzt_units; + this->DimInfo = hdr->dim_info; + strncpy(this->Descrip, hdr->descrip, sizeof(hdr->descrip)); + strncpy(this->AuxFile, hdr->aux_file, sizeof(hdr->aux_file)); + if (isnifti) + { + this->QFormCode = hdr->qform_code; + this->SFormCode = hdr->sform_code; + this->QuaternB = hdr->quatern_b; + this->QuaternC = hdr->quatern_c; + this->QuaternD = hdr->quatern_d; + this->QOffsetX = hdr->qoffset_x; + this->QOffsetY = hdr->qoffset_y; + this->QOffsetZ = hdr->qoffset_z; + for (int i = 0; i < 4; i++) + { + this->SRowX[i] = hdr->srow_x[i]; + this->SRowY[i] = hdr->srow_y[i]; + this->SRowZ[i] = hdr->srow_z[i]; + } + } +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::GetHeader(nifti_1_header *hdr) +{ + hdr->sizeof_hdr = NIFTI1HeaderSize; + memcpy(hdr->magic, this->Magic, sizeof(hdr->magic)); + memset(hdr->data_type, '\0', 10); + memset(hdr->db_name, '\0', 18); + hdr->extents = 0; + hdr->session_error = 0; + hdr->regular = 0; + hdr->dim_info = this->DimInfo; + hdr->intent_p1 = vtkNIFTINormalizeFloat(this->IntentP1); + hdr->intent_p2 = vtkNIFTINormalizeFloat(this->IntentP2); + hdr->intent_p3 = vtkNIFTINormalizeFloat(this->IntentP3); + hdr->intent_code = static_cast<short>(this->IntentCode); + hdr->datatype = static_cast<short>(this->DataType); + hdr->bitpix = static_cast<short>(this->BitPix); + hdr->slice_start = this->SliceStart; + for (int i = 0; i < 8; i++) + { + hdr->dim[i] = static_cast<short>(this->Dim[i]); + hdr->pixdim[i] = vtkNIFTINormalizeFloat(this->PixDim[i]); + } + hdr->vox_offset = static_cast<float>(this->VoxOffset); + strncpy(hdr->intent_name, this->IntentName, sizeof(hdr->intent_name)); + hdr->scl_slope = vtkNIFTINormalizeFloat(this->SclSlope); + hdr->scl_inter = vtkNIFTINormalizeFloat(this->SclInter); + hdr->cal_min = vtkNIFTINormalizeFloat(this->CalMin); + hdr->cal_max = vtkNIFTINormalizeFloat(this->CalMax); + hdr->slice_duration = vtkNIFTINormalizeFloat(this->SliceDuration); + hdr->toffset = vtkNIFTINormalizeFloat(this->TOffset); + hdr->glmax = 0; + hdr->glmin = 0; + hdr->slice_end = this->SliceEnd; + hdr->slice_code = this->SliceCode; + hdr->xyzt_units = this->XYZTUnits; + strncpy(hdr->descrip, this->Descrip, sizeof(hdr->descrip)); + strncpy(hdr->aux_file, this->AuxFile, sizeof(hdr->aux_file)); + hdr->qform_code = static_cast<short>(this->QFormCode); + hdr->sform_code = static_cast<short>(this->SFormCode); + hdr->quatern_b = vtkNIFTINormalizeFloat(this->QuaternB); + hdr->quatern_c = vtkNIFTINormalizeFloat(this->QuaternC); + hdr->quatern_d = vtkNIFTINormalizeFloat(this->QuaternD); + hdr->qoffset_x = vtkNIFTINormalizeFloat(this->QOffsetX); + hdr->qoffset_y = vtkNIFTINormalizeFloat(this->QOffsetY); + hdr->qoffset_z = vtkNIFTINormalizeFloat(this->QOffsetZ); + for (int i = 0; i < 4; i++) + { + hdr->srow_x[i] = vtkNIFTINormalizeFloat(this->SRowX[i]); + hdr->srow_y[i] = vtkNIFTINormalizeFloat(this->SRowY[i]); + hdr->srow_z[i] = vtkNIFTINormalizeFloat(this->SRowZ[i]); + } +} + + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::SetHeader(const nifti_2_header *hdr) +{ + memcpy(this->Magic, hdr->magic, sizeof(hdr->magic)); + this->VoxOffset = hdr->vox_offset; + this->DataType = hdr->datatype; + this->BitPix = hdr->bitpix; + for (int i = 0; i < 8; i++) + { + this->Dim[i] = hdr->dim[i]; + this->PixDim[i] = hdr->pixdim[i]; + } + this->IntentCode = hdr->intent_code; + strncpy(this->IntentName, hdr->intent_name, sizeof(hdr->intent_name)); + this->IntentP1 = hdr->intent_p1; + this->IntentP2 = hdr->intent_p2; + this->IntentP3 = hdr->intent_p3; + this->SclSlope = hdr->scl_slope; + this->SclInter = hdr->scl_inter; + this->CalMin = hdr->cal_min; + this->CalMax = hdr->cal_max; + this->SliceDuration = hdr->slice_duration; + this->TOffset = hdr->toffset; + this->SliceStart = hdr->slice_start; + this->SliceEnd = hdr->slice_end; + this->SliceCode = hdr->slice_code; + this->XYZTUnits = hdr->xyzt_units; + this->DimInfo = hdr->dim_info; + strncpy(this->Descrip, hdr->descrip, sizeof(hdr->descrip)); + strncpy(this->AuxFile, hdr->aux_file, sizeof(hdr->aux_file)); + this->QFormCode = hdr->qform_code; + this->SFormCode = hdr->sform_code; + this->QuaternB = hdr->quatern_b; + this->QuaternC = hdr->quatern_c; + this->QuaternD = hdr->quatern_d; + this->QOffsetX = hdr->qoffset_x; + this->QOffsetY = hdr->qoffset_y; + this->QOffsetZ = hdr->qoffset_z; + for (int i = 0; i < 4; i++) + { + this->SRowX[i] = hdr->srow_x[i]; + this->SRowY[i] = hdr->srow_y[i]; + this->SRowZ[i] = hdr->srow_z[i]; + } +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::GetHeader(nifti_2_header *hdr) +{ + hdr->sizeof_hdr = NIFTI2HeaderSize; + memcpy(hdr->magic, this->Magic, sizeof(hdr->magic)); + hdr->datatype = static_cast<short>(this->DataType); + hdr->bitpix = static_cast<short>(this->BitPix); + for (int i = 0; i < 8; i++) + { + hdr->dim[i] = static_cast<short>(this->Dim[i]); + hdr->pixdim[i] = vtkNIFTINormalizeDouble(this->PixDim[i]); + } + hdr->intent_p1 = vtkNIFTINormalizeDouble(this->IntentP1); + hdr->intent_p2 = vtkNIFTINormalizeDouble(this->IntentP2); + hdr->intent_p3 = vtkNIFTINormalizeDouble(this->IntentP3); + hdr->vox_offset = this->VoxOffset; + hdr->scl_slope = vtkNIFTINormalizeDouble(this->SclSlope); + hdr->scl_inter = vtkNIFTINormalizeDouble(this->SclInter); + hdr->cal_min = vtkNIFTINormalizeDouble(this->CalMin); + hdr->cal_max = vtkNIFTINormalizeDouble(this->CalMax); + hdr->slice_duration = vtkNIFTINormalizeDouble(this->SliceDuration); + hdr->toffset = vtkNIFTINormalizeDouble(this->TOffset); + hdr->slice_start = this->SliceStart; + hdr->slice_end = this->SliceEnd; + strncpy(hdr->descrip, this->Descrip, sizeof(hdr->descrip)); + strncpy(hdr->aux_file, this->AuxFile, sizeof(hdr->aux_file)); + hdr->qform_code = static_cast<short>(this->QFormCode); + hdr->sform_code = static_cast<short>(this->SFormCode); + hdr->quatern_b = vtkNIFTINormalizeDouble(this->QuaternB); + hdr->quatern_c = vtkNIFTINormalizeDouble(this->QuaternC); + hdr->quatern_d = vtkNIFTINormalizeDouble(this->QuaternD); + hdr->qoffset_x = vtkNIFTINormalizeDouble(this->QOffsetX); + hdr->qoffset_y = vtkNIFTINormalizeDouble(this->QOffsetY); + hdr->qoffset_z = vtkNIFTINormalizeDouble(this->QOffsetZ); + for (int i = 0; i < 4; i++) + { + hdr->srow_x[i] = vtkNIFTINormalizeDouble(this->SRowX[i]); + hdr->srow_y[i] = vtkNIFTINormalizeDouble(this->SRowY[i]); + hdr->srow_z[i] = vtkNIFTINormalizeDouble(this->SRowZ[i]); + } + hdr->slice_code = this->SliceCode; + hdr->xyzt_units = this->XYZTUnits; + hdr->intent_code = static_cast<short>(this->IntentCode); + strncpy(hdr->intent_name, this->IntentName, sizeof(hdr->intent_name)); + hdr->dim_info = static_cast<char>(this->DimInfo); + memset(hdr->unused_str, '\0', 15); +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::DeepCopy(vtkNIFTIImageHeader *o) +{ + if (o) + { + nifti_2_header hdr; + o->GetHeader(&hdr); + this->SetHeader(&hdr); + } + else + { + this->Initialize(); + } +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::PrintSelf(ostream& os, vtkIndent indent) +{ + this->Superclass::PrintSelf(os, indent); + + os.setf(std::ios::hex, std::ios::basefield); + os << indent << "DimInfo: 0x" << this->DimInfo << "\n"; + os.unsetf(std::ios::hex); + os << indent << "Dim:"; + for (int i = 0; i < 8; i++) + { + os << " " << this->Dim[i]; + } + os << indent << "\n"; + os << indent << "PixDim:"; + for (int i = 0; i < 8; i++) + { + os << " " << this->PixDim[i]; + } + os << indent << "\n"; + os << indent << "VoxOffset:" << this->VoxOffset << "\n"; + os << indent << "IntentP1: " << this->IntentP1 << "\n"; + os << indent << "IntentP2: " << this->IntentP2 << "\n"; + os << indent << "IntentP3: " << this->IntentP3 << "\n"; + os << indent << "IntentCode: " << this->IntentCode << "\n"; + os << indent << "DataType: " << this->DataType << "\n"; + os << indent << "BitPix: " << this->BitPix << "\n"; + os << indent << "SliceStart: " << this->SliceStart << "\n"; + os << indent << "SclSlope: " << this->SclSlope << "\n"; + os << indent << "SclInter: " << this->SclInter << "\n"; + os << indent << "SliceEnd: " << this->SliceEnd << "\n"; + os << indent << "SliceCode: " << static_cast<int>(this->SliceCode) << "\n"; + os.setf(std::ios::hex, std::ios::basefield); + os << indent << "XYZTUnits: 0x" << static_cast<int>(this->XYZTUnits) << "\n"; + os.unsetf(std::ios::hex); + os << indent << "CalMax: " << this->CalMax << "\n"; + os << indent << "CalMin: " << this->CalMin << "\n"; + os << indent << "SliceDuration: " << this->SliceDuration << "\n"; + os << indent << "TOffset: " << this->TOffset << "\n"; + os << indent << "Descrip: \""; + for (size_t j = 0; j < 80 && this->Descrip[j] != '\0'; j++) + { + os << (isprint(this->Descrip[j]) ? this->Descrip[j] : '?'); + } + os << "\"\n"; + os << indent << "AuxFile: \""; + for (size_t j = 0; j < 24 && this->AuxFile[j] != '\0'; j++) + { + os << (isprint(this->AuxFile[j]) ? this->AuxFile[j] : '?'); + } + os << "\"\n"; + os << indent << "QFormCode: " << this->QFormCode << "\n"; + os << indent << "SFormCode: " << this->SFormCode << "\n"; + os << indent << "QuaternB: " << this->QuaternB << "\n"; + os << indent << "QuaternC: " << this->QuaternC << "\n"; + os << indent << "QuaternD: " << this->QuaternD << "\n"; + os << indent << "QOffsetX: " << this->QOffsetX << "\n"; + os << indent << "QOffsetY: " << this->QOffsetY << "\n"; + os << indent << "QOffsetZ: " << this->QOffsetZ << "\n"; + os << indent << "SRowX:"; + for (int i = 0; i < 4; i++) + { + os << " " << this->SRowX[i]; + } + os << "\n"; + os << indent << "SRowY:"; + for (int i = 0; i < 4; i++) + { + os << " " << this->SRowY[i]; + } + os << "\n"; + os << indent << "SRowZ:"; + for (int i = 0; i < 4; i++) + { + os << " " << this->SRowZ[i]; + } + os << "\n"; + os << indent << "IntentName: \""; + for (size_t j = 0; j < 16 && this->IntentName[j] != '\0'; j++) + { + os << (isprint(this->IntentName[j]) ? this->IntentName[j] : '?'); + } + os << "\"\n"; + os << indent << "Magic: \""; + for (size_t j = 0; j < 4 && this->Magic[j] != '\0'; j++) + { + os << (isprint(this->Magic[j]) ? this->Magic[j] : '?'); + } + os << "\"\n"; +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::SetStringValue(char *x, const char *y, size_t n) +{ + if (y == 0) + { + y = ""; + } + if (strncmp(x, y, n) != 0) + { + strncpy(x, y, n); + x[n] = '\0'; + this->Modified(); + } +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::SetIntentName(const char *val) +{ + this->SetStringValue(this->IntentName, val, 16); +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::SetDescrip(const char *val) +{ + this->SetStringValue(this->Descrip, val, 80); +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageHeader::SetAuxFile(const char *val) +{ + this->SetStringValue(this->AuxFile, val, 24); +} diff --git a/IO/Image/vtkNIFTIImageHeader.h b/IO/Image/vtkNIFTIImageHeader.h new file mode 100644 index 0000000000000000000000000000000000000000..c7a7f85d1466bc23f12afabecdc1e61398aae25f --- /dev/null +++ b/IO/Image/vtkNIFTIImageHeader.h @@ -0,0 +1,383 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: vtkNIFTIImageHeader.h + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ +// .NAME vtkNIFTIImageHeader - Store NIfTI header information. +// .SECTION Description +// This class stores the header of a NIfTI file in a VTK-friendly format. +// By using this class, it is possible to specify the header information +// that will be stored in a file written by the vtkNIFTIImageWriter. Note +// that the SForm and QForm orientation information in this class will be +// ignored by the writer if an SForm and QForm have been explicitly set +// via the writer's SetSForm and SetQForm methods. Also note that all +// info like Dim, PixDim, DataType, etc. will be ignored by the writer +// because this information must instead be taken from the vtkImageData +// information. Finally, note that the vtkNIFTIImageWriter will ignore the +// Descrip field, since it has its own SetDescription method. +// .SECTION Thanks +// This class was contributed to VTK by the Calgary Image Processing and +// Analysis Centre (CIPAC). +// .SECTION See Also +// vtkNIFTIImageReader, vtkNIFTIImageWriter + +#ifndef __vtkNIFTIImageHeader_h +#define __vtkNIFTIImageHeader_h + +#include "vtkIOImageModule.h" // For export macro +#include "vtkObject.h" + +struct nifti_1_header; +struct nifti_2_header; + +//---------------------------------------------------------------------------- +class VTKIOIMAGE_EXPORT vtkNIFTIImageHeader : public vtkObject +{ +public: + + // Description: + // NIFTI intent codes. + enum IntentCodeEnum { + IntentNone = 0, + IntentCorrel = 2, + IntentTTest = 3, + IntentFTest = 4, + IntentZScore = 5, + IntentChiSQ = 6, + IntentBeta = 7, + IntentBinom = 8, + IntentGamma = 9, + IntentPoisson = 10, + IntentNormal = 11, + IntentFTestNonc = 12, + IntentChiSQNonc = 13, + IntentLogistic = 14, + IntentLaplace = 15, + IntentUniform = 16, + IntentTTestNonc = 17, + IntentWeibull = 18, + IntentChi = 19, + IntentInvGauss = 20, + IntentExtVal = 21, + IntentPVal = 22, + IntentLogPVal = 23, + IntentLog10PVal = 24, + IntentEstimate = 1001, + IntentLabel = 1002, + IntentNeuroName = 1003, + IntentGenMatrix = 1004, + IntentSymMatrix = 1005, + IntentDispVect = 1006, + IntentVector = 1007, + IntentPointSet = 1008, + IntentTriangle = 1009, + IntentQuaternion = 1010, + IntentDimless = 1011, + IntentTimeSeries = 2001, + IntentNodeIndex = 2002, + IntentRGBVector = 2003, + IntentRGBAVector = 2004, + IntentShape = 2005 + }; + + // Description: + // NIFTI transform codes. + enum XFormCodeEnum { + XFormUnkown = 0, + XFormScannerAnat = 1, + XFormAlignedAnat = 2, + XFormTalairach = 3, + XFormMNI152 = 4 + }; + + // Description: + // NIFTI slice codes. + enum SliceCodeEnum { + SliceUnknown = 0, + SliceSeqInc = 1, + SliceSeqDec = 2, + SliceAltInc = 3, + SliceAltDec = 4, + SliceAltInc2 = 5, + SliceAltDec2 = 6 + }; + + // Description: + // NIFTI unit codes. + enum UnitsXYZTEnum { + UnitsUnknown = 0, + UnitsMeter = 1, + UnitsMM = 2, + UnitsMicron = 3, + UnitsSpace = 7, + UnitsSec = 8, + UnitsMSec = 16, + UnitsUSec = 24, + UnitsHz = 32, + UnitsPPM = 40, + UnitsRads = 48, + UnitsTime = 56 + }; + + // Description: + // NIFTI data types. + // Types RGB24 and RGB32 are represented in VTK as a multi-component + // unsigned char image. Complex values are represented as two-component + // images. The NIFTI types Float128 and Complex256 are not supported. + enum DataTypeEnum { + TypeUInt8 = 2, + TypeInt16 = 4, + TypeInt32 = 8, + TypeFloat32 = 16, + TypeComplex64 = 32, + TypeFloat64 = 64, + TypeRGB24 = 128, + TypeInt8 = 256, + TypeUInt16 = 512, + TypeUInt32 = 768, + TypeInt64 = 1024, + TypeUInt64 = 1280, + TypeFloat128 = 1536, + TypeComplex128 = 1792, + TypeComplex256 = 2048, + TypeRGBA32 = 2304 + }; + + // Description: + // NIFTI header sizes. + enum HeaderSizeEnum { + NIFTI1HeaderSize = 348, + NIFTI2HeaderSize = 540 + }; + + // Description: + // Static method for construction. + static vtkNIFTIImageHeader *New(); + vtkTypeMacro(vtkNIFTIImageHeader, vtkObject); + + // Description: + // Print information about this object. + virtual void PrintSelf(ostream& os, vtkIndent indent); + + // Description: + // Get the magic number for the NIFTI file as a null-terminated string. + const char *GetMagic() { return this->Magic; } + + // Description: + // Get the offset to the pixel data within the file. + vtkTypeInt64 GetVoxOffset() { return this->VoxOffset; } + + // Description: + // Get the data type. + int GetDataType() { return this->DataType; } + + // Description: + // Get the number of bits per pixel. + int GetBitPix() { return this->BitPix; } + + // Description: + // Get the nth dimension of the data, where GetDim(0) returns the + // number of dimensions that are defined for the file. + vtkTypeInt64 GetDim(int i) { + return (i < 0 || i > 7 ? 0 : this->Dim[i]); } + + // Description: + // Get the sample spacing in the nth dimension. If GetPixDim(0) is + // negative, then the quaternion for the qform describes the correct + // orientation only after the slice ordering has been reversed. + double GetPixDim(int i) { + return (i < 0 || i > 7 ? 0.0 : this->PixDim[i]); } + + // Description: + // Get the NIFTI intent code. This is an enumerated value in the NIFTI + // header that states what the data is intended to be used for. + vtkSetMacro(IntentCode, int); + int GetIntentCode() { return this->IntentCode; } + + // Description: + // Get the intent name. This should match the intent code. + void SetIntentName(const char *name); + const char *GetIntentName() { return this->IntentName; } + + // Description: + // Get one of the NIFTI intent parameters. The definition of these + // parameters varies according to the IntentCode. + vtkSetMacro(IntentP1, double); + double GetIntentP1() { return this->IntentP1; } + vtkSetMacro(IntentP2, double); + double GetIntentP2() { return this->IntentP2; } + vtkSetMacro(IntentP3, double); + double GetIntentP3() { return this->IntentP3; } + + // Description: + // Get the scale and slope to apply to the data in order to get + // the real-valued data values. + vtkSetMacro(SclSlope, double); + double GetSclSlope() { return this->SclSlope; } + vtkSetMacro(SclInter, double); + double GetSclInter() { return this->SclInter; } + + // Description: + // Get the calibrated range of the data, i.e. the values stored in the + // cal_min and cal_max fields in the header. + vtkSetMacro(CalMin, double); + double GetCalMin() { return this->CalMin; } + vtkSetMacro(CalMax, double); + double GetCalMax() { return this->CalMax; } + + // Description: + // Get the slice_duration and toffset from the header. + vtkSetMacro(SliceDuration, double); + double GetSliceDuration() { return this->SliceDuration; } + vtkSetMacro(TOffset, double); + double GetTOffset() { return this->TOffset; } + + // Description: + // Get the slice range for the data. + vtkSetMacro(SliceStart, vtkTypeInt64); + vtkTypeInt64 GetSliceStart() { return this->SliceStart; } + vtkSetMacro(SliceEnd, vtkTypeInt64); + vtkTypeInt64 GetSliceEnd() { return this->SliceEnd; } + + // Description: + // Get the slice code for the data. + vtkSetMacro(SliceCode, int); + int GetSliceCode() { return this->SliceCode; } + + // Description: + // Get a bitfield that describes the units for the first 4 dims. + vtkSetMacro(XYZTUnits, int); + int GetXYZTUnits() { return this->XYZTUnits; } + + // Description: + // Get a bitfield with extra information about the dimensions, it + // states which dimensions are the phase encode, frequency encode, + // and slice encode dimensions for MRI acquisitions. + vtkSetMacro(DimInfo, int); + int GetDimInfo() { return this->DimInfo; } + + // Description: + // Get a null-terminated file descriptor, this usually gives the + // name of the software that wrote the file. It will have a maximum + // length of 80 characters. Use ASCII to ensure compatibility with + // all NIFTI software, the NIFTI standard itself does not specify what + // encodings are permitted. + void SetDescrip(const char *descrip); + const char *GetDescrip() { return this->Descrip; } + + // Description: + // Get an auxilliary file, e.g. a color table, that is associated + // with this data. The length of the filename must be a maximum of + // 24 characters, and it will be assumed to be in the same directory + // as the NIFTI file. + void SetAuxFile(const char *auxfile); + const char *GetAuxFile() { return this->AuxFile; } + + // Description: + // Get the QForm or SForm code. + vtkSetMacro(QFormCode, int); + int GetQFormCode() { return this->QFormCode; } + vtkSetMacro(SFormCode, int); + int GetSFormCode() { return this->SFormCode; } + + // Description: + // Get information about the quaternion transformation. Note that + // the vtkNIFTIImageWriter ignores this part of the header if a quaternion + // has been set via vtkNIFTIImageWriter::SetQFormMatrix(). + vtkSetMacro(QuaternB, double); + double GetQuaternB() { return this->QuaternB; } + vtkSetMacro(QuaternC, double); + double GetQuaternC() { return this->QuaternC; } + vtkSetMacro(QuaternD, double); + double GetQuaternD() { return this->QuaternD; } + vtkSetMacro(QOffsetX, double); + double GetQOffsetX() { return this->QOffsetX; } + vtkSetMacro(QOffsetY, double); + double GetQOffsetY() { return this->QOffsetY; } + vtkSetMacro(QOffsetZ, double); + double GetQOffsetZ() { return this->QOffsetZ; } + + // Description: + // Get information about the matrix transformation. Note that + // the vtkNIFTIImageWriter ignores this part of the header if a matrix + // has been set via vtkNIFTIImageWriter::SetSFormMatrix(). + vtkSetVector4Macro(SRowX, double); + vtkGetVector4Macro(SRowX, double); + vtkSetVector4Macro(SRowY, double); + vtkGetVector4Macro(SRowY, double); + vtkSetVector4Macro(SRowZ, double); + vtkGetVector4Macro(SRowZ, double); + + // Description: + // Initialize the header to default values. + void Initialize(); + + // Description: + // Make a copy of the header. + void DeepCopy(vtkNIFTIImageHeader *o); + + // Description: + // Set the values from an existing nifti struct, or store + // the values in an existing nifti struct. + void SetHeader(const nifti_1_header *hdr); + void GetHeader(nifti_1_header *hdr); + void SetHeader(const nifti_2_header *hdr); + void GetHeader(nifti_2_header *hdr); + +protected: + vtkNIFTIImageHeader(); + ~vtkNIFTIImageHeader(); + + char Magic[12]; + vtkTypeInt64 VoxOffset; + int DataType; + int BitPix; + vtkTypeInt64 Dim[8]; + double PixDim[8]; + int IntentCode; + char IntentName[18]; + double IntentP1; + double IntentP2; + double IntentP3; + double SclSlope; + double SclInter; + double CalMin; + double CalMax; + double SliceDuration; + double TOffset; + vtkTypeInt64 SliceStart; + vtkTypeInt64 SliceEnd; + int SliceCode; + int XYZTUnits; + int DimInfo; + char Descrip[82]; + char AuxFile[26]; + int QFormCode; + int SFormCode; + double QuaternB; + double QuaternC; + double QuaternD; + double QOffsetX; + double QOffsetY; + double QOffsetZ; + double SRowX[4]; + double SRowY[4]; + double SRowZ[4]; + + void SetStringValue(char *x, const char *y, size_t n); + +private: + vtkNIFTIImageHeader(const vtkNIFTIImageHeader&); // Not implemented. + void operator=(const vtkNIFTIImageHeader&); // Not implemented. +}; + +#endif // __vtkNIFTIImageHeader_h diff --git a/IO/Image/vtkNIFTIImagePrivate.h b/IO/Image/vtkNIFTIImagePrivate.h new file mode 100644 index 0000000000000000000000000000000000000000..d7074ac99c8b3af5bbaa0814d5079fabf1a6de39 --- /dev/null +++ b/IO/Image/vtkNIFTIImagePrivate.h @@ -0,0 +1,277 @@ +#ifndef _NIFTI_HEADER_ +#define _NIFTI_HEADER_ + +/***************************************************************************** + ** This file defines the "NIFTI-1" header format. ** + ** It is derived from 2 meetings at the NIH (31 Mar 2003 and ** + ** 02 Sep 2003) of the Data Format Working Group (DFWG), ** + ** chartered by the NIfTI (Neuroimaging Informatics Technology ** + ** Initiative) at the National Institutes of Health (NIH). ** + **--------------------------------------------------------------** + ** Neither the National Institutes of Health (NIH), the DFWG, ** + ** nor any of the members or employees of these institutions ** + ** imply any warranty of usefulness of this material for any ** + ** purpose, and do not assume any liability for damages, ** + ** incidental or otherwise, caused by any use of this document. ** + ** If these conditions are not acceptable, do not use this! ** + **--------------------------------------------------------------** + ** Author: Robert W Cox (NIMH, Bethesda) ** + ** Advisors: John Ashburner (FIL, London), ** + ** Stephen Smith (FMRIB, Oxford), ** + ** Mark Jenkinson (FMRIB, Oxford) ** +******************************************************************************/ + +/*=================*/ +#ifdef __cplusplus +extern "C" { +#endif +/*=================*/ + +/*! \struct nifti_1_header + \brief Data structure defining the fields in the nifti1 header. + This binary header should be found at the beginning of a valid + NIFTI-1 header file. + */ + /*************************/ /************/ +struct nifti_1_header { /* NIFTI-1 usage */ /* offset */ + /*************************/ /************/ + + int sizeof_hdr; /*!< MUST be 348 */ /* 0 */ + char data_type[10]; /*!< ++UNUSED++ */ /* 4 */ + char db_name[18]; /*!< ++UNUSED++ */ /* 14 */ + int extents; /*!< ++UNUSED++ */ /* 32 */ + short session_error; /*!< ++UNUSED++ */ /* 36 */ + char regular; /*!< ++UNUSED++ */ /* 38 */ + char dim_info; /*!< MRI slice ordering. */ /* 39 */ + short dim[8]; /*!< Data array dimensions.*/ /* 40 */ + float intent_p1; /*!< 1st intent parameter. */ /* 56 */ + float intent_p2; /*!< 2nd intent parameter. */ /* 60 */ + float intent_p3; /*!< 3rd intent parameter. */ /* 64 */ + short intent_code; /*!< NIFTI_INTENT_* code. */ /* 68 */ + short datatype; /*!< Defines data type! */ /* 70 */ + short bitpix; /*!< Number bits/voxel. */ /* 72 */ + short slice_start; /*!< First slice index. */ /* 74 */ + float pixdim[8]; /*!< Grid spacings. */ /* 76 */ + float vox_offset; /*!< Offset into .nii file */ /* 108 */ + float scl_slope; /*!< Data scaling: slope. */ /* 112 */ + float scl_inter; /*!< Data scaling: offset. */ /* 116 */ + short slice_end; /*!< Last slice index. */ /* 120 */ + char slice_code; /*!< Slice timing order. */ /* 122 */ + char xyzt_units; /*!< Units of pixdim[1..4] */ /* 123 */ + float cal_max; /*!< Max display intensity */ /* 124 */ + float cal_min; /*!< Min display intensity */ /* 128 */ + float slice_duration;/*!< Time for 1 slice. */ /* 132 */ + float toffset; /*!< Time axis shift. */ /* 136 */ + int glmax; /*!< ++UNUSED++ */ /* 140 */ + int glmin; /*!< ++UNUSED++ */ /* 144 */ + char descrip[80]; /*!< any text you like. */ /* 148 */ + char aux_file[24]; /*!< auxiliary filename. */ /* 228 */ + short qform_code; /*!< NIFTI_XFORM_* code. */ /* 252 */ + short sform_code; /*!< NIFTI_XFORM_* code. */ /* 254 */ + float quatern_b; /*!< Quaternion b param. */ /* 256 */ + float quatern_c; /*!< Quaternion c param. */ /* 260 */ + float quatern_d; /*!< Quaternion d param. */ /* 264 */ + float qoffset_x; /*!< Quaternion x shift. */ /* 268 */ + float qoffset_y; /*!< Quaternion y shift. */ /* 272 */ + float qoffset_z; /*!< Quaternion z shift. */ /* 276 */ + float srow_x[4]; /*!< 1st row affine transform. */ /* 280 */ + float srow_y[4]; /*!< 2nd row affine transform. */ /* 296 */ + float srow_z[4]; /*!< 3rd row affine transform. */ /* 312 */ + char intent_name[16];/*!< 'name' or meaning of data. */ /* 328 */ + char magic[4]; /*!< MUST be "ni1\0" or "n+1\0". */ /* 344 */ + +}; /**** 348 bytes total ****/ + +typedef struct nifti_1_header nifti_1_header; + +/*---------------------------------------------------------------------------*/ +/* TYPE OF DATA (acceptable values for datatype field): + --------------------------------------------------- + Values of datatype smaller than 256 are ANALYZE 7.5 compatible. + Larger values are NIFTI-1 additions. These are all multiples of 256, so + that no bits below position 8 are set in datatype. But there is no need + to use only powers-of-2, as the original ANALYZE 7.5 datatype codes do. + + The additional codes are intended to include a complete list of basic + scalar types, including signed and unsigned integers from 8 to 64 bits, + floats from 32 to 128 bits, and complex (float pairs) from 64 to 256 bits. + + Note that most programs will support only a few of these datatypes! + A NIFTI-1 program should fail gracefully (e.g., print a warning message) + when it encounters a dataset with a type it doesn't like. +-----------------------------------------------------------------------------*/ + +/*! \defgroup NIFTI1_DATATYPE_ALIASES + \brief aliases for the nifti1 datatype codes + @{ + */ + /*! unsigned char. */ +#define NIFTI_TYPE_UINT8 2 + /*! signed short. */ +#define NIFTI_TYPE_INT16 4 + /*! signed int. */ +#define NIFTI_TYPE_INT32 8 + /*! 32 bit float. */ +#define NIFTI_TYPE_FLOAT32 16 + /*! 64 bit complex = 2 32 bit floats. */ +#define NIFTI_TYPE_COMPLEX64 32 + /*! 64 bit float = double. */ +#define NIFTI_TYPE_FLOAT64 64 + /*! 3 8 bit bytes. */ +#define NIFTI_TYPE_RGB24 128 + /*! signed char. */ +#define NIFTI_TYPE_INT8 256 + /*! unsigned short. */ +#define NIFTI_TYPE_UINT16 512 + /*! unsigned int. */ +#define NIFTI_TYPE_UINT32 768 + /*! signed long long. */ +#define NIFTI_TYPE_INT64 1024 + /*! unsigned long long. */ +#define NIFTI_TYPE_UINT64 1280 + /*! 128 bit float = long double. */ +#define NIFTI_TYPE_FLOAT128 1536 + /*! 128 bit complex = 2 64 bit floats. */ +#define NIFTI_TYPE_COMPLEX128 1792 + /*! 256 bit complex = 2 128 bit floats */ +#define NIFTI_TYPE_COMPLEX256 2048 + /*! 4 8 bit bytes. */ +#define NIFTI_TYPE_RGBA32 2304 +/* @} */ + + +/*---------------------------------------------------------------------------*/ +/* MISCELLANEOUS C MACROS +-----------------------------------------------------------------------------*/ + +/*.................*/ +/*! Given a nifti_1_header struct, check if it has a good magic number. + Returns NIFTI version number (1..9) if magic is good, 0 if it is not. */ + +#define NIFTI_VERSION(h) \ + ( ( (h).magic[0]=='n' && (h).magic[3]=='\0' && \ + ( (h).magic[1]=='i' || (h).magic[1]=='+' ) && \ + ( (h).magic[2]>='1' && (h).magic[2]<='9' ) ) \ + ? (h).magic[2]-'0' : 0 ) + +/*.................*/ +/*! Check if a nifti_1_header struct says if the data is stored in the + same file or in a separate file. Returns 1 if the data is in the same + file as the header, 0 if it is not. */ + +#define NIFTI_ONEFILE(h) ( (h).magic[1] == '+' ) + +/*.................*/ +/*! Check if a nifti_1_header struct needs to be byte swapped. + Returns 1 if it needs to be swapped, 0 if it does not. */ + +#define NIFTI_NEEDS_SWAP(h) ( (h).dim[0] < 0 || (h).dim[0] > 7 ) + + +/*=================*/ +#ifdef __cplusplus +} +#endif +/*=================*/ + +#endif /* _NIFTI_HEADER_ */ + +#ifndef __NIFTI2_HEADER +#define __NIFTI2_HEADER + +/*---------------------------------------------------------------------------*/ +/* Changes to the header from NIFTI-1 to NIFTI-2 are intended to allow for + larger and more accurate fields. The changes are as follows: + + - short dim[8] -> int64_t dim[8] + - float intent_p1,2,3 -> double intent_p1,2,3 (3 fields) + - float pixdim[8] -> double pixdim[8] + - float vox_offset -> int64_t vox_offset + - float scl_slope -> double scl_slope + - float scl_inter -> double scl_inter + - float cal_max -> double cal_max + - float cal_min -> double cal_min + - float slice_duration -> double slice_duration + - float toffset -> double toffset + - short slice_start -> int64_t slice_start + - short slice_end -> int64_t slice_end + - char slice_code -> int32_t slice_code + - char xyzt_units -> int32_t xyzt_units + - short intent_code -> int32_t intent_code + - short qform_code -> int32_t qform_code + - short sform_code -> int32_t sform_code + - float quatern_b,c,d -> double quatern_b,c,d (3 fields) + - float srow_x,y,z[4] -> double srow_x,y,z[4] (3 fields) + - char magic[4] -> char magic[8] + - char unused_str[15] -> padding added at the end of the header + + - previously unused fields have been removed: + data_type, db_name, extents, session_error, regular, glmax, glmin + + - the field ordering has been changed +-----------------------------------------------------------------------------*/ + +/*=================*/ +#ifdef __cplusplus +extern "C" { +#endif +/*=================*/ + +/*! \struct nifti_2_header + \brief Data structure defining the fields in the nifti2 header. + This binary header should be found at the beginning of a valid + NIFTI-2 header file. + */ + + /*************************/ /************/ +struct nifti_2_header { /* NIFTI-2 usage */ /* offset */ + /*************************/ /************/ + int sizeof_hdr; /*!< MUST be 540 */ /* 0 */ + char magic[8]; /*!< MUST be valid signature. */ /* 4 */ + short datatype; /*!< Defines data type! */ /* 12 */ + short bitpix; /*!< Number bits/voxel. */ /* 14 */ + long long dim[8]; /*!< Data array dimensions.*/ /* 16 */ + double intent_p1; /*!< 1st intent parameter. */ /* 80 */ + double intent_p2; /*!< 2nd intent parameter. */ /* 88 */ + double intent_p3; /*!< 3rd intent parameter. */ /* 96 */ + double pixdim[8]; /*!< Grid spacings. */ /* 104 */ + long long vox_offset; /*!< Offset into .nii file */ /* 168 */ + double scl_slope; /*!< Data scaling: slope. */ /* 176 */ + double scl_inter; /*!< Data scaling: offset. */ /* 184 */ + double cal_max; /*!< Max display intensity */ /* 192 */ + double cal_min; /*!< Min display intensity */ /* 200 */ + double slice_duration;/*!< Time for 1 slice. */ /* 208 */ + double toffset; /*!< Time axis shift. */ /* 216 */ + long long slice_start;/*!< First slice index. */ /* 224 */ + long long slice_end; /*!< Last slice index. */ /* 232 */ + char descrip[80]; /*!< any text you like. */ /* 240 */ + char aux_file[24]; /*!< auxiliary filename. */ /* 320 */ + int qform_code; /*!< NIFTI_XFORM_* code. */ /* 344 */ + int sform_code; /*!< NIFTI_XFORM_* code. */ /* 348 */ + double quatern_b; /*!< Quaternion b param. */ /* 352 */ + double quatern_c; /*!< Quaternion c param. */ /* 360 */ + double quatern_d; /*!< Quaternion d param. */ /* 368 */ + double qoffset_x; /*!< Quaternion x shift. */ /* 376 */ + double qoffset_y; /*!< Quaternion y shift. */ /* 384 */ + double qoffset_z; /*!< Quaternion z shift. */ /* 392 */ + double srow_x[4]; /*!< 1st row affine transform. */ /* 400 */ + double srow_y[4]; /*!< 2nd row affine transform. */ /* 432 */ + double srow_z[4]; /*!< 3rd row affine transform. */ /* 464 */ + int slice_code; /*!< Slice timing order. */ /* 496 */ + int xyzt_units; /*!< Units of pixdim[1..4] */ /* 500 */ + int intent_code; /*!< NIFTI_INTENT_* code. */ /* 504 */ + char intent_name[16]; /*!< 'name' or meaning of data. */ /* 508 */ + char dim_info; /*!< MRI slice ordering. */ /* 524 */ + char unused_str[15]; /*!< unused, filled with \0 */ /* 525 */ +}; /**** 540 bytes total ****/ + +typedef struct nifti_2_header nifti_2_header; + +/*=================*/ +#ifdef __cplusplus +} +#endif +/*=================*/ + +#endif /* __NIFTI2_HEADER */ +// VTK-HeaderTest-Exclude: vtkNIFTIImagePrivate.h diff --git a/IO/Image/vtkNIFTIImageReader.cxx b/IO/Image/vtkNIFTIImageReader.cxx new file mode 100644 index 0000000000000000000000000000000000000000..d9313d1f2f55ec59b51b65e3bf9e61fcbd8836f2 --- /dev/null +++ b/IO/Image/vtkNIFTIImageReader.cxx @@ -0,0 +1,1337 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: vtkNIFTIImageReader.cxx + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ + +#include "vtkNIFTIImageReader.h" +#include "vtkObjectFactory.h" +#include "vtkImageData.h" +#include "vtkPointData.h" +#include "vtkDataArray.h" +#include "vtkByteSwap.h" +#include "vtkMatrix4x4.h" +#include "vtkMath.h" +#include "vtkCommand.h" +#include "vtkErrorCode.h" +#include "vtkInformation.h" +#include "vtkInformationVector.h" +#include "vtkStreamingDemandDrivenPipeline.h" +#include "vtkStringArray.h" +#include "vtkVersion.h" + +#include "vtksys/SystemTools.hxx" +#include "vtksys/ios/sstream" + +// Header for NIFTI +#include "vtkNIFTIImageHeader.h" +#include "vtkNIFTIImagePrivate.h" + +// Header for zlib +#include "vtk_zlib.h" + +#include <ctype.h> +#include <string.h> +#include <string> + +vtkStandardNewMacro(vtkNIFTIImageReader); + +//---------------------------------------------------------------------------- +vtkNIFTIImageReader::vtkNIFTIImageReader() +{ + for (int i = 0; i < 8; i++) + { + this->Dim[i] = 0; + } + for (int i = 0; i < 8; i++) + { + this->PixDim[i] = 1.0; + } + this->TimeAsVector = false; + this->RescaleSlope = 1.0; + this->RescaleIntercept = 0.0; + this->QFac = 1.0; + this->QFormMatrix = 0; + this->SFormMatrix = 0; + this->NIFTIHeader = 0; +} + +//---------------------------------------------------------------------------- +vtkNIFTIImageReader::~vtkNIFTIImageReader() +{ + if (this->QFormMatrix) + { + this->QFormMatrix->Delete(); + } + if (this->SFormMatrix) + { + this->SFormMatrix->Delete(); + } + if (this->NIFTIHeader) + { + this->NIFTIHeader->Delete(); + } +} + +//---------------------------------------------------------------------------- +namespace { // anonymous namespace + +void vtkNIFTIImageReaderSwapHeader(nifti_1_header *hdr) +{ + // Common to NIFTI and Analyze 7.5 + vtkByteSwap::SwapVoidRange(&hdr->sizeof_hdr, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->extents, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->session_error, 1, 2); + vtkByteSwap::SwapVoidRange(hdr->dim, 8, 2); + vtkByteSwap::SwapVoidRange(&hdr->intent_p1, 1, 4); // unused in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->intent_p2, 1, 4); // unused in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->intent_p3, 1, 4); // unused in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->intent_code, 1, 2); // unused in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->datatype, 1, 2); + vtkByteSwap::SwapVoidRange(&hdr->bitpix, 1, 2); + vtkByteSwap::SwapVoidRange(&hdr->slice_start, 1, 2); // dim_un0 in 7.5 + vtkByteSwap::SwapVoidRange(hdr->pixdim, 8, 4); + vtkByteSwap::SwapVoidRange(&hdr->vox_offset, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->scl_slope, 1, 4); // unused in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->scl_inter, 1, 4); // unused in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->slice_end, 1, 2); // unused in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->cal_max, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->cal_min, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->slice_duration,1, 4); // compressed in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->toffset, 1, 4); // verified in 7.5 + vtkByteSwap::SwapVoidRange(&hdr->glmax, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->glmin, 1, 4); + + // All NIFTI-specific (meaning is totally different in Analyze 7.5) + if (strncmp(hdr->magic, "ni1", 4) == 0 || + strncmp(hdr->magic, "n+1", 4) == 0) + { + vtkByteSwap::SwapVoidRange(&hdr->qform_code, 1, 2); + vtkByteSwap::SwapVoidRange(&hdr->sform_code, 1, 2); + vtkByteSwap::SwapVoidRange(&hdr->quatern_b, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->quatern_c, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->quatern_d, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->qoffset_x, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->qoffset_y, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->qoffset_z, 1, 4); + vtkByteSwap::SwapVoidRange(hdr->srow_x, 4, 4); + vtkByteSwap::SwapVoidRange(hdr->srow_y, 4, 4); + vtkByteSwap::SwapVoidRange(hdr->srow_z, 4, 4); + } +} + +void vtkNIFTIImageReaderSwapHeader(nifti_2_header *hdr) +{ + vtkByteSwap::SwapVoidRange(&hdr->sizeof_hdr, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->datatype, 1, 2); + vtkByteSwap::SwapVoidRange(&hdr->bitpix, 1, 2); + vtkByteSwap::SwapVoidRange(hdr->dim, 8, 8); + vtkByteSwap::SwapVoidRange(&hdr->intent_p1, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->intent_p2, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->intent_p3, 1, 8); + vtkByteSwap::SwapVoidRange(hdr->pixdim, 8, 8); + vtkByteSwap::SwapVoidRange(&hdr->vox_offset, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->scl_slope, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->scl_inter, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->cal_max, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->cal_min, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->slice_duration,1, 8); + vtkByteSwap::SwapVoidRange(&hdr->toffset, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->slice_start, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->slice_end, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->qform_code, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->sform_code, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->quatern_b, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->quatern_c, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->quatern_d, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->qoffset_x, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->qoffset_y, 1, 8); + vtkByteSwap::SwapVoidRange(&hdr->qoffset_z, 1, 8); + vtkByteSwap::SwapVoidRange(hdr->srow_x, 4, 8); + vtkByteSwap::SwapVoidRange(hdr->srow_y, 4, 8); + vtkByteSwap::SwapVoidRange(hdr->srow_z, 4, 8); + vtkByteSwap::SwapVoidRange(&hdr->slice_code, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->xyzt_units, 1, 4); + vtkByteSwap::SwapVoidRange(&hdr->intent_code, 1, 4); +} + +} // end anonymous namespace + +//---------------------------------------------------------------------------- +vtkNIFTIImageHeader *vtkNIFTIImageReader::GetNIFTIHeader() +{ + if (!this->NIFTIHeader) + { + this->NIFTIHeader = vtkNIFTIImageHeader::New(); + } + return this->NIFTIHeader; +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageReader::PrintSelf(ostream& os, vtkIndent indent) +{ + this->Superclass::PrintSelf(os, indent); + + os << indent << "TimeAsVector: " + << (this->TimeAsVector ? "On\n" : "Off\n"); + os << indent << "TimeDimension: " << this->GetTimeDimension() << "\n"; + os << indent << "TimeSpacing: " << this->GetTimeSpacing() << "\n"; + os << indent << "RescaleSlope: " << this->RescaleSlope << "\n"; + os << indent << "RescaleIntercept: " << this->RescaleIntercept << "\n"; + os << indent << "QFac: " << this->QFac << "\n"; + + os << indent << "QFormMatrix:"; + if (this->QFormMatrix) + { + double mat[16]; + vtkMatrix4x4::DeepCopy(mat, this->QFormMatrix); + for (int i = 0; i < 16; i++) + { + os << " " << mat[i]; + } + os << "\n"; + } + else + { + os << " (none)\n"; + } + + os << indent << "SFormMatrix:"; + if (this->SFormMatrix) + { + double mat[16]; + vtkMatrix4x4::DeepCopy(mat, this->SFormMatrix); + for (int i = 0; i < 16; i++) + { + os << " " << mat[i]; + } + os << "\n"; + } + else + { + os << " (none)\n"; + } + + os << indent << "NIFTIHeader:" << (this->NIFTIHeader ? "\n" : " (none)\n"); +} + +//---------------------------------------------------------------------------- +bool vtkNIFTIImageReader::CheckExtension( + const char *filename, const char *ext) +{ + if (strlen(ext) == 4 && ext[0] == '.') + { + size_t n = strlen(filename); + if (n > 2 && filename[n-3] == '.' && + tolower(filename[n-2]) == 'g' && + tolower(filename[n-1]) == 'z') + { + n -= 3; + } + if (n > 3 && filename[n-4] == '.' && + tolower(filename[n-3]) == tolower(ext[1]) && + tolower(filename[n-2]) == tolower(ext[2]) && + tolower(filename[n-1]) == tolower(ext[3])) + { + return true; + } + } + return false; +} + +//---------------------------------------------------------------------------- +char *vtkNIFTIImageReader::ReplaceExtension( + const char *filename, const char *ext1, const char *ext2) +{ + char *newname = 0; + + if (strlen(ext1) == 4 && ext1[0] == '.' && + strlen(ext2) == 4 && ext2[0] == '.') + { + size_t n = strlen(filename); + size_t m = n; + newname = new char[n+4]; + strcpy(newname, filename); + + // check for trailing .gz + if (n > 2 && filename[n-3] == '.' && + tolower(filename[n-2]) == 'g' && + tolower(filename[n-1]) == 'z') + { + m = n - 3; + } + + if (vtkNIFTIImageReader::CheckExtension(filename, ext1)) + { + // replace the extension + if (isupper(filename[m-3])) + { + newname[m-3] = toupper(ext2[1]); + newname[m-2] = toupper(ext2[2]); + newname[m-1] = toupper(ext2[3]); + } + else + { + newname[m-3] = tolower(ext2[1]); + newname[m-2] = tolower(ext2[2]); + newname[m-1] = tolower(ext2[3]); + } + } + + // existence of file + for (int i = 0; i < 2; i++) + { + if (vtksys::SystemTools::FileExists(newname)) + { + return newname; + } + if (i == 0) + { + if (m < n) + { + // try again without the ".gz" + newname[m] = '\0'; + n = m; + } + else + { + // try again with the ".gz" + newname[m] = '.'; + newname[m+1] = (isupper(newname[m-3]) ? 'G' : 'g'); + newname[m+2] = (isupper(newname[m-3]) ? 'Z' : 'z'); + newname[m+3] = '\0'; + } + } + } + + delete [] newname; + newname = 0; + } + + return newname; +} + +//---------------------------------------------------------------------------- +int vtkNIFTIImageReader::CheckNIFTIVersion(const nifti_1_header *hdr) +{ + int version = 0; + + // Check for NIFTIv2. The NIFTIv2 magic number is stored where + // the data_type appears in the NIFTIv1 header. + if (hdr->data_type[0] == 'n' && + (hdr->data_type[1] == '+' || hdr->data_type[1] == 'i') && + (hdr->data_type[2] >= '2' && hdr->data_type[2] <= '9') && + hdr->data_type[3] == '\0') + { + version = (hdr->data_type[2] - '0'); + + if (hdr->data_type[4] != '\r' || + hdr->data_type[5] != '\n' || + hdr->data_type[6] != '\032' || + hdr->data_type[7] != '\n') + { + // Indicate that file was corrupted by newline conversion + version = -version; + } + } + // Check for NIFTIv1 + else if (hdr->magic[0] == 'n' && + (hdr->magic[1] == '+' || hdr->magic[1] == 'i') && + hdr->magic[2] == '1' && + hdr->magic[3] == '\0') + { + version = 1; + } + + return version; +} + +//---------------------------------------------------------------------------- +bool vtkNIFTIImageReader::CheckAnalyzeHeader(const nifti_1_header *hdr) +{ + if (hdr->sizeof_hdr == 348 || // Analyze 7.5 header size + hdr->sizeof_hdr == 1543569408) // byte-swapped 348 + { + return true; + } + return false; +} + +//---------------------------------------------------------------------------- +int vtkNIFTIImageReader::CanReadFile(const char *filename) +{ + vtkDebugMacro("Opening NIFTI file " << filename); + + char *hdrname = vtkNIFTIImageReader::ReplaceExtension( + filename, ".img", ".hdr"); + + if (hdrname == 0) + { + return 0; + } + + // try opening file + gzFile file = gzopen(hdrname, "rb"); + + delete [] hdrname; + + if (!file) + { + return 0; + } + + // read and check the header + bool canRead = false; + nifti_1_header hdr; + int hsize = vtkNIFTIImageHeader::NIFTI1HeaderSize; // nifti_1 header size + int rsize = gzread(file, &hdr, hsize); + if (rsize == hsize) + { + int version = vtkNIFTIImageReader::CheckNIFTIVersion(&hdr); + if (version > 0) + { + // NIFTI file + canRead = true; + } + else if (version == 0) + { + // Analyze 7.5 file + canRead = vtkNIFTIImageReader::CheckAnalyzeHeader(&hdr); + } + } + + gzclose(file); + + return canRead; +} + +//---------------------------------------------------------------------------- +int vtkNIFTIImageReader::RequestInformation( + vtkInformation* vtkNotUsed(request), + vtkInformationVector** vtkNotUsed(inputVector), + vtkInformationVector* outputVector) +{ + // Clear the error indicator. + this->SetErrorCode(vtkErrorCode::NoError); + + // Create the header object + if (!this->NIFTIHeader) + { + this->NIFTIHeader = vtkNIFTIImageHeader::New(); + } + + // default byte order is native byte order +#ifdef VTK_WORDS_BIGENDIAN + bool isLittleEndian = false; +#else + bool isLittleEndian = true; +#endif + + const char *filename = 0; + char *hdrname = 0; + + if (this->FileNames) + { + vtkIdType n = this->FileNames->GetNumberOfValues(); + int headers = 0; + for (vtkIdType i = 0; i < n; i++) + { + filename = this->FileNames->GetValue(i); + // this checks for .hdr and .hdr.gz, case insensitive + if (vtkNIFTIImageReader::CheckExtension(filename, ".hdr")) + { + headers++; + hdrname = new char[strlen(filename) + 1]; + strcpy(hdrname, filename); + } + } + if (n != 2 || headers != 1) + { + vtkErrorMacro("There must be two files and one must be a .hdr file."); + return 0; + } + } + else + { + filename = this->GetFileName(); + } + + if (filename == 0) + { + vtkErrorMacro("A FileName must be provided"); + this->SetErrorCode(vtkErrorCode::NoFileNameError); + return 0; + } + + if (hdrname == 0) + { + hdrname = vtkNIFTIImageReader::ReplaceExtension( + filename, ".img", ".hdr"); + } + + if (hdrname == 0) + { + vtkErrorMacro("Unable to locate header for file " << filename); + this->SetErrorCode(vtkErrorCode::CannotOpenFileError); + return 0; + } + + vtkDebugMacro("Opening NIFTI file " << hdrname); + + // try opening file + gzFile file = gzopen(hdrname, "rb"); + + if (!file) + { + vtkErrorMacro("Cannot open file " << hdrname); + delete [] hdrname; + this->SetErrorCode(vtkErrorCode::CannotOpenFileError); + return 0; + } + + // read and check the header + bool canRead = false; + int niftiVersion = 0; + nifti_1_header *hdr1 = new nifti_1_header; + nifti_2_header hdr2obj; + nifti_2_header *hdr2 = &hdr2obj; + const int hsize = vtkNIFTIImageHeader::NIFTI1HeaderSize; + int rsize = gzread(file, hdr1, hsize); + if (rsize == hsize) + { + niftiVersion = vtkNIFTIImageReader::CheckNIFTIVersion(hdr1); + if (niftiVersion >= 2) + { + // the header was a NIFTIv2 header + const int h2size = vtkNIFTIImageHeader::NIFTI2HeaderSize; + // copy what was read into the NIFTIv1 header + memcpy(hdr2, hdr1, hsize); + // read the remainder of the NIFTIv2 header + rsize = gzread(file, reinterpret_cast<char *>(hdr2)+hsize, h2size-hsize); + if (rsize == h2size-hsize) + { + canRead = true; + } + } + else if (niftiVersion == 1) + { + // the header was a NIFTIv1 header + canRead = true; + } + else if (niftiVersion == 0) + { + // Analyze 7.5 file + canRead = vtkNIFTIImageReader::CheckAnalyzeHeader(hdr1); + } + } + + if (canRead) + { + if (niftiVersion >= 2) + { + if (NIFTI_NEEDS_SWAP(*hdr2)) + { + vtkNIFTIImageReaderSwapHeader(hdr2); + isLittleEndian = !isLittleEndian; + } + this->NIFTIHeader->SetHeader(hdr2); + } + else + { + if (NIFTI_NEEDS_SWAP(*hdr1)) + { + vtkNIFTIImageReaderSwapHeader(hdr1); + isLittleEndian = !isLittleEndian; + } + // convert NIFTIv1 header into NIFTIv2 + this->NIFTIHeader->SetHeader(hdr1); + this->NIFTIHeader->GetHeader(hdr2); + } + } + + gzclose(file); + + // delete the NIFTIv1 header, use the NIFTIv2 header + delete hdr1; + hdr1 = 0; + + if (!canRead) + { + const char *message = (niftiVersion <= -2 ? + "NIfTI header has newline corruption " : + "Bad NIfTI header in file "); + vtkErrorMacro(<< message << hdrname); + this->SetErrorCode(vtkErrorCode::UnrecognizedFileTypeError); + delete [] hdrname; + return 0; + } + + delete [] hdrname; + + // number of dimensions + int ndim = hdr2->dim[0]; + if (ndim < 0 || ndim > 7) + { + vtkErrorMacro("NIfTI image has illegal ndim of " << ndim); + this->SetErrorCode(vtkErrorCode::FileFormatError); + return 0; + } + + // sanity checks + for (int i = 0; i < 8; i++) + { + // voxel spacing cannot be zero + if (hdr2->pixdim[i] == 0) + { + hdr2->pixdim[i] = 1.0; + } + if (i > ndim) + { + // dimensions greater than ndim have size of 1 + hdr2->dim[i] = 1; + } + else if (hdr2->dim[i] < 0) + { + vtkErrorMacro("NIfTI image dimension " << i << " is negative"); + this->SetErrorCode(vtkErrorCode::FileFormatError); + return 0; + } + else if ((hdr2->dim[i] & 0x7fffffff) != hdr2->dim[i]) + { + // dimension does not fit in signed int + vtkErrorMacro("NIfTI image dimension " << i << " is larger than int32"); + this->SetErrorCode(vtkErrorCode::FileFormatError); + return 0; + } + } + + if (niftiVersion > 0) + { + // pass rescale info to user (do not rescale in the reader) + this->RescaleSlope = hdr2->scl_slope; + this->RescaleIntercept = hdr2->scl_inter; + } + else + { + // rescale information not available for Analyze 7.5 + this->RescaleSlope = 1.0; + this->RescaleIntercept = 0.0; + } + + // header might be extended, vox_offset says where data starts + this->SetHeaderSize(static_cast<unsigned long>(hdr2->vox_offset)); + + // endianness of data + if (isLittleEndian) + { + this->SetDataByteOrderToLittleEndian(); + } + else + { + this->SetDataByteOrderToBigEndian(); + } + + // NIFTI images are stored in a single file, not one file per slice + this->SetFileDimensionality(3); + + // NIFTI uses a lower-left-hand origin + this->FileLowerLeftOn(); + + // dim + this->SetDataExtent(0, hdr2->dim[1]-1, + 0, hdr2->dim[2]-1, + 0, hdr2->dim[3]-1); + + // pixdim + this->SetDataSpacing(hdr2->pixdim[1], + hdr2->pixdim[2], + hdr2->pixdim[3]); + + // offset is part of the transform, so set origin to zero + this->SetDataOrigin(0.0, 0.0, 0.0); + + // map the NIFTI type to a VTK type and number of components + static const int typeMap[][3] = { + { NIFTI_TYPE_INT8, VTK_TYPE_INT8, 1}, + { NIFTI_TYPE_UINT8, VTK_TYPE_UINT8, 1 }, + { NIFTI_TYPE_INT16, VTK_TYPE_INT16, 1 }, + { NIFTI_TYPE_UINT16, VTK_TYPE_UINT16, 1 }, + { NIFTI_TYPE_INT32, VTK_TYPE_INT32, 1 }, + { NIFTI_TYPE_UINT32, VTK_TYPE_UINT32, 1 }, + { NIFTI_TYPE_INT64, VTK_TYPE_INT64, 1 }, + { NIFTI_TYPE_UINT64, VTK_TYPE_UINT64, 1 }, + { NIFTI_TYPE_FLOAT32, VTK_TYPE_FLOAT32, 1 }, + { NIFTI_TYPE_FLOAT64, VTK_TYPE_FLOAT64, 1 }, + { NIFTI_TYPE_COMPLEX64, VTK_TYPE_FLOAT32, 2 }, + { NIFTI_TYPE_COMPLEX128, VTK_TYPE_FLOAT64, 2 }, + { NIFTI_TYPE_RGB24, VTK_TYPE_UINT8, 3 }, + { NIFTI_TYPE_RGBA32, VTK_TYPE_UINT8, 4 }, + { 0, 0, 0 } + }; + + int scalarType = 0; + int numComponents = 0; + + for (int i = 0; typeMap[2] != 0; i++) + { + if (hdr2->datatype == typeMap[i][0]) + { + scalarType = typeMap[i][1]; + numComponents = typeMap[i][2]; + break; + } + } + + // if loop finished without finding a match + if (numComponents == 0) + { + vtkErrorMacro("Unrecognized NIFTI data type: " << hdr2->datatype); + this->SetErrorCode(vtkErrorCode::FileFormatError); + return 0; + } + + // vector planes become vector components + if (ndim >= 5) + { + numComponents *= hdr2->dim[5]; + } + if (ndim >= 4 && this->TimeAsVector) + { + numComponents *= hdr2->dim[4]; + } + + this->SetDataScalarType(scalarType); + this->SetNumberOfScalarComponents(numComponents); + + // Set the output information. + vtkInformation* outInfo = outputVector->GetInformationObject(0); + vtkDataObject::SetPointDataActiveScalarInfo( + outInfo, this->DataScalarType, this->NumberOfScalarComponents); + + outInfo->Set(vtkDataObject::SPACING(), this->DataSpacing, 3); + outInfo->Set(vtkDataObject::ORIGIN(), this->DataOrigin, 3); + + outInfo->Set(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT(), + this->DataExtent, 6); + + // copy dim for when RequestData is called + for (int j = 0; j < 8; j++) + { + this->Dim[j] = hdr2->dim[j]; + this->PixDim[j] = hdr2->pixdim[j]; + } + + // === Image Orientation in NIfTI files === + // + // The vtkImageData class does not provide a way of storing image + // orientation. So when we read a NIFTI file, we should also provide + // the user with a 4x4 matrix that can transform VTK's data coordinates + // into NIFTI's intended coordinate system for the image. NIFTI defines + // these coordinate systems as: + // 1) NIFTI_XFORM_SCANNER_ANAT - coordinate system of the imaging device + // 2) NIFTI_XFORM_ALIGNED_ANAT - result of registration to another image + // 3) NIFTI_XFORM_TALAIRACH - a brain-specific coordinate system + // 4) NIFTI_XFORM_MNI_152 - a similar brain-specific coordinate system + // + // NIFTI images can store orientation in two ways: + // 1) via a quaternion (orientation and offset, i.e. rigid-body) + // 2) via a matrix (used to store e.g. the results of registration) + // + // A NIFTI file can have both a quaternion (qform) and matrix (xform) + // stored in the same file. The NIFTI documentation recommends that + // the qform be used to record the "scanner anatomical" coordinates + // and that the sform, if present, be used to define a secondary + // coordinate system, e.g. a coordinate system derived through + // registration to a template. + // + // -- Quaternion Representation -- + // + // If the "quaternion" form is used, then the following equation + // defines the transformation from voxel indices to NIFTI's world + // coordinates, where R is the rotation matrix computed from the + // quaternion components: + // + // [ x ] [ R11 R12 R13 ] [ pixdim[1] * i ] [ qoffset_x ] + // [ y ] = [ R21 R22 R23 ] [ pixdim[2] * j ] + [ qoffset_y ] + // [ z ] [ R31 R32 R33 ] [ pixdim[3] * k * qfac ] [ qoffset_z ] + // + // qfac is stored in pixdim[0], if it is equal to -1 then the slices + // are stacked in reverse: VTK will have to reorder the slices in order + // to maintain a right-handed coordinate transformation between indices + // and coordinates. + // + // Let's call VTK data coordinates X,Y,Z to distinguish them from + // the NIFTI coordinates x,y,z. The relationship between X,Y,Z and + // x,y,z is expressed by a 4x4 matrix M: + // + // [ x ] [ M11 M12 M13 M14 ] [ X ] + // [ y ] = [ M21 M22 M23 M24 ] [ Y ] + // [ z ] [ M31 M32 M33 M34 ] [ Z ] + // [ 1 ] [ 0 0 0 1 ] [ 1 ] + // + // where the VTK data coordinates X,Y,Z are related to the + // VTK structured coordinates IJK (i.e. point indices) by: + // + // X = I*Spacing[0] + Origin[0] + // Y = J*Spacing[1] + Origin[1] + // Z = K*Spacing[2] + Origin[2] + // + // Now let's consider: when we read a NIFTI image, how should we set + // the Spacing, the Origin, and the matrix M? Let's consider the + // cases: + // + // 1) If there is no qform, then R is identity and qoffset is zero, + // and qfac will be 1 (never -1). So: + // I,J,K = i,j,k, Spacing = pixdim, Origin = 0, M = Identity + // + // 2) If there is a qform, and qfac is 1, then: + // + // I,J,K = i,j,k (i.e. voxel order in VTK same as in NIFTI) + // + // Spacing[0] = pixdim[1] + // Spacing[1] = pixdim[2] + // Spacing[2] = pixdim[3] + // + // Origin[0] = 0.0 + // Origin[1] = 0.0 + // Origin[2] = 0.0 + // + // [ R11 R12 R13 qoffset_x ] + // M = [ R21 R22 R23 qoffset_y ] + // [ R31 R32 R33 qoffset_z ] + // [ 0 0 0 1 ] + // + // Note that we cannot store qoffset in the origin. That would + // be mathematically incorrect. It would only give us the right + // offset when R is the identity matrix. + // + // 3) If there is a qform and qfac is -1, then the situation is more + // compilcated. We have three choices, each of which is a compromise: + // a) we can use Spacing[2] = qfac*pixdim[3], i.e. use a negative + // slice spacing, which might cause some VTK algorithms to + // misbehave (the VTK tests only use images with positive spacing). + // b) we can use M13 = -R13, M23 = -R23, M33 = -R33 i.e. introduce + // a flip into the matrix, which is very bad for VTK rendering + // algorithms and should definitely be avoided. + // c) we can reverse the order of the slices in VTK relative to + // NIFTI, which allows us to preserve positive spacing and retain + // a well-behaved rotation matrix, by using these equations: + // + // J = number_of_slices - j - 1 + // + // M14 = qoffset_x - (number_of_slices - 1)*pixdim[3]*R13 + // M24 = qoffset_y - (number_of_slices - 1)*pixdim[3]*R23 + // M34 = qoffset_z - (number_of_slices - 1)*pixdim[3]*R33 + // + // This will give us data that will be well-behaved in VTK, at + // the expense of making VTK slice numbers not match with + // the original NIFTI slice numbers. NIFTY slice 0 will become + // VTK slice N-1, and the order will be reversed. + // + // -- Matrix Representation -- + // + // If the "matrix" form is used, then pixdim[] is ignored, and the + // voxel spacing is implicitly stored in the matrix. In addition, + // the matrix may have a negative determinant, there is no "qfac" + // flip-factor as there is in the quaternion representation. + // + // Let S be the matrix stored in the NIFTI header, and let M be our + // desired coordinate tranformation from VTK data coordinates X,Y,Z + // to NIFTI data coordinates x,y,z (see discussion above for more + // information). Let's consider the cases where the determinant + // is positive, or negative. + // + // 1) If the determinant is positive, we will factor the spacing + // (but not the origin) out of the matrix. + // + // Spacing[0] = pixdim[1] + // Spacing[1] = pixdim[2] + // Spacing[2] = pixdim[3] + // + // Origin[0] = 0.0 + // Origin[1] = 0.0 + // Origin[2] = 0.0 + // + // [ S11/pixdim[1] S12/pixdim[2] S13/pixdim[3] S14 ] + // M = [ S21/pixdim[1] S22/pixdim[2] S23/pixdim[3] S24 ] + // [ S31/pixdim[1] S32/pixdim[2] S33/pixdim[3] S34 ] + // [ 0 0 0 1 ] + // + // 2) If the determinant is negative, then we face the same choices + // as when qfac is -1 for the quaternion transformation. We can: + // a) use a negative Z spacing and multiply the 3rd column of M by -1 + // b) keep the matrix as is (with a negative determinant) + // c) reorder the slices, multiply the 3rd column by -1, and adjust + // the 4th column of the matrix: + // + // M14 = S14 - (number_of_slices - 1)*S13 + // M24 = S24 - (number_of_slices - 1)*S23 + // M34 = S34 - (number_of_slices - 1)*S33 + // + // The third choice will provide a VTK image that has positive + // spacing and a matrix with a positive determinant. + // + // -- Analyze 7.5 Orientation -- + // + // This reader provides only bare-bones backwards compatibility with + // the Analyze 7.5 file header. We do not orient these files. + + // Initialize + this->QFac = 1.0; + if (this->QFormMatrix) + { + this->QFormMatrix->Delete(); + this->QFormMatrix = NULL; + } + if (this->SFormMatrix) + { + this->SFormMatrix->Delete(); + this->SFormMatrix = NULL; + } + + // Set the QFormMatrix from the quaternion data in the header. + // See the long discussion above for more information. + if (niftiVersion > 0 && hdr2->qform_code > 0) + { + double mmat[16]; + double rmat[3][3]; + double quat[4]; + + quat[1] = hdr2->quatern_b; + quat[2] = hdr2->quatern_c; + quat[3] = hdr2->quatern_d; + + quat[0] = 1.0 - quat[1]*quat[1] - quat[2]*quat[2] - quat[3]*quat[3]; + if (quat[0] > 0.0) + { + quat[0] = sqrt(quat[0]); + } + else + { + quat[0] = 0.0; + } + + vtkMath::QuaternionToMatrix3x3(quat, rmat); + + // If any matrix values are close to zero, then they should actually + // be zero but aren't due to limited numerical precision in the + // quaternion-to-matrix conversion. + const double tol = 2.384185791015625e-07; // 2**-22 + for (int i = 0; i < 3; i++) + { + for (int j = 0; j < 3; j++) + { + if (fabs(rmat[i][j]) < tol) + { + rmat[i][j] = 0.0; + } + } + vtkMath::Normalize(rmat[i]); + } + + // first row + mmat[0] = rmat[0][0]; + mmat[1] = rmat[0][1]; + mmat[2] = rmat[0][2]; + mmat[3] = hdr2->qoffset_x; + + // second row + mmat[4] = rmat[1][0]; + mmat[5] = rmat[1][1]; + mmat[6] = rmat[1][2]; + mmat[7] = hdr2->qoffset_y; + + // third row + mmat[8] = rmat[2][0]; + mmat[9] = rmat[2][1]; + mmat[10] = rmat[2][2]; + mmat[11] = hdr2->qoffset_z; + + mmat[12] = 0.0; + mmat[13] = 0.0; + mmat[14] = 0.0; + mmat[15] = 1.0; + + this->QFac = ((hdr2->pixdim[0] < 0) ? -1.0 : 1.0); + + if (this->QFac < 0) + { + // We will be reversing the order of the slices, so the first VTK + // slice will be at the position of the last NIfTI slice, and we + // must adjust the offset to compensate for this. + mmat[3] -= rmat[0][2]*hdr2->pixdim[3]*(hdr2->dim[3] - 1); + mmat[7] -= rmat[1][2]*hdr2->pixdim[3]*(hdr2->dim[3] - 1); + mmat[11] -= rmat[2][2]*hdr2->pixdim[3]*(hdr2->dim[3] - 1); + } + + this->QFormMatrix = vtkMatrix4x4::New(); + this->QFormMatrix->DeepCopy(mmat); + } + + // Set the SFormMatrix from the matrix information in the header. + // See the long discussion above for more information. + if (niftiVersion > 0 && hdr2->sform_code > 0) + { + double mmat[16]; + + // first row + mmat[0] = hdr2->srow_x[0]/hdr2->pixdim[1]; + mmat[1] = hdr2->srow_x[1]/hdr2->pixdim[2]; + mmat[2] = hdr2->srow_x[2]/hdr2->pixdim[3]; + mmat[3] = hdr2->srow_x[3]; + + // second row + mmat[4] = hdr2->srow_y[0]/hdr2->pixdim[1]; + mmat[5] = hdr2->srow_y[1]/hdr2->pixdim[2]; + mmat[6] = hdr2->srow_y[2]/hdr2->pixdim[3]; + mmat[7] = hdr2->srow_y[3]; + + // third row + mmat[8] = hdr2->srow_z[0]/hdr2->pixdim[1]; + mmat[9] = hdr2->srow_z[1]/hdr2->pixdim[2]; + mmat[10] = hdr2->srow_z[2]/hdr2->pixdim[3]; + mmat[11] = hdr2->srow_z[3]; + + mmat[12] = 0.0; + mmat[13] = 0.0; + mmat[14] = 0.0; + mmat[15] = 1.0; + + // Set QFac to -1 if the determinant is negative, unless QFac + // has already been set by the qform information. + if (vtkMatrix4x4::Determinant(mmat) < 0 && hdr2->qform_code == 0) + { + this->QFac = -1.0; + } + + if (this->QFac < 0) + { + // If QFac is set to -1 then the slices will be reversed, and we must + // reverse the slice orientation vector (the third column of the matrix) + // to compensate. + + // reverse the slice orientation vector + mmat[2] = -mmat[2]; + mmat[6] = -mmat[6]; + mmat[10] = -mmat[10]; + + // adjust the offset to compensate for changed slice ordering + mmat[3] -= hdr2->srow_x[2]*(hdr2->dim[3] - 1); + mmat[7] -= hdr2->srow_y[2]*(hdr2->dim[3] - 1); + mmat[11] -= hdr2->srow_z[2]*(hdr2->dim[3] - 1); + } + + this->SFormMatrix = vtkMatrix4x4::New(); + this->SFormMatrix->DeepCopy(mmat); + + if (this->SFormMatrix->Determinant() < 0) + { + vtkWarningMacro("SFormMatrix is flipped compared to QFormMatrix"); + } + } + + return 1; +} + +//---------------------------------------------------------------------------- +int vtkNIFTIImageReader::RequestData( + vtkInformation* request, + vtkInformationVector** vtkNotUsed(inputVector), + vtkInformationVector* outputVector) +{ + // check whether the reader is in an error state + if (this->GetErrorCode() != vtkErrorCode::NoError) + { + return 0; + } + + // which output port did the request come from + int outputPort = + request->Get(vtkDemandDrivenPipeline::FROM_OUTPUT_PORT()); + + // for now, this reader has only one output + if (outputPort > 0) + { + return 1; + } + + vtkInformation* outInfo = outputVector->GetInformationObject(0); + + int extent[6]; + outInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_EXTENT(), extent); + + // get the data object, allocate memory + vtkImageData *data = + static_cast<vtkImageData *>(outInfo->Get(vtkDataObject::DATA_OBJECT())); +#if VTK_MAJOR_VERSION >= 6 + this->AllocateOutputData(data, outInfo, extent); +#else + this->AllocateOutputData(data, extent); +#endif + + data->GetPointData()->GetScalars()->SetName("NIFTI"); + + const char *filename = 0; + char *imgname = 0; + + if (this->FileNames) + { + vtkIdType n = this->FileNames->GetNumberOfValues(); + int headers = 0; + for (vtkIdType i = 0; i < n; i++) + { + filename = this->FileNames->GetValue(i); + // this checks for .hdr and .hdr.gz, case insensitive + if (vtkNIFTIImageReader::CheckExtension(filename, ".hdr")) + { + headers++; + } + else + { + imgname = new char[strlen(filename) + 1]; + strcpy(imgname, filename); + } + } + if (n != 2 || headers != 1) + { + vtkErrorMacro("There must be two files and one must be a .hdr file."); + return 0; + } + } + else + { + filename = this->GetFileName(); + } + + if (filename == 0) + { + vtkErrorMacro("A FileName must be provided"); + return 0; + } + + if (imgname == 0) + { + imgname = vtkNIFTIImageReader::ReplaceExtension(filename, ".hdr", ".img"); + } + + if (imgname == 0) + { + vtkErrorMacro("Unable to locate image for file " << filename); + return 0; + } + + vtkDebugMacro("Opening NIFTI file " << imgname); + + data->GetPointData()->GetScalars()->SetName("NIFTI"); + + unsigned char *dataPtr = + static_cast<unsigned char *>(data->GetScalarPointer()); + + gzFile file = gzopen(imgname, "rb"); + + delete [] imgname; + + if (!file) + { + return 0; + } + + int swapBytes = this->GetSwapBytes(); + int scalarSize = data->GetScalarSize(); + int numComponents = data->GetNumberOfScalarComponents(); + int timeDim = (this->Dim[0] >= 4 ? this->Dim[4] : 1); + int vectorDim = (this->Dim[0] >= 5 ? this->Dim[5] : 1); + if (this->TimeAsVector) + { + vectorDim *= timeDim; + } + + int outSizeX = extent[1] - extent[0] + 1; + int outSizeY = extent[3] - extent[2] + 1; + int outSizeZ = extent[5] - extent[4] + 1; + + z_off_t fileVoxelIncr = scalarSize*numComponents/vectorDim; + z_off_t fileRowIncr = fileVoxelIncr*this->Dim[1]; + z_off_t fileSliceIncr = fileRowIncr*this->Dim[2]; + z_off_t fileTimeIncr = fileSliceIncr*this->Dim[3]; + z_off_t fileVectorIncr = fileTimeIncr*this->Dim[4]; + if (this->TimeAsVector) + { + fileVectorIncr = fileTimeIncr; + } + + // add a buffer for planar-vector to packed-vector conversion + unsigned char *rowBuffer = 0; + if (vectorDim > 1) + { + rowBuffer = new unsigned char[outSizeX*fileVoxelIncr]; + } + + // special increment to reverse the slices if needed + vtkIdType sliceOffset = 0; + + if (this->GetQFac() < 0) + { + // put slices in reverse order + sliceOffset = scalarSize*numComponents; + sliceOffset *= outSizeX; + sliceOffset *= outSizeY; + dataPtr += sliceOffset*(outSizeZ - 1); + } + + // report progress every 2% of the way to completion + this->InvokeEvent(vtkCommand::StartEvent); + this->UpdateProgress(0.0); + vtkIdType target = + static_cast<vtkIdType>(0.02*outSizeY*outSizeZ*vectorDim) + 1; + vtkIdType count = 0; + + // seek to the start of the data + z_off_t offset = static_cast<z_off_t>(this->GetHeaderSize()); + offset += extent[0]*fileVoxelIncr; + offset += extent[2]*fileRowIncr; + offset += extent[4]*fileSliceIncr; + + // read the data one row at a time, do planar-to-packed conversion + // of vector components if NIFTI file has a vector dimension + int rowSize = numComponents/vectorDim*outSizeX; + int t = 0; // counter for time + int c = 0; // counter for vector components + int j = 0; // counter for rows + int k = 0; // counter for slices + unsigned char *ptr = dataPtr; + + int errorCode = 0; + + while (!this->AbortExecute) + { + if (offset) + { + int rval = gzseek(file, offset, SEEK_CUR); + if (rval == -1) + { + errorCode = vtkErrorCode::FileFormatError; + if (gzeof(file)) + { + errorCode = vtkErrorCode::PrematureEndOfFileError; + } + break; + } + } + + if (vectorDim == 1) + { + // read directly into the output instead of into a buffer + rowBuffer = ptr; + } + + int code = gzread(file, rowBuffer, rowSize*scalarSize); + if (code != rowSize*scalarSize) + { + errorCode = vtkErrorCode::FileFormatError; + if (gzeof(file)) + { + errorCode = vtkErrorCode::PrematureEndOfFileError; + } + break; + } + + if (swapBytes != 0 && scalarSize > 1) + { + vtkByteSwap::SwapVoidRange(rowBuffer, rowSize, scalarSize); + } + + if (vectorDim == 1) + { + // advance the pointer to the next row + ptr += outSizeX*numComponents*scalarSize; + rowBuffer = 0; + } + else + { + // write vector plane to packed vector component + unsigned char *tmpPtr = rowBuffer; + z_off_t skipOther = scalarSize*numComponents - fileVoxelIncr; + for (int i = 0; i < outSizeX; i++) + { + // write one vector component of one voxel + z_off_t n = fileVoxelIncr; + do { *ptr++ = *tmpPtr++; } while (--n); + // skip past the other components + ptr += skipOther; + } + } + + if (++count % target == 0) + { + this->UpdateProgress(0.02*count/target); + } + + // offset to skip unread sections of the file, for when + // the update extent is less than the whole extent + offset = fileRowIncr - outSizeX*fileVoxelIncr; + if (++j == outSizeY) + { + j = 0; + offset += fileSliceIncr - outSizeY*fileRowIncr; + ptr -= 2*sliceOffset; // for reverse slice order + if (++k == outSizeZ) + { + k = 0; + offset += fileVectorIncr - outSizeZ*fileSliceIncr; + if (++t == timeDim) + { + t = 0; + } + if (++c == vectorDim) + { + break; + } + // back up the ptr to the beginning of the image, + // then increment to the next vector component + ptr = dataPtr + c*fileVoxelIncr; + + if (this->TimeAsVector) + { + // if timeDim is included in the vectorDim (and hence in the + // VTK scalar components) then we have to make sure that + // the vector components are packed before the time steps + ptr = dataPtr + (c + t*(vectorDim - 1))/timeDim*fileVoxelIncr; + } + } + } + } + + if (vectorDim > 1) + { + delete [] rowBuffer; + } + + gzclose(file); + + if (errorCode) + { + const char *errorText = "Error in NIFTI file, cannot read."; + if (errorCode == vtkErrorCode::PrematureEndOfFileError) + { + errorText = "NIFTI file is truncated, some data is missing."; + } + this->SetErrorCode(errorCode); + vtkErrorMacro(<< errorText); + return 0; + } + + this->UpdateProgress(1.0); + this->InvokeEvent(vtkCommand::EndEvent); + + return 1; +} diff --git a/IO/Image/vtkNIFTIImageReader.h b/IO/Image/vtkNIFTIImageReader.h new file mode 100644 index 0000000000000000000000000000000000000000..9aaede4bfa4c2ef5495ca5fd1a64accb5188804b --- /dev/null +++ b/IO/Image/vtkNIFTIImageReader.h @@ -0,0 +1,212 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: vtkNIFTIImageReader.h + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ +// .NAME vtkNIFTIImageReader - Read NIfTI-1 and NIfTI-2 medical image files +// .SECTION Description +// This class reads NIFTI files, either in .nii format or as separate +// .img and .hdr files. If two files are used, then they can be passed +// by using SetFileNames() instead of SetFileName(). Files ending in .gz +// are decompressed on-the-fly while they are being read. Files with +// complex numbers or vector dimensions will be read as multi-component +// images. If a NIFTI file has a time dimension, then by default only the +// first image in the time series will be read, but the TimeAsVector +// flag can be set to read the time steps as vector components. Files in +// Analyze 7.5 format are also supported by this reader. +// .SECTION Thanks +// This class was contributed to VTK by the Calgary Image Processing and +// Analysis Centre (CIPAC). +// .SECTION See Also +// vtkNIFTIImageWriter, vtkNIFTIImageHeader + +#ifndef __vtkNIFTIImageReader_h +#define __vtkNIFTIImageReader_h + +#include "vtkIOImageModule.h" // For export macro +#include "vtkImageReader2.h" + +class vtkNIFTIImageHeader; +class vtkMatrix4x4; + +struct nifti_1_header; + +//---------------------------------------------------------------------------- +class VTKIOIMAGE_EXPORT vtkNIFTIImageReader : public vtkImageReader2 +{ +public: + // Description: + // Static method for construction. + static vtkNIFTIImageReader *New(); + vtkTypeMacro(vtkNIFTIImageReader, vtkImageReader2); + + // Description: + // Print information about this object. + virtual void PrintSelf(ostream& os, vtkIndent indent); + + // Description: + // Valid extensions for this file type. + virtual const char* GetFileExtensions() { + return ".nii .nii.gz .img .img.gz .hdr .hdr.gz"; } + + // Description: + // Return a descriptive name that might be useful in a GUI. + virtual const char* GetDescriptiveName() { + return "NIfTI"; } + + // Description: + // Return true if this reader can read the given file. + int CanReadFile(const char* filename); + + // Description: + // Read the time dimension as scalar components (default: Off). + // If this is on, then each time point will be stored as a component in + // the image data. If the file has both a time dimension and a vector + // dimension, then the number of components will be the product of these + // two dimensions, i.e. the components will store a sequence of vectors. + vtkGetMacro(TimeAsVector, bool); + vtkSetMacro(TimeAsVector, bool); + vtkBooleanMacro(TimeAsVector, bool); + + // Description: + // Get the time dimension that was stored in the NIFTI header. + int GetTimeDimension() { return this->Dim[4]; } + double GetTimeSpacing() { return this->PixDim[4]; } + + // Description: + // Get the slope and intercept for rescaling the scalar values. + // These values allow calibration of the data to real values. + // Use the equation v = u*RescaleSlope + RescaleIntercept. + // This directly returns the values stored in the scl_slope and + // scl_inter fields in the NIFTI header. + double GetRescaleSlope() { return this->RescaleSlope; } + double GetRescaleIntercept() { return this->RescaleIntercept; } + + // Description: + // QFac gives the slice order in the NIFTI file versus the VTK image. + // If QFac is -1, then the VTK slice index J is related to the NIFTI + // slice index j by the equation J = (num_slices - j - 1). VTK requires + // the slices to be ordered so that the voxel indices (I,J,K) provide a + // right-handed coordinate system, whereas NIFTI does not. Instead, + // NIFTI stores a factor called "qfac" in the header to signal when the + // (i,j,k) indices form a left-handed coordinate system. QFac will only + // ever have values of +1 or -1. + double GetQFac() { return this->QFac; } + + // Description: + // Get a matrix that gives the "qform" orientation and offset for the data. + // If no qform matrix was stored in the file, the return value is NULL. + // This matrix will transform VTK data coordinates into the NIFTI oriented + // data coordinates, where +X points right, +Y points anterior (toward the + // front), and +Z points superior (toward the head). The qform matrix will + // always have a positive determinant. The offset that is stored in the + // matrix gives the position of the first pixel in the first slice of the + // VTK image data. Note that if QFac is -1, then the first slice in the + // VTK image data is the last slice in the NIFTI file, and the Z offset + // will automatically be adjusted to compensate for this. + vtkMatrix4x4 *GetQFormMatrix() { return this->QFormMatrix; } + + // Description: + // Get a matrix that gives the "sform" orientation and offset for the data. + // If no sform matrix was stored in the file, the return value is NULL. + // Like the qform matrix, this matrix will transform VTK data coordinates + // into a NIFTI coordinate system. Unlike the qform matrix, the sform + // matrix can contain scaling information and can even (rarely) have + // a negative determinant, i.e. a flip. This matrix is modified slightly + // as compared to the sform matrix stored in the NIFTI header: the pixdim + // pixel spacing is factored out. Also, if QFac is -1, then the VTK slices + // are in reverse order as compared to the NIFTI slices, hence as compared + // to the sform matrix stored in the header, the third column of this matrix + // is multiplied by -1 and the Z offset is shifted to compensate for the + // fact that the last slice has become the first. + vtkMatrix4x4 *GetSFormMatrix() { return this->SFormMatrix; } + + // Description: + // Get the raw header information from the NIfTI file. + vtkNIFTIImageHeader *GetNIFTIHeader(); + +protected: + vtkNIFTIImageReader(); + ~vtkNIFTIImageReader(); + + // Description: + // Read the header information. + virtual int RequestInformation( + vtkInformation* request, vtkInformationVector** inputVector, + vtkInformationVector* outputVector); + + // Description: + // Read the voxel data. + virtual int RequestData( + vtkInformation* request, vtkInformationVector** inputVector, + vtkInformationVector* outputVector); + + // Description: + // Do a case-insensitive check for the given extension. + // The check will succeed if the filename ends in ".gz", and if the + // extension matches after removing the ".gz". + static bool CheckExtension(const char *fname, const char *ext); + + // Description: + // Make a new filename by replacing extension "ext1" with "ext2". + // The extensions must include a period, must be three characters + // long, and must be lower case. This method also verifies that + // the file exists, and adds or subtracts a ".gz" as necessary + // If the file exists, a new string is returned that must be + // deleted by the caller. Otherwise, the return value is NULL. + static char *ReplaceExtension( + const char *fname, const char *ext1, const char *ext2); + + // Description: + // Check the version of the header. + static int CheckNIFTIVersion(const nifti_1_header *hdr); + + // Description: + // Return true if an Analyze 7.5 header was found. + static bool CheckAnalyzeHeader(const nifti_1_header *hdr); + + // Description: + // Read the time dimension as if it was a vector dimension. + bool TimeAsVector; + + // Description: + // Information for rescaling data to quantitative units. + double RescaleIntercept; + double RescaleSlope; + + // Description: + // Is -1 if VTK slice order is opposite to NIFTI slice order, +1 otherwise. + double QFac; + + // Description: + // The orientation matrices for the NIFTI file. + vtkMatrix4x4 *QFormMatrix; + vtkMatrix4x4 *SFormMatrix; + + // Description: + // The dimensions of the NIFTI file. + int Dim[8]; + + // Description: + // The spacings in the NIFTI file. + double PixDim[8]; + + // Description: + // A copy of the header from the file that was most recently read. + vtkNIFTIImageHeader *NIFTIHeader; + +private: + vtkNIFTIImageReader(const vtkNIFTIImageReader&); // Not implemented. + void operator=(const vtkNIFTIImageReader&); // Not implemented. +}; + +#endif // __vtkNIFTIImageReader_h diff --git a/IO/Image/vtkNIFTIImageWriter.cxx b/IO/Image/vtkNIFTIImageWriter.cxx new file mode 100644 index 0000000000000000000000000000000000000000..93d9e92c5e04074a4b480e6d3ebe7c407e1be4a3 --- /dev/null +++ b/IO/Image/vtkNIFTIImageWriter.cxx @@ -0,0 +1,930 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: vtkNIFTIImageWriter.cxx + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ + +#include "vtkNIFTIImageWriter.h" +#include "vtkObjectFactory.h" +#include "vtkNIFTIImageReader.h" +#include "vtkImageData.h" +#include "vtkPointData.h" +#include "vtkInformation.h" +#include "vtkInformationVector.h" +#include "vtkStreamingDemandDrivenPipeline.h" +#include "vtkErrorCode.h" +#include "vtkByteSwap.h" +#include "vtkMatrix4x4.h" +#include "vtkMath.h" +#include "vtkCommand.h" +#include "vtkVersion.h" + +#include "vtksys/SystemTools.hxx" +#include "vtksys/ios/sstream" + +// Header for NIFTI +#include "vtkNIFTIImageHeader.h" +#include "vtkNIFTIImagePrivate.h" + +// Header for zlib +#include "vtk_zlib.h" + +#include <stdio.h> +#include <string.h> +#include <float.h> +#include <math.h> + +vtkStandardNewMacro(vtkNIFTIImageWriter); +vtkCxxSetObjectMacro(vtkNIFTIImageWriter,QFormMatrix,vtkMatrix4x4); +vtkCxxSetObjectMacro(vtkNIFTIImageWriter,SFormMatrix,vtkMatrix4x4); +vtkCxxSetObjectMacro(vtkNIFTIImageWriter,NIFTIHeader,vtkNIFTIImageHeader); + +//---------------------------------------------------------------------------- +vtkNIFTIImageWriter::vtkNIFTIImageWriter() +{ + this->FileLowerLeft = 1; + this->FileDimensionality = 3; + this->TimeDimension = 0; + this->TimeSpacing = 1.0; + // If slope,inter are 0,0 then default slope,inter of 1,0 is used + this->RescaleSlope = 0.0; + this->RescaleIntercept = 0.0; + this->QFac = 0.0; + this->QFormMatrix = 0; + this->SFormMatrix = 0; + this->OwnHeader = 0; + this->NIFTIHeader = 0; + this->NIFTIVersion = 0; + // Default description is "VTKX.Y.Z" + const char *version = vtkVersion::GetVTKVersion(); + size_t l = strlen(version); + this->Description = new char[l + 4]; + strncpy(this->Description, "VTK", 3); + strncpy(&this->Description[3], version, l); + this->Description[l + 3] = '\0'; +} + +//---------------------------------------------------------------------------- +vtkNIFTIImageWriter::~vtkNIFTIImageWriter() +{ + if (this->QFormMatrix) + { + this->QFormMatrix->Delete(); + } + if (this->SFormMatrix) + { + this->SFormMatrix->Delete(); + } + if (this->OwnHeader) + { + this->OwnHeader->Delete(); + } + if (this->NIFTIHeader) + { + this->NIFTIHeader->Delete(); + } + delete [] this->Description; +} + +//---------------------------------------------------------------------------- +vtkNIFTIImageHeader *vtkNIFTIImageWriter::GetNIFTIHeader() +{ + if (!this->NIFTIHeader) + { + this->NIFTIHeader = vtkNIFTIImageHeader::New(); + } + return this->NIFTIHeader; +} + +//---------------------------------------------------------------------------- +void vtkNIFTIImageWriter::PrintSelf(ostream& os, vtkIndent indent) +{ + this->Superclass::PrintSelf(os, indent); + + os << indent << "Description: " << this->Description << "\n"; + os << indent << "TimeDimension: " << this->TimeDimension << "\n"; + os << indent << "TimeSpacing: " << this->TimeSpacing << "\n"; + os << indent << "RescaleSlope: " << this->RescaleSlope << "\n"; + os << indent << "RescaleIntercept: " << this->RescaleIntercept << "\n"; + os << indent << "QFac: " << this->QFac << "\n"; + + os << indent << "QFormMatrix:"; + if (this->QFormMatrix) + { + double mat[16]; + vtkMatrix4x4::DeepCopy(mat, this->QFormMatrix); + for (int i = 0; i < 16; i++) + { + os << " " << mat[i]; + } + os << "\n"; + } + else + { + os << " (none)\n"; + } + + os << indent << "SFormMatrix:"; + if (this->SFormMatrix) + { + double mat[16]; + vtkMatrix4x4::DeepCopy(mat, this->SFormMatrix); + for (int i = 0; i < 16; i++) + { + os << " " << mat[i]; + } + os << "\n"; + } + else + { + os << " (none)\n"; + } + + os << indent << "NIFTIHeader: "; + if (this->NIFTIHeader) + { + os << this->NIFTIHeader << "\n"; + } + else + { + os << "(none)\n"; + } + os << indent << "NIFTIVersion: " << this->NIFTIVersion << "\n"; +} + +//---------------------------------------------------------------------------- +char *vtkNIFTIImageWriter::ReplaceExtension( + const char *filename, const char *ext1, const char *ext2) +{ + size_t n = strlen(filename); + size_t m = n; + char *newname = new char[n+4]; + strcpy(newname, filename); + + if (n > 2 && filename[n-3] == '.' && + tolower(filename[n-2]) == 'g' && + tolower(filename[n-1]) == 'z') + { + m -= 3; + } + if (m > 3 && filename[m-4] == '.' && + tolower(filename[m-3]) == tolower(ext1[1]) && + tolower(filename[m-2]) == tolower(ext1[2]) && + tolower(filename[m-1]) == tolower(ext1[3])) + { + if (isupper(filename[m-3])) + { + newname[m-3] = toupper(ext2[1]); + newname[m-2] = toupper(ext2[2]); + newname[m-1] = toupper(ext2[3]); + } + else + { + newname[m-3] = tolower(ext2[1]); + newname[m-2] = tolower(ext2[2]); + newname[m-1] = tolower(ext2[3]); + } + } + + return newname; +} + +//---------------------------------------------------------------------------- +namespace { + +// Initialize the NIFTI header with only the most basic information: +// - NIFTI data type is set from VTK data type +// - NIFTI pixdim set from VTK spacing +// - dimensionality is: +// - 5 if number of components is greater than one +// - 2 if Z dimension is one and number of components is one +// - 3 if Z dimension is greater than one and number of components is one +// - units are NIFTI_UNITS_UNKNOWN +// - intent is NIFTI_INTENT_NONE +// - magic is "n+1" (i.e. a .nii file, header+image in one file) +// - vox_offset is set to the header size plus 64-bit alignment padding +// - everything else is initialized to zero +// After initialization, the following should be set: +// - if file is ".hdr", set magic to "ni1" and vox_offset to zero +// - intent should be set, if known +// - units should be set, if known +// - qform and sform should be set, if known +// - pixdim[0] should be set to qfac (1 or -1) if qform is known +// - slope and inter should be set, if known +// - descrip and intent_name should be set, if known +// - for RGB and RGBA images, header should be modified as necessary +// - for complex images, header should be modified as necessary + +void vtkNIFTIImageWriterSetInformation( + nifti_2_header *hdr, + vtkInformation *info) +{ + // get the scalar information + vtkInformation *scalarInfo = vtkDataObject::GetActiveFieldInformation( + info, vtkDataObject::FIELD_ASSOCIATION_POINTS, + vtkDataSetAttributes::SCALARS); + + int extent[6]; + info->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT(), extent); + + double spacing[3]; + info->Get(vtkDataObject::SPACING(), spacing); + + int scalarType = scalarInfo->Get(vtkDataObject::FIELD_ARRAY_TYPE()); + int numComponents = scalarInfo->Get( + vtkDataObject::FIELD_NUMBER_OF_COMPONENTS()); + + // map VTK type to NIFTI type and bits + static const int typeMap[][3] = { +#if VTK_TYPE_CHAR_IS_SIGNED + { VTK_CHAR, NIFTI_TYPE_INT8, 8 }, +#else + { VTK_CHAR, NIFTI_TYPE_UINT8, 8 }, +#endif + { VTK_SIGNED_CHAR, NIFTI_TYPE_INT8, 8 }, + { VTK_UNSIGNED_CHAR, NIFTI_TYPE_UINT8, 8 }, + { VTK_SHORT, NIFTI_TYPE_INT16, 16 }, + { VTK_UNSIGNED_SHORT, NIFTI_TYPE_UINT16, 16 }, + { VTK_INT, NIFTI_TYPE_INT32, 32 }, + { VTK_UNSIGNED_INT, NIFTI_TYPE_UINT32, 32 }, +#if VTK_SIZEOF_LONG == 4 + { VTK_LONG, NIFTI_TYPE_INT32, 32 }, + { VTK_UNSIGNED_LONG, NIFTI_TYPE_UINT32, 32 }, +#else + { VTK_LONG, NIFTI_TYPE_INT64, 64 }, + { VTK_UNSIGNED_LONG, NIFTI_TYPE_UINT64, 64 }, +#endif + { VTK_LONG_LONG, NIFTI_TYPE_INT64, 64 }, + { VTK_UNSIGNED_LONG_LONG, NIFTI_TYPE_UINT64, 64 }, + { VTK___INT64, NIFTI_TYPE_INT64, 64 }, + { VTK_UNSIGNED___INT64, NIFTI_TYPE_UINT64, 64 }, + { VTK_FLOAT, NIFTI_TYPE_FLOAT32, 32 }, + { VTK_DOUBLE, NIFTI_TYPE_FLOAT64, 64 }, + { 0, 0, 0 } + }; + + short datatype = 0; + short databits = 0; + + // the end of the typemap has been reached when typeMap[2] is 0 + for (int i = 0; typeMap[2] != 0; i++) + { + if (scalarType == typeMap[i][0]) + { + datatype = typeMap[i][1]; + databits = typeMap[i][2]; + break; + } + } + + // number of spatial dimensions + int spaceDim = (extent[4] == extent[5] ? 2 : 3); + + hdr->dim[0] = (numComponents == 1 ? spaceDim : 5); + hdr->dim[1] = extent[1] - extent[0] + 1; + hdr->dim[2] = extent[3] - extent[2] + 1; + hdr->dim[3] = extent[5] - extent[4] + 1; + hdr->dim[4] = 1; + hdr->dim[5] = numComponents; + hdr->dim[6] = 1; + hdr->dim[7] = 1; + + hdr->datatype = datatype; + hdr->bitpix = databits; + + hdr->slice_start = 0; + hdr->pixdim[0] = 0.0; + hdr->pixdim[1] = spacing[0]; + hdr->pixdim[2] = spacing[1]; + hdr->pixdim[3] = spacing[2]; + hdr->pixdim[4] = 1.0; + hdr->pixdim[5] = 1.0; + hdr->pixdim[6] = 1.0; + hdr->pixdim[7] = 1.0; +} + +// Set the QForm from a 4x4 matrix +void vtkNIFTIImageWriterSetQForm( + nifti_2_header *hdr, double mmat[16], double qfac) +{ + double rmat[3][3]; + rmat[0][0] = mmat[0]; + rmat[0][1] = mmat[1]; + rmat[0][2] = mmat[2]; + rmat[1][0] = mmat[4]; + rmat[1][1] = mmat[5]; + rmat[1][2] = mmat[6]; + rmat[2][0] = mmat[8]; + rmat[2][1] = mmat[9]; + rmat[2][2] = mmat[10]; + + double quat[4]; + vtkMath::Matrix3x3ToQuaternion(rmat, quat); + if (quat[0] < 0) + { + quat[0] = -quat[0]; + quat[1] = -quat[1]; + quat[2] = -quat[2]; + quat[3] = -quat[3]; + } + + if (qfac < 0) + { + // We will be reversing the order of the slices, so the first VTK + // slice will be at the position of the last NIfTI slice, and we + // must adjust the offset to compensate for this. + mmat[3] -= rmat[0][2]*hdr->pixdim[3]*(hdr->dim[3] - 1); + mmat[7] -= rmat[1][2]*hdr->pixdim[3]*(hdr->dim[3] - 1); + mmat[11] -= rmat[2][2]*hdr->pixdim[3]*(hdr->dim[3] - 1); + } + + hdr->pixdim[0] = qfac; + hdr->quatern_b = quat[1]; + hdr->quatern_c = quat[2]; + hdr->quatern_d = quat[3]; + hdr->qoffset_x = mmat[3]; + hdr->qoffset_y = mmat[7]; + hdr->qoffset_z = mmat[11]; +} + +// Set the SForm from a 4x4 matrix +void vtkNIFTIImageWriterSetSForm( + nifti_2_header *hdr, double mmat[16], double qfac) +{ + if (qfac < 0) + { + // If QFac is set to -1 (which only occurs if qform_code was set) + // then the slices will be reversed, and we must reverse the slice + // orientation vector (the third column of the matrix) to compensate. + + // adjust the offset to compensate for changed slice ordering + mmat[3] -= mmat[2]*(hdr->dim[3] - 1); + mmat[7] -= mmat[6]*(hdr->dim[3] - 1); + mmat[11] -= mmat[10]*(hdr->dim[3] - 1); + + // reverse the slice orientation vector + mmat[2] = -mmat[2]; + mmat[6] = -mmat[6]; + mmat[10] = -mmat[10]; + } + + // first row + hdr->srow_x[0] = mmat[0] * hdr->pixdim[1]; + hdr->srow_x[1] = mmat[1] * hdr->pixdim[2]; + hdr->srow_x[2] = mmat[2] * hdr->pixdim[3]; + hdr->srow_x[3] = mmat[3]; + + // second row + hdr->srow_y[0] = mmat[4] * hdr->pixdim[1]; + hdr->srow_y[1] = mmat[5] * hdr->pixdim[2]; + hdr->srow_y[2] = mmat[6] * hdr->pixdim[3]; + hdr->srow_y[3] = mmat[7]; + + // third row + hdr->srow_z[0] = mmat[8] * hdr->pixdim[1]; + hdr->srow_z[1] = mmat[9] * hdr->pixdim[2]; + hdr->srow_z[2] = mmat[10] * hdr->pixdim[3]; + hdr->srow_z[3] = mmat[11]; +} + +void vtkNIFTIImageWriterMatrix( + double mmat[16], vtkMatrix4x4 *matrix, const double origin[3]) +{ + // find new offset by multiplying the origin by the matrix + double offset[4]; + offset[0] = origin[0]; + offset[1] = origin[1]; + offset[2] = origin[2]; + offset[3] = 1.0; + + if (matrix) + { + matrix->MultiplyPoint(offset, offset); + vtkMatrix4x4::DeepCopy(mmat, matrix); + } + else + { + vtkMatrix4x4::Identity(mmat); + } + + mmat[3] = offset[0]; + mmat[7] = offset[1]; + mmat[11] = offset[2]; +} + +} // end anonymous namespace + +//---------------------------------------------------------------------------- +int vtkNIFTIImageWriter::GenerateHeader(vtkInformation *info, bool singleFile) +{ + // create the header + nifti_2_header hdr; + int version = 0; + if (this->OwnHeader == 0) + { + this->OwnHeader = vtkNIFTIImageHeader::New(); + } + else + { + this->OwnHeader->Initialize(); + } + if (this->NIFTIHeader) + { + // use the header supplied by SetNIFTIHeader() + this->NIFTIHeader->GetHeader(&hdr); + version = hdr.magic[2] - '0'; + if (version > 2) + { + version = 2; + } + } + else + { + // start with a blank header + this->OwnHeader->GetHeader(&hdr); + hdr.scl_slope = 1.0; + } + + // copy the image information into the header + vtkNIFTIImageWriterSetInformation(&hdr, info); + if (hdr.datatype == 0) + { + vtkErrorMacro("Illegal data type for NIFTI file."); + return 0; + } + + // override the version if set via SetNIFTIVersion + if (this->NIFTIVersion != 0) + { + version = this->NIFTIVersion; + } + + // set the rescale slope/intercept if not (0.0,0.0) + if (this->RescaleSlope != 0.0 || this->RescaleIntercept != 0.0) + { + hdr.scl_slope = this->RescaleSlope; + hdr.scl_inter = this->RescaleIntercept; + } + + // set the header size + hdr.sizeof_hdr = (version == 2 ? + vtkNIFTIImageHeader::NIFTI2HeaderSize : + vtkNIFTIImageHeader::NIFTI1HeaderSize); + + // modify magic number and voxel offset for .img files + if (!singleFile) + { + strncpy(hdr.magic, (version == 2 ? "ni2" : "ni1"), 4); + hdr.vox_offset = 0; + } + else + { + strncpy(hdr.magic, (version == 2 ? "n+2" : "n+1"), 4); + hdr.vox_offset = (version == 2 ? 544 : 352); + } + if (version == 2) + { + // version 2 has four bytes for newline transfer checks + strncpy(&hdr.magic[4], "\r\n\032\n", 4); + } + + // set the description + if (this->Description) + { + strncpy(hdr.descrip, this->Description, 80); + } + + // qfac dictates the slice ordering in the file + double qfac = (this->QFac < 0 ? -1.0 : 1.0); + + // origin must be incorporated into qform and sform + double origin[3]; + info->Get(vtkDataObject::ORIGIN(), origin); + + if (this->QFormMatrix || + (origin[0] != 0 || origin[1] != 0 || origin[2] != 0)) + { + hdr.qform_code = 1; // SCANNER_ANAT + double mat16[16]; + vtkNIFTIImageWriterMatrix(mat16, this->QFormMatrix, origin); + vtkNIFTIImageWriterSetQForm(&hdr, mat16, qfac); + } + + if (this->SFormMatrix) + { + hdr.sform_code = 2; // ALIGNED_ANAT + double mat16[16]; + vtkNIFTIImageWriterMatrix(mat16, this->SFormMatrix, origin); + vtkNIFTIImageWriterSetSForm(&hdr, mat16, qfac); + } + + // base dimension not counting vector dimension + int basedim = (hdr.dim[3] == 1 ? 2 : 3); + + if (this->TimeDimension) + { + int tdim = this->TimeDimension; + if (hdr.dim[5] % tdim != 0) + { + vtkErrorMacro("Number of components in the image data must be " + "divisible by the TimeDimension"); + return 0; + } + hdr.pixdim[4] = this->TimeSpacing; + hdr.dim[4] = tdim; + hdr.dim[5] /= tdim; + hdr.dim[0] = (hdr.dim[5] > 1 ? 5 : 4); + basedim = 4; + } + + if (hdr.dim[5] == 2 && hdr.datatype == NIFTI_TYPE_FLOAT32) + { + // float with 2 components becomes COMPLEX64 + hdr.datatype = NIFTI_TYPE_COMPLEX64; + hdr.bitpix = 64; + hdr.dim[0] = basedim; + hdr.dim[5] = 1; + } + else if (hdr.dim[5] == 2 && hdr.datatype == NIFTI_TYPE_FLOAT64) + { + // double with 2 components becomes COMPLEX128 + hdr.datatype = NIFTI_TYPE_COMPLEX128; + hdr.bitpix = 32; + hdr.dim[0] = basedim; + hdr.dim[5] = 1; + } + else if (hdr.dim[5] == 3 && hdr.datatype == NIFTI_TYPE_UINT8) + { + // unsigned char with 3 components becomes RGB24 + hdr.datatype = NIFTI_TYPE_RGB24; + hdr.bitpix = 24; + hdr.dim[0] = basedim; + hdr.dim[5] = 1; + } + else if (hdr.dim[5] == 4 && hdr.datatype == NIFTI_TYPE_UINT8) + { + // unsigned char with 4 components becomes RGBA32 + hdr.datatype = NIFTI_TYPE_RGBA32; + hdr.bitpix = 32; + hdr.dim[0] = basedim; + hdr.dim[5] = 1; + } + + this->OwnHeader->SetHeader(&hdr); + return 1; +} + +//---------------------------------------------------------------------------- +int vtkNIFTIImageWriter::RequestData( + vtkInformation* vtkNotUsed(request), + vtkInformationVector** inputVector, + vtkInformationVector* vtkNotUsed(outputVector)) +{ + this->SetErrorCode(vtkErrorCode::NoError); + + vtkInformation *info = inputVector[0]->GetInformationObject(0); + vtkImageData *data = + vtkImageData::SafeDownCast(info->Get(vtkDataObject::DATA_OBJECT())); + + if (data == NULL) + { + vtkErrorMacro("No input provided!"); + return 0; + } + + const char *filename = this->GetFileName(); + if (filename == NULL) + { + vtkErrorMacro("A FileName must be provided"); + this->SetErrorCode(vtkErrorCode::NoFileNameError); + return 0; + } + + int extent[6]; + info->Get(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT(), extent); + + // use compression if name ends in .gz + bool isCompressed = false; + size_t n = strlen(filename); + size_t m = n; + if (n > 2 && filename[n-3] == '.' && + tolower(filename[n-2]) == 'g' && + tolower(filename[n-1]) == 'z') + { + m = n - 3; + isCompressed = true; + } + + // after the optional ".gz" is removed, is it a ".img/.hdr" file? + bool singleFile = true; + if (m > 4 && filename[m-4] == '.' && + ((tolower(filename[m-3]) == 'h' && + tolower(filename[m-2]) == 'd' && + tolower(filename[m-1]) == 'r') || + (tolower(filename[m-3]) == 'i' && + tolower(filename[m-2]) == 'm' && + tolower(filename[m-1]) == 'g'))) + { + singleFile = false; + } + + // generate the header information + if (this->GenerateHeader(info, singleFile) == 0) + { + return 0; + } + + // if file is not .nii, then get .hdr and .img filenames + char *hdrname = vtkNIFTIImageWriter::ReplaceExtension( + filename, ".img", ".hdr"); + char *imgname = vtkNIFTIImageWriter::ReplaceExtension( + filename, ".hdr", ".img"); + + vtkDebugMacro(<< "Writing NIFTI file " << hdrname); + + // get either a NIFTIv1 or a NIFTIv2 header + nifti_1_header hdr1; + nifti_2_header hdr2; + void *hdrptr = 0; + size_t hdrsize = 0; + int version = this->OwnHeader->GetMagic()[2] - '0'; + if (version == 2) + { + this->OwnHeader->GetHeader(&hdr2); + hdrptr = &hdr2; + hdrsize = hdr2.sizeof_hdr; + } + else + { + this->OwnHeader->GetHeader(&hdr1); + hdrptr = &hdr1; + hdrsize = hdr1.sizeof_hdr; + if (extent[1] - extent[0] + 1 > VTK_SHORT_MAX || + extent[3] - extent[2] + 1 > VTK_SHORT_MAX || + extent[5] - extent[4] + 1 > VTK_SHORT_MAX) + { + vtkErrorMacro("Image too large to store in NIFTI-1 format"); + delete [] hdrname; + delete [] imgname; + return 0; + } + } + + // try opening file + gzFile file = 0; + FILE *ufile = 0; + if (isCompressed) + { + file = gzopen(hdrname, "wb"); + } + else + { + ufile = fopen(hdrname, "wb"); + } + + if (!file && !ufile) + { + vtkErrorMacro("Cannot open file " << hdrname); + delete [] hdrname; + delete [] imgname; + this->SetErrorCode(vtkErrorCode::CannotOpenFileError); + return 0; + } + + this->InvokeEvent(vtkCommand::StartEvent); + this->UpdateProgress(0.0); + + // write the header + size_t bytesWritten = 0; + if (isCompressed) + { + unsigned int hsize = static_cast<unsigned int>(hdrsize); + int code = gzwrite(file, hdrptr, hsize); + bytesWritten = (code < 0 ? 0 : code); + } + else + { + bytesWritten = fwrite(hdrptr, 1, hdrsize, ufile); + } + if (bytesWritten < hdrsize) + { + this->SetErrorCode(vtkErrorCode::OutOfDiskSpaceError); + } + + if (singleFile && !this->ErrorCode) + { + // write the padding between the header and the image to the .nii file + size_t padsize = (static_cast<size_t>(this->OwnHeader->GetVoxOffset()) - + hdrsize); + char *padding = new char[padsize]; + memset(padding, '\0', padsize); + if (isCompressed) + { + int code = gzwrite(file, padding, static_cast<unsigned int>(padsize)); + bytesWritten = (code < 0 ? 0 : code); + } + else + { + bytesWritten = fwrite(padding, 1, padsize, ufile); + } + delete [] padding; + if (bytesWritten < padsize) + { + this->SetErrorCode(vtkErrorCode::OutOfDiskSpaceError); + } + } + else if (!this->ErrorCode) + { + // close the .hdr file and open the .img file + if (isCompressed) + { + gzclose(file); + file = gzopen(imgname, "wb"); + } + else + { + fclose(ufile); + ufile = fopen(imgname, "wb"); + } + } + + if (!file && !ufile) + { + vtkErrorMacro("Cannot open file " << imgname); + this->SetErrorCode(vtkErrorCode::CannotOpenFileError); + } + + // write the image + unsigned char *dataPtr = + static_cast<unsigned char *>(data->GetScalarPointer()); + + int swapBytes = 0; + int scalarSize = data->GetScalarSize(); + int numComponents = data->GetNumberOfScalarComponents(); + int outSizeX = static_cast<int>(this->OwnHeader->GetDim(1)); + int outSizeY = static_cast<int>(this->OwnHeader->GetDim(2)); + int outSizeZ = static_cast<int>(this->OwnHeader->GetDim(3)); + int timeDim = static_cast<int>(this->OwnHeader->GetDim(4)); + int vectorDim = static_cast<int>(this->OwnHeader->GetDim(5)); + + // for counting, include timeDim in vectorDim + vectorDim *= timeDim; + + z_off_t fileVoxelIncr = scalarSize*numComponents/vectorDim; + + // add a buffer for planar-vector to packed-vector conversion + unsigned char *rowBuffer = 0; + if (vectorDim > 1 || swapBytes) + { + rowBuffer = new unsigned char[outSizeX*fileVoxelIncr]; + } + + // special increment to reverse the slices if needed + vtkIdType sliceOffset = 0; + + if (this->QFac < 0) + { + // put slices in reverse order + sliceOffset = scalarSize*numComponents; + sliceOffset *= outSizeX; + sliceOffset *= outSizeY; + dataPtr += sliceOffset*(outSizeZ - 1); + } + + // report progress every 2% of the way to completion + vtkIdType target = + static_cast<vtkIdType>(0.02*outSizeY*outSizeZ*vectorDim) + 1; + vtkIdType count = 0; + + // write the data one row at a time, do planar-to-packed conversion + // of vector components if NIFTI file has a vector dimension + int rowSize = numComponents/vectorDim*outSizeX; + int c = 0; // counter for vector components + int j = 0; // counter for rows + int k = 0; // counter for slices + int t = 0; // counter for time + + unsigned char *ptr = dataPtr; + + while (!this->AbortExecute && !this->ErrorCode) + { + if (vectorDim == 1 && swapBytes == 0) + { + // write directly from input, instead of using a buffer + rowBuffer = ptr; + ptr += outSizeX*numComponents*scalarSize; + } + else + { + // create a vector plane from packed vector components + unsigned char *tmpPtr = rowBuffer; + z_off_t skipOther = scalarSize*numComponents - fileVoxelIncr; + for (int i = 0; i < outSizeX; i++) + { + // write one vector component of one voxel + z_off_t nn = fileVoxelIncr; + do { *tmpPtr++ = *ptr++; } while (--nn); + // skip past the other components + ptr += skipOther; + } + } + + if (swapBytes != 0 && scalarSize > 1) + { + vtkByteSwap::SwapVoidRange(rowBuffer, rowSize, scalarSize); + } + + if (isCompressed) + { + int code = gzwrite(file, rowBuffer, rowSize*scalarSize); + bytesWritten = (code < 0 ? 0 : code); + } + else + { + bytesWritten = fwrite(rowBuffer, scalarSize, rowSize, ufile)*scalarSize; + } + if (bytesWritten < static_cast<size_t>(rowSize*scalarSize)) + { + this->SetErrorCode(vtkErrorCode::OutOfDiskSpaceError); + break; + } + + if (++count % target == 0) + { + this->UpdateProgress(0.02*count/target); + } + + if (++j == outSizeY) + { + j = 0; + ptr -= 2*sliceOffset; // for reverse slice order + if (++k == outSizeZ) + { + k = 0; + if (++t == timeDim) + { + t = 0; + } + if (++c == vectorDim) + { + break; + } + // back up the ptr to the beginning of the image, + // then increment to the next vector component + ptr = dataPtr + c*fileVoxelIncr; + + if (timeDim > 1) + { + // if timeDim is included in the vectorDim (and hence in the + // VTK scalar components) then we have to make sure that + // the vector components are packed before the time steps + ptr = dataPtr + (c + t*(vectorDim - 1))/timeDim*fileVoxelIncr; + } + } + } + } + + // only delete this if it was alloced (if it was not alloced, it + // would have been set directly to a row out the output image) + if (vectorDim > 1 || swapBytes) + { + delete [] rowBuffer; + } + + if (isCompressed) + { + gzclose(file); + } + else + { + fclose(ufile); + } + + if (this->ErrorCode == vtkErrorCode::OutOfDiskSpaceError) + { + // erase the file, rather than leave a corrupt file on disk + vtkErrorMacro("Out of disk space, removing incomplete file " << imgname); + vtksys::SystemTools::RemoveFile(imgname); + if (!singleFile) + { + vtksys::SystemTools::RemoveFile(hdrname); + } + } + + this->UpdateProgress(1.0); + this->InvokeEvent(vtkCommand::EndEvent); + + delete [] hdrname; + delete [] imgname; + + return 1; +} diff --git a/IO/Image/vtkNIFTIImageWriter.h b/IO/Image/vtkNIFTIImageWriter.h new file mode 100644 index 0000000000000000000000000000000000000000..9729d9f424abd8a616ceb51630e689709d79cf38 --- /dev/null +++ b/IO/Image/vtkNIFTIImageWriter.h @@ -0,0 +1,177 @@ +/*========================================================================= + + Program: Visualization Toolkit + Module: vtkNIFTIImageWriter.h + + Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen + All rights reserved. + See Copyright.txt or http://www.kitware.com/Copyright.htm for details. + + This software is distributed WITHOUT ANY WARRANTY; without even + the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR + PURPOSE. See the above copyright notice for more information. + +=========================================================================*/ +// .NAME vtkNIFTIImageWriter - Write NIfTI-1 and NIfTI-2 medical image files +// .SECTION Description +// This class writes NIFTI files, either in .nii format or as separate +// .img and .hdr files. If told to write a file that ends in ".gz", +// then the writer will automatically compress the file with zlib. +// Images of type unsigned char that have 3 or 4 scalar components +// will automatically be written as RGB or RGBA respectively. Images +// of type float or double that have 2 components will automatically be +// written as complex values. +// .SECTION Thanks +// This class was contributed to VTK by the Calgary Image Processing and +// Analysis Centre (CIPAC). +// .SECTION See Also +// vtkNIFTIImageReader + +#ifndef __vtkNIFTIImageWriter_h +#define __vtkNIFTIImageWriter_h + +#include "vtkIOImageModule.h" // For export macro +#include "vtkImageWriter.h" + +class vtkMatrix4x4; +class vtkNIFTIImageHeader; + +class VTKIOIMAGE_EXPORT vtkNIFTIImageWriter : public vtkImageWriter +{ +public: + // Description: + // Static method for construction. + static vtkNIFTIImageWriter *New(); + vtkTypeMacro(vtkNIFTIImageWriter, vtkImageWriter); + + // Description: + // Print information about this object. + virtual void PrintSelf(ostream& os, vtkIndent indent); + + // Description: + // Set the version number for the NIfTI file format to use. + // This can be 1, 2, or 0 (the default). If set to zero, then it + // will save as NIfTI version 1 unless SetNIFTIHeader() provided + // header information from a NIfTI version 2 file. + vtkSetMacro(NIFTIVersion, int); + vtkGetMacro(NIFTIVersion, int); + + // Description: + // Set a short description (max 80 chars) of how the file was produced. + // The default description is "VTKX.Y" where X.Y is the VTK version. + vtkSetStringMacro(Description); + vtkGetStringMacro(Description); + + // Description: + // Set the time dimension to use in the NIFTI file (or zero if none). + // The number of components of the input data must be divisible by the time + // dimension if the time dimension is not set to zero. The vector dimension + // will be set to the number of components divided by the time dimension. + vtkGetMacro(TimeDimension, int); + vtkSetMacro(TimeDimension, int); + vtkGetMacro(TimeSpacing, double); + vtkSetMacro(TimeSpacing, double); + + // Description: + // Set the slope and intercept for calibrating the scalar values. + // Other programs that read the NIFTI file can use the equation + // v = u*RescaleSlope + RescaleIntercept to rescale the data to + // real values. If both the slope and the intercept are zero, + // then the SclSlope and SclIntercept in the header info provided + // via SetNIFTIHeader() are used instead. + vtkSetMacro(RescaleSlope, double); + vtkGetMacro(RescaleSlope, double); + vtkSetMacro(RescaleIntercept, double); + vtkGetMacro(RescaleIntercept, double); + + // Description: + // The QFac sets the ordering of the slices in the NIFTI file. + // If QFac is -1, then the slice ordering in the file will be reversed + // as compared to VTK. Use with caution. + vtkSetMacro(QFac, double); + vtkGetMacro(QFac, double); + + // Description: + // Set the "qform" orientation and offset for the image data. + // The 3x3 portion of the matrix must be orthonormal and have a + // positive determinant, it will be used to compute the quaternion. + // The last column of the matrix will be used for the offset. + // In the NIFTI header, the qform_code will be set to 1. + void SetQFormMatrix(vtkMatrix4x4 *); + vtkMatrix4x4 *GetQFormMatrix() { return this->QFormMatrix; } + + // Description: + // Set a matrix for the "sform" transformation stored in the file. + // Unlike the qform matrix, the sform matrix can contain scaling + // information. Before being stored in the NIFTI header, the + // first three columns of the matrix will be multipled by the voxel + // spacing. In the NIFTI header, the sform_code will be set to 2. + void SetSFormMatrix(vtkMatrix4x4 *); + vtkMatrix4x4 *GetSFormMatrix() { return this->SFormMatrix; } + + // Description: + // Set the NIFTI header information to use when writing the file. + // The data dimensions and pixdim from the supplied header will be + // ignored. Likewise, the QForm and SForm information in the supplied + // header will be ignored if you have called SetQFormMatrix() or + // SetSFormMatrix() to provide the orientation information for the file. + void SetNIFTIHeader(vtkNIFTIImageHeader *hdr); + vtkNIFTIImageHeader *GetNIFTIHeader(); + +protected: + vtkNIFTIImageWriter(); + ~vtkNIFTIImageWriter(); + + // Description: + // Generate the header information for the file. + int GenerateHeader(vtkInformation *info, bool singleFile); + + // Description: + // The main execution method, which writes the file. + virtual int RequestData(vtkInformation *request, + vtkInformationVector** inputVector, + vtkInformationVector* outputVector); + + // Description: + // Make a new filename by replacing extension "ext1" with "ext2". + // The extensions must include a period, must be three characters + // long, and must be lower case. A new string is returned that must + // be deleted by the caller. + static char *ReplaceExtension( + const char *fname, const char *ext1, const char *ext2); + + // Description: + // The size and spacing of the Time dimension to use in the file. + int TimeDimension; + double TimeSpacing; + + // Description: + // Information for rescaling data to quantitative units. + double RescaleIntercept; + double RescaleSlope; + + // Description: + // Is -1 if VTK slice order is opposite to NIFTI slice order, +1 otherwise. + double QFac; + + // Description: + // The orientation matrices for the NIFTI file. + vtkMatrix4x4 *QFormMatrix; + vtkMatrix4x4 *SFormMatrix; + + // Description + // A description of how the file was produced. + char *Description; + + // Description: + // The header information. + vtkNIFTIImageHeader *NIFTIHeader; + vtkNIFTIImageHeader *OwnHeader; + int NIFTIVersion; + +private: + vtkNIFTIImageWriter(const vtkNIFTIImageWriter&); // Not implemented. + void operator=(const vtkNIFTIImageWriter&); // Not implemented. +}; + +#endif // __vtkNIFTIImageWriter_h diff --git a/Testing/Data/ANALYZE.HDR.md5 b/Testing/Data/ANALYZE.HDR.md5 new file mode 100644 index 0000000000000000000000000000000000000000..d8c34f7c19dc68f93e2581818e4c23fbf28ffdcd --- /dev/null +++ b/Testing/Data/ANALYZE.HDR.md5 @@ -0,0 +1 @@ +5b2b4124474f98f661ac1fe226795372 diff --git a/Testing/Data/ANALYZE.IMG.GZ.md5 b/Testing/Data/ANALYZE.IMG.GZ.md5 new file mode 100644 index 0000000000000000000000000000000000000000..727d8abb3122e84980b5d23c7720182fda20c7e7 --- /dev/null +++ b/Testing/Data/ANALYZE.IMG.GZ.md5 @@ -0,0 +1 @@ +c911fa47ecb18f42c1105efd159b5c45 diff --git a/Testing/Data/avg152T1_RL_nifti2.nii.gz.md5 b/Testing/Data/avg152T1_RL_nifti2.nii.gz.md5 new file mode 100644 index 0000000000000000000000000000000000000000..5274fc32da2cdbacc87572e3594a067942f12432 --- /dev/null +++ b/Testing/Data/avg152T1_RL_nifti2.nii.gz.md5 @@ -0,0 +1 @@ +f03354cbec08d2fdb406284250ced1dc diff --git a/Testing/Data/filtered_func_data.nii.gz.md5 b/Testing/Data/filtered_func_data.nii.gz.md5 new file mode 100644 index 0000000000000000000000000000000000000000..9921bd95c5506656d487ffcf79d398c32c14e2ae --- /dev/null +++ b/Testing/Data/filtered_func_data.nii.gz.md5 @@ -0,0 +1 @@ +8a41f1e82d2bdafa821b7eb3b5b02b8c diff --git a/Testing/Data/nifti_rgb.nii.gz.md5 b/Testing/Data/nifti_rgb.nii.gz.md5 new file mode 100644 index 0000000000000000000000000000000000000000..e97c19e175394880ac39ed29d347f4e3d1c91bb1 --- /dev/null +++ b/Testing/Data/nifti_rgb.nii.gz.md5 @@ -0,0 +1 @@ +0e22bda40ebfbdcd51052ddd05d29194