Commit e7aa61da authored by David Gobbi's avatar David Gobbi

Add reader/writer for NIFTIv1 and NIFTIv2 files.

This patch adds support for reading and writing NIFTI files, as well
as full access to the NIFTI header.  Versions 1 and 2 of the NIFTI
header are supported, and the reader can also be used to read Analyze 7.5
files, though the only Analyze header elements that are supported are the
ones that overlap with NIFTI.  Automatic compression and expansion for
gzipped files is also provided whenever the filename ends with a .gz
extension.

Change-Id: Iead3afb1a65d8c01aa52f3ca55c6da2f08703375
parent c13f410a
......@@ -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})
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
......
/*=========================================================================
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);
}
/*=========================================================================
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);
}
This diff is collapsed.
......@@ -6,6 +6,7 @@ vtk_add_test_python(
TestTIFFReader.py
dem.py
TestMetaImage2D.py
TestNIFTIReaderWriter.py
TestSetFileNames.py
TestImageJSONWriter.py,NO_VALID
)
#! /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)
......@@ -20,6 +20,7 @@ vtk_module(vtkIOImage
vtkTestingCore
vtkImagingSources
vtkImagingMath
vtkInteractionStyle
vtkRenderingContext2D
vtkRenderingCore
vtkTestingCore
......
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
This diff is collapsed.
/*=========================================================================
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