Cloned library of VTK-5.0.0 with extra build files for internal package management.
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
 
 
 
 
 
 

733 lines
22 KiB

/*=========================================================================
Program: Visualization Toolkit
Module: $RCSfile: vtkImageReader.cxx,v $
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 "vtkImageReader.h"
#include "vtkByteSwap.h"
#include "vtkImageData.h"
#include "vtkInformation.h"
#include "vtkInformationVector.h"
#include "vtkObjectFactory.h"
#include "vtkPointData.h"
#include "vtkStreamingDemandDrivenPipeline.h"
#include "vtkTransform.h"
vtkCxxRevisionMacro(vtkImageReader, "$Revision: 1.121 $");
vtkStandardNewMacro(vtkImageReader);
vtkCxxSetObjectMacro(vtkImageReader,Transform,vtkTransform);
#ifdef read
#undef read
#endif
#ifdef close
#undef close
#endif
//----------------------------------------------------------------------------
vtkImageReader::vtkImageReader()
{
int idx;
for (idx = 0; idx < 3; ++idx)
{
this->DataVOI[idx*2] = this->DataVOI[idx*2 + 1] = 0;
}
// Left over from short reader
this->DataMask = 0xffff;
this->Transform = NULL;
this->ScalarArrayName = NULL;
this->SetScalarArrayName("ImageFile");
}
//----------------------------------------------------------------------------
vtkImageReader::~vtkImageReader()
{
this->SetTransform(NULL);
this->SetScalarArrayName(NULL);
}
//----------------------------------------------------------------------------
void vtkImageReader::PrintSelf(ostream& os, vtkIndent indent)
{
int idx;
this->Superclass::PrintSelf(os,indent);
os << indent << "Data Mask: " << this->DataMask << "\n";
os << indent << "DataVOI: (" << this->DataVOI[0];
for (idx = 1; idx < 6; ++idx)
{
os << ", " << this->DataVOI[idx];
}
os << ")\n";
if ( this->Transform )
{
os << indent << "Transform: " << this->Transform << "\n";
}
else
{
os << indent << "Transform: (none)\n";
}
os << indent << "ScalarArrayName: "
<< (this->ScalarArrayName ? this->ScalarArrayName : "(none)") << endl;
}
// This method returns the largest data that can be generated.
int vtkImageReader::RequestInformation (
vtkInformation * vtkNotUsed( request ),
vtkInformationVector** vtkNotUsed( inputVector ),
vtkInformationVector * outputVector)
{
// call the old method to help with backwards compatiblity
this->ExecuteInformation();
// get the info objects
vtkInformation* outInfo = outputVector->GetInformationObject(0);
double spacing[3];
int extent[6];
double origin[3];
// set the extent, if the VOI has not been set then default to
// the DataExtent
if (this->DataVOI[0] || this->DataVOI[1] ||
this->DataVOI[2] || this->DataVOI[3] ||
this->DataVOI[4] || this->DataVOI[5])
{
this->ComputeTransformedExtent(this->DataVOI,extent);
}
else
{
this->ComputeTransformedExtent(this->DataExtent,extent);
}
outInfo->Set(vtkStreamingDemandDrivenPipeline::WHOLE_EXTENT(),extent, 6);
// set the spacing
this->ComputeTransformedSpacing(spacing);
outInfo->Set(vtkDataObject::SPACING(), this->DataSpacing, 3);
// set the origin.
this->ComputeTransformedOrigin(origin);
outInfo->Set(vtkDataObject::ORIGIN(), this->DataOrigin, 3);
vtkDataObject::SetPointDataActiveScalarInfo(outInfo, this->DataScalarType,
this->NumberOfScalarComponents);
return 1;
}
int vtkImageReader::OpenAndSeekFile(int dataExtent[6], int idx)
{
unsigned long streamStart;
if (!this->FileName && !this->FilePattern)
{
vtkErrorMacro(<<"Either a FileName or FilePattern must be specified.");
return 0;
}
this->ComputeInternalFileName(idx);
this->OpenFile();
if ( !this->File )
{
return 0;
}
// convert data extent into constants that can be used to seek.
streamStart =
(dataExtent[0] - this->DataExtent[0]) * this->DataIncrements[0];
if (this->FileLowerLeft)
{
streamStart = streamStart +
(dataExtent[2] - this->DataExtent[2]) * this->DataIncrements[1];
}
else
{
streamStart = streamStart +
(this->DataExtent[3] - this->DataExtent[2] - dataExtent[2]) *
this->DataIncrements[1];
}
// handle three and four dimensional files
if (this->GetFileDimensionality() >= 3)
{
streamStart = streamStart +
(dataExtent[4] - this->DataExtent[4]) * this->DataIncrements[2];
}
streamStart += this->GetHeaderSize(idx);
// error checking
this->File->seekg(static_cast<long>(streamStart), ios::beg);
if (this->File->fail())
{
vtkErrorMacro(<< "File operation failed: " << streamStart << ", ext: "
<< dataExtent[0] << ", " << dataExtent[1] << ", "
<< dataExtent[2] << ", " << dataExtent[3] << ", "
<< dataExtent[4] << ", " << dataExtent[5]);
vtkErrorMacro(<< "Header size: " << this->GetHeaderSize(idx) << ", file ext: "
<< this->DataExtent[0] << ", " << this->DataExtent[1] << ", "
<< this->DataExtent[2] << ", " << this->DataExtent[3] << ", "
<< this->DataExtent[4] << ", " << this->DataExtent[5]);
return 0;
}
return 1;
}
//----------------------------------------------------------------------------
// This function reads in one data of data.
// templated to handle different data types.
template <class IT, class OT>
void vtkImageReaderUpdate2(vtkImageReader *self, vtkImageData *data,
IT *inPtr, OT *outPtr)
{
vtkIdType inIncr[3], outIncr[3];
OT *outPtr0, *outPtr1, *outPtr2;
long streamSkip0, streamSkip1;
unsigned long streamRead;
int idx0, idx1, idx2, pixelRead;
unsigned char *buf;
int inExtent[6];
int dataExtent[6];
int comp, pixelSkip;
long filePos, correction = 0;
unsigned long count = 0;
unsigned short DataMask;
unsigned long target;
// Get the requested extents.
data->GetExtent(inExtent);
// Convert them into to the extent needed from the file.
self->ComputeInverseTransformedExtent(inExtent,dataExtent);
// get and transform the increments
data->GetIncrements(inIncr);
self->ComputeInverseTransformedIncrements(inIncr,outIncr);
DataMask = self->GetDataMask();
// compute outPtr2
outPtr2 = outPtr;
if (outIncr[0] < 0)
{
outPtr2 = outPtr2 - outIncr[0]*(dataExtent[1] - dataExtent[0]);
}
if (outIncr[1] < 0)
{
outPtr2 = outPtr2 - outIncr[1]*(dataExtent[3] - dataExtent[2]);
}
if (outIncr[2] < 0)
{
outPtr2 = outPtr2 - outIncr[2]*(dataExtent[5] - dataExtent[4]);
}
// length of a row, num pixels read at a time
pixelRead = dataExtent[1] - dataExtent[0] + 1;
streamRead = static_cast<unsigned long>(pixelRead *
self->GetDataIncrements()[0]);
streamSkip0 = (long)(self->GetDataIncrements()[1] - streamRead);
streamSkip1 = (long)(self->GetDataIncrements()[2] -
(dataExtent[3] - dataExtent[2] + 1)* self->GetDataIncrements()[1]);
pixelSkip = data->GetNumberOfScalarComponents();
// read from the bottom up
if (!self->GetFileLowerLeft())
{
streamSkip0 = (long)(-static_cast<long>(streamRead)
- self->GetDataIncrements()[1]);
streamSkip1 = (long)(self->GetDataIncrements()[2] +
(dataExtent[3] - dataExtent[2] + 1)* self->GetDataIncrements()[1]);
}
// create a buffer to hold a row of the data
buf = new unsigned char[streamRead];
target = (unsigned long)((dataExtent[5]-dataExtent[4]+1)*
(dataExtent[3]-dataExtent[2]+1)/50.0);
target++;
// read the data row by row
if (self->GetFileDimensionality() == 3)
{
if ( !self->OpenAndSeekFile(dataExtent,0) )
{
delete [] buf;
return;
}
}
for (idx2 = dataExtent[4]; idx2 <= dataExtent[5]; ++idx2)
{
if (self->GetFileDimensionality() == 2)
{
if ( !self->OpenAndSeekFile(dataExtent,idx2) )
{
delete [] buf;
return;
}
}
outPtr1 = outPtr2;
for (idx1 = dataExtent[2];
!self->AbortExecute && idx1 <= dataExtent[3]; ++idx1)
{
if (!(count%target))
{
self->UpdateProgress(count/(50.0*target));
}
count++;
outPtr0 = outPtr1;
// read the row.
self->GetFile()->read((char *)buf, streamRead);
#ifdef __APPLE_CC__
if (static_cast<unsigned long>(self->GetFile()->gcount()) != streamRead)
// Apple's gcc3 returns fail when reading _to_ eof
#else
if ( static_cast<unsigned long>(self->GetFile()->gcount()) !=
streamRead || self->GetFile()->fail())
#endif
{
vtkGenericWarningMacro("File operation failed. row = " << idx1
<< ", Tried to Read = " << streamRead
<< ", Read = " << self->GetFile()->gcount()
<< ", Skip0 = " << streamSkip0
<< ", Skip1 = " << streamSkip1
<< ", FilePos = " << static_cast<vtkIdType>(self->GetFile()->tellg()));
delete [] buf;
return;
}
// handle swapping
if (self->GetSwapBytes())
{
// pixelSkip is the number of components in data
vtkByteSwap::SwapVoidRange(buf, pixelRead*pixelSkip, sizeof(IT));
}
// copy the bytes into the typed data
inPtr = (IT *)(buf);
for (idx0 = dataExtent[0]; idx0 <= dataExtent[1]; ++idx0)
{
// Copy pixel into the output.
if (DataMask == 0xffff)
{
for (comp = 0; comp < pixelSkip; comp++)
{
outPtr0[comp] = (OT)(inPtr[comp]);
}
}
else
{
// left over from short reader (what about other types.
for (comp = 0; comp < pixelSkip; comp++)
{
outPtr0[comp] = (OT)((short)(inPtr[comp]) & DataMask);
}
}
// move to next pixel
inPtr += pixelSkip;
outPtr0 += outIncr[0];
}
// move to the next row in the file and data
filePos = self->GetFile()->tellg();
/* Unfortunately this doesn't work as a fix, but I'll leave it here for a bit
to provoke ideas.
if (filePos == -1)
{
self->GetFile()->clear(self->GetFile()->rdstate() & ~ios::eofbit);
self->GetFile()->clear(self->GetFile()->rdstate() & ~ios::failbit);
filePos = self->GetFile()->tellg();
}
*/
#if defined(VTK_USE_ANSI_STDLIB ) && ((defined(__sgi) && !defined(__GNUC__)) \
|| (__BORLANDC__>=0x0560) || defined (__APPLE_CC__))
// this check is required for SGI's when vtk is build with VTK_USE_ANSI_STDLIB
// seems that after a read that just reaches EOF, tellg reports a -1.
// clear() does not work, so we have to reopen the file.
// NB: Also Borland CBuilder 6 suffers from same trouble
// As does Apple's gcc3
if (filePos == -1)
{
self->OpenFile();
self->GetFile()->seekg(0,ios::end);
filePos = self->GetFile()->tellg();
}
#endif
// watch for case where we might rewind too much
// if that happens, store the value in correction and apply later
if (filePos + streamSkip0 >= 0)
{
self->GetFile()->seekg(static_cast<long>(self->GetFile()->tellg()) + streamSkip0, ios::beg);
correction = 0;
}
else
{
correction = streamSkip0;
}
outPtr1 += outIncr[1];
}
// move to the next image in the file and data
self->GetFile()->seekg(static_cast<long>(self->GetFile()->tellg()) + streamSkip1 + correction,
ios::beg);
outPtr2 += outIncr[2];
}
// delete the temporary buffer
delete [] buf;
}
//----------------------------------------------------------------------------
// This function reads in one data of one slice.
// templated to handle different data types.
template <class T>
void vtkImageReaderUpdate1(vtkImageReader *self, vtkImageData *data, T *inPtr)
{
void *outPtr;
// Call the correct templated function for the input
outPtr = data->GetScalarPointer();
switch (data->GetScalarType())
{
vtkTemplateMacro(vtkImageReaderUpdate2(self, data, inPtr,
(VTK_TT *)(outPtr)));
default:
vtkGenericWarningMacro("Update1: Unknown data type\n");
}
}
//----------------------------------------------------------------------------
// This function reads a data from a file. The datas extent/axes
// are assumed to be the same as the file extent/order.
void vtkImageReader::ExecuteData(vtkDataObject *output)
{
vtkImageData *data = this->AllocateOutputData(output);
void *ptr = NULL;
int *ext;
if (!this->FileName && !this->FilePattern)
{
vtkErrorMacro("Either a valid FileName or FilePattern must be specified.");
return;
}
ext = data->GetExtent();
if (!data->GetPointData()->GetScalars())
{
return;
}
data->GetPointData()->GetScalars()->SetName(this->ScalarArrayName);
vtkDebugMacro("Reading extent: " << ext[0] << ", " << ext[1] << ", "
<< ext[2] << ", " << ext[3] << ", " << ext[4] << ", " << ext[5]);
this->ComputeDataIncrements();
// Call the correct templated function for the output
switch (this->GetDataScalarType())
{
vtkTemplateMacro(vtkImageReaderUpdate1(this, data, (VTK_TT *)(ptr)));
default:
vtkErrorMacro(<< "UpdateFromFile: Unknown data type");
}
}
void vtkImageReader::ComputeTransformedSpacing (double Spacing[3])
{
if (!this->Transform)
{
memcpy (Spacing, this->DataSpacing, 3 * sizeof (double));
}
else
{
double transformedSpacing[3];
memcpy (transformedSpacing, this->DataSpacing, 3 * sizeof (double));
this->Transform->TransformVector(transformedSpacing, transformedSpacing);
for (int i = 0; i < 3; i++)
{
Spacing[i] = fabs(transformedSpacing[i]);
}
vtkDebugMacro("Transformed Spacing " << Spacing[0] << ", " << Spacing[1] << ", " << Spacing[2]);
}
}
// if the spacing is negative we need to tranlate the origin
// basically O' = O + spacing*(dim-1) for any axis that would
// have a negative spaing
void vtkImageReader::ComputeTransformedOrigin (double origin[3])
{
if (!this->Transform)
{
memcpy (origin, this->DataOrigin, 3 * sizeof (double));
}
else
{
double transformedOrigin[3];
double transformedSpacing[3];
int transformedExtent[6];
memcpy (transformedSpacing, this->DataSpacing, 3 * sizeof (double));
this->Transform->TransformVector(transformedSpacing, transformedSpacing);
memcpy (transformedOrigin, this->DataOrigin, 3 * sizeof (double));
this->Transform->TransformPoint(transformedOrigin, transformedOrigin);
this->ComputeTransformedExtent(this->DataExtent,transformedExtent);
for (int i = 0; i < 3; i++)
{
if (transformedSpacing[i] < 0)
{
origin[i] = transformedOrigin[i] + transformedSpacing[i]*
(transformedExtent[i*2+1] - transformedExtent[i*2] + 1);
}
else
{
origin[i] = transformedOrigin[i];
}
}
vtkDebugMacro("Transformed Origin " << origin[0] << ", " << origin[1] << ", " << origin[2]);
}
}
void vtkImageReader::ComputeTransformedExtent(int inExtent[6],
int outExtent[6])
{
double transformedExtent[3];
int temp;
int idx;
int dataExtent[6];
if (!this->Transform)
{
memcpy (outExtent, inExtent, 6 * sizeof (int));
memcpy (dataExtent, this->DataExtent, 6 * sizeof(int));
}
else
{
// need to know how far to translate to start at 000
// first transform the data extent
transformedExtent[0] = this->DataExtent[0];
transformedExtent[1] = this->DataExtent[2];
transformedExtent[2] = this->DataExtent[4];
this->Transform->TransformPoint(transformedExtent, transformedExtent);
dataExtent[0] = (int) transformedExtent[0];
dataExtent[2] = (int) transformedExtent[1];
dataExtent[4] = (int) transformedExtent[2];
transformedExtent[0] = this->DataExtent[1];
transformedExtent[1] = this->DataExtent[3];
transformedExtent[2] = this->DataExtent[5];
this->Transform->TransformPoint(transformedExtent, transformedExtent);
dataExtent[1] = (int) transformedExtent[0];
dataExtent[3] = (int) transformedExtent[1];
dataExtent[5] = (int) transformedExtent[2];
for (idx = 0; idx < 6; idx += 2)
{
if (dataExtent[idx] > dataExtent[idx+1])
{
temp = dataExtent[idx];
dataExtent[idx] = dataExtent[idx+1];
dataExtent[idx+1] = temp;
}
}
// now transform the inExtent
transformedExtent[0] = inExtent[0];
transformedExtent[1] = inExtent[2];
transformedExtent[2] = inExtent[4];
this->Transform->TransformPoint(transformedExtent, transformedExtent);
outExtent[0] = (int) transformedExtent[0];
outExtent[2] = (int) transformedExtent[1];
outExtent[4] = (int) transformedExtent[2];
transformedExtent[0] = inExtent[1];
transformedExtent[1] = inExtent[3];
transformedExtent[2] = inExtent[5];
this->Transform->TransformPoint(transformedExtent, transformedExtent);
outExtent[1] = (int) transformedExtent[0];
outExtent[3] = (int) transformedExtent[1];
outExtent[5] = (int) transformedExtent[2];
}
for (idx = 0; idx < 6; idx += 2)
{
if (outExtent[idx] > outExtent[idx+1])
{
temp = outExtent[idx];
outExtent[idx] = outExtent[idx+1];
outExtent[idx+1] = temp;
}
// do the slide to 000 origin by subtracting the minimum extent
outExtent[idx] -= dataExtent[idx];
outExtent[idx+1] -= dataExtent[idx];
}
vtkDebugMacro(<< "Transformed extent are:"
<< outExtent[0] << ", " << outExtent[1] << ", "
<< outExtent[2] << ", " << outExtent[3] << ", "
<< outExtent[4] << ", " << outExtent[5]);
}
void vtkImageReader::ComputeInverseTransformedExtent(int inExtent[6],
int outExtent[6])
{
double transformedExtent[3];
int temp;
int idx;
if (!this->Transform)
{
memcpy (outExtent, inExtent, 6 * sizeof (int));
for (idx = 0; idx < 6; idx += 2)
{
// do the slide to 000 origin by subtracting the minimum extent
outExtent[idx] += this->DataExtent[idx];
outExtent[idx+1] += this->DataExtent[idx];
}
}
else
{
// need to know how far to translate to start at 000
int dataExtent[6];
// first transform the data extent
transformedExtent[0] = this->DataExtent[0];
transformedExtent[1] = this->DataExtent[2];
transformedExtent[2] = this->DataExtent[4];
this->Transform->TransformPoint(transformedExtent, transformedExtent);
dataExtent[0] = (int) transformedExtent[0];
dataExtent[2] = (int) transformedExtent[1];
dataExtent[4] = (int) transformedExtent[2];
transformedExtent[0] = this->DataExtent[1];
transformedExtent[1] = this->DataExtent[3];
transformedExtent[2] = this->DataExtent[5];
this->Transform->TransformPoint(transformedExtent, transformedExtent);
dataExtent[1] = (int) transformedExtent[0];
dataExtent[3] = (int) transformedExtent[1];
dataExtent[5] = (int) transformedExtent[2];
for (idx = 0; idx < 6; idx += 2)
{
if (dataExtent[idx] > dataExtent[idx+1])
{
temp = dataExtent[idx];
dataExtent[idx] = dataExtent[idx+1];
dataExtent[idx+1] = temp;
}
}
for (idx = 0; idx < 6; idx += 2)
{
// do the slide to 000 origin by subtracting the minimum extent
inExtent[idx] += dataExtent[idx];
inExtent[idx+1] += dataExtent[idx];
}
transformedExtent[0] = inExtent[0];
transformedExtent[1] = inExtent[2];
transformedExtent[2] = inExtent[4];
this->Transform->GetLinearInverse()->TransformPoint(transformedExtent,
transformedExtent);
outExtent[0] = (int) transformedExtent[0];
outExtent[2] = (int) transformedExtent[1];
outExtent[4] = (int) transformedExtent[2];
transformedExtent[0] = inExtent[1];
transformedExtent[1] = inExtent[3];
transformedExtent[2] = inExtent[5];
this->Transform->GetLinearInverse()->TransformPoint(transformedExtent,
transformedExtent);
outExtent[1] = (int) transformedExtent[0];
outExtent[3] = (int) transformedExtent[1];
outExtent[5] = (int) transformedExtent[2];
for (idx = 0; idx < 6; idx += 2)
{
if (outExtent[idx] > outExtent[idx+1])
{
temp = outExtent[idx];
outExtent[idx] = outExtent[idx+1];
outExtent[idx+1] = temp;
}
}
}
vtkDebugMacro(<< "Inverse Transformed extent are:"
<< outExtent[0] << ", " << outExtent[1] << ", "
<< outExtent[2] << ", " << outExtent[3] << ", "
<< outExtent[4] << ", " << outExtent[5]);
}
void vtkImageReader::ComputeTransformedIncrements(vtkIdType inIncr[3],
vtkIdType outIncr[3])
{
double transformedIncr[3];
if (!this->Transform)
{
memcpy (outIncr, inIncr, 3 * sizeof (vtkIdType));
}
else
{
transformedIncr[0] = inIncr[0];
transformedIncr[1] = inIncr[1];
transformedIncr[2] = inIncr[2];
this->Transform->TransformVector(transformedIncr, transformedIncr);
outIncr[0] = (vtkIdType) transformedIncr[0];
outIncr[1] = (vtkIdType) transformedIncr[1];
outIncr[2] = (vtkIdType) transformedIncr[2];
vtkDebugMacro(<< "Transformed Incr are:"
<< outIncr[0] << ", " << outIncr[1] << ", " << outIncr[2]);
}
}
void vtkImageReader::ComputeInverseTransformedIncrements(vtkIdType inIncr[3],
vtkIdType outIncr[3])
{
double transformedIncr[3];
if (!this->Transform)
{
memcpy (outIncr, inIncr, 3 * sizeof (vtkIdType));
}
else
{
transformedIncr[0] = inIncr[0];
transformedIncr[1] = inIncr[1];
transformedIncr[2] = inIncr[2];
this->Transform->GetLinearInverse()->TransformVector(transformedIncr,
transformedIncr);
outIncr[0] = (vtkIdType) transformedIncr[0];
outIncr[1] = (vtkIdType) transformedIncr[1];
outIncr[2] = (vtkIdType) transformedIncr[2];
vtkDebugMacro(<< "Inverse Transformed Incr are:"
<< outIncr[0] << ", " << outIncr[1] << ", " << outIncr[2]);
}
}