From ce8987fa603ce1591d6d595273761d092467575d Mon Sep 17 00:00:00 2001 From: Simon Rit Date: Tue, 13 Nov 2012 18:21:13 +0100 Subject: [PATCH] Copied ESRF Edf processing from rtk and use rtk file names to emphasize that they are comming from rtk. --- common/CMakeLists.txt | 2 + common/clitkIO.cxx | 2 + common/rtkEdfImageIO.cxx | 286 ++++++++++++++++++++++++++++++++ common/rtkEdfImageIO.h | 157 ++++++++++++++++++ common/rtkEdfImageIOFactory.cxx | 29 ++++ common/rtkEdfImageIOFactory.h | 77 +++++++++ 6 files changed, 553 insertions(+) create mode 100644 common/rtkEdfImageIO.cxx create mode 100644 common/rtkEdfImageIO.h create mode 100644 common/rtkEdfImageIOFactory.cxx create mode 100644 common/rtkEdfImageIOFactory.h diff --git a/common/CMakeLists.txt b/common/CMakeLists.txt index d23e0b3..8ad60a0 100644 --- a/common/CMakeLists.txt +++ b/common/CMakeLists.txt @@ -22,6 +22,8 @@ SET(clitkCommon_SRC rtkHisImageIOFactory.cxx rtkHndImageIO.cxx rtkHndImageIOFactory.cxx + rtkEdfImageIO.cxx + rtkEdfImageIOFactory.cxx clitkDicomRTDoseIO.cxx clitkDicomRTDoseIOFactory.cxx clitkOrientation.cxx diff --git a/common/clitkIO.cxx b/common/clitkIO.cxx index 46e0580..271d5a0 100644 --- a/common/clitkIO.cxx +++ b/common/clitkIO.cxx @@ -32,6 +32,7 @@ #include "clitkXdrImageIOFactory.h" #include "rtkHisImageIOFactory.h" #include "rtkHndImageIOFactory.h" +#include "rtkEdfImageIOFactory.h" #include "clitkGateAsciiImageIOFactory.h" #include "clitkConfiguration.h" #if CLITK_PRIVATE_FEATURES @@ -59,5 +60,6 @@ void clitk::RegisterClitkFactories() clitk::XdrImageIOFactory::RegisterOneFactory(); rtk::HisImageIOFactory::RegisterOneFactory(); rtk::HndImageIOFactory::RegisterOneFactory(); + rtk::EdfImageIOFactory::RegisterOneFactory(); } //// diff --git a/common/rtkEdfImageIO.cxx b/common/rtkEdfImageIO.cxx new file mode 100644 index 0000000..ed3cf10 --- /dev/null +++ b/common/rtkEdfImageIO.cxx @@ -0,0 +1,286 @@ +/*========================================================================= + * + * Copyright RTK Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ + +#include "rtkEdfImageIO.h" + +// itk include (for itkReadRawBytesAfterSwappingMacro) +#include +#include + +//-------------------------------------------------------------------- +/* Find value_ptr as pointer to the parameter of the given key in the header. + * Returns NULL on success. + */ +char* +rtk::EdfImageIO::edf_findInHeader( char* header, const char* key ) +{ + char *value_ptr = strstr( header, key ); + + if (!value_ptr) return NULL; + /* an edf line is "key = value ;" */ + value_ptr = 1 + strchr( value_ptr + strlen(key), '=' ); + while (isspace(*value_ptr) ) value_ptr++; + return value_ptr; +} + +//-------------------------------------------------------------------- +// Read Image Information +void rtk::EdfImageIO::ReadImageInformation() +{ + int k; + char * header = NULL; + int header_size = 0; + char * p; + gzFile inp; + + inp = gzopen(m_FileName.c_str(), "rb"); + if (!inp) + itkGenericExceptionMacro(<< "Cannot open input file " << m_FileName); + + // read header: it is a multiple of 512 B ending by "}\n" + while (header_size == 0 || strncmp(&header[header_size-2],"}\n",2) ) { + int header_size_prev = header_size; + header_size += 512; + if (!header) + header = (char*)malloc(header_size+1); + else + header = (char*)realloc(header, header_size+1); + header[header_size_prev] = 0; /* protection against empty file */ + // fread(header+header_size_prev, 512, 1, fp); + k = gzread(inp, header+header_size_prev, 512); + if (k < 512) { /* protection against infinite loop */ + gzclose(inp); + free(header); + itkGenericExceptionMacro(<< "Damaged EDF header of " + << m_FileName + << ": not multiple of 512 B."); + } + header[header_size] = 0; /* end of string: protection against strstr later + on */ + } + + // parse the header + int dim1 = -1, dim2 = -1, datalen = -1; + char *otherfile_name = 0; // this file, or another file with the data (EDF vs + // EHF formats) + int otherfile_skip = 0; + + if ( (p = edf_findInHeader(header, "EDF_BinaryFileName") ) ) { + int plen = strcspn(p, " ;\n"); + otherfile_name = (char*)realloc(otherfile_name, plen+1); + strncpy(otherfile_name, p, plen); + otherfile_name[plen] = '\0'; + if ( (p = edf_findInHeader(header, "EDF_BinaryFilePosition") ) ) + otherfile_skip = atoi(p); + } + + if ( (p = edf_findInHeader(header, "Dim_1") ) ) + dim1 = atoi(p); + if ( (p = edf_findInHeader(header, "Dim_2") ) ) + dim2 = atoi(p); + +// int orig1 = -1, orig2 = -1; +// if ((p = edf_findInHeader(header, "row_beg"))) +// orig1 = atoi(p); +// if ((p = edf_findInHeader(header, "col_beg"))) +// orig2 = atoi(p); + + static const struct table3 edf_datatype_table[] = + { + { "UnsignedByte", U_CHAR_DATATYPE, 1 }, + { "SignedByte", CHAR_DATATYPE, 1 }, + { "UnsignedShort", U_SHORT_DATATYPE, 2 }, + { "SignedShort", SHORT_DATATYPE, 2 }, + { "UnsignedInteger", U_INT_DATATYPE, 4 }, + { "SignedInteger", INT_DATATYPE, 4 }, + { "UnsignedLong", U_L_INT_DATATYPE, 4 }, + { "SignedLong", L_INT_DATATYPE, 4 }, + { "FloatValue", FLOAT_DATATYPE, 4 }, + { "DoubleValue", DOUBLE_DATATYPE, 8 }, + { "Float", FLOAT_DATATYPE, 4 }, // Float and FloatValue + // are synonyms + { "Double", DOUBLE_DATATYPE, 8 }, // Double and DoubleValue + // are synonyms + { NULL, -1, -1 } + }; + if ( (p = edf_findInHeader(header, "DataType") ) ) { + k = lookup_table3_nth(edf_datatype_table, p); + if (k < 0) { // unknown EDF DataType + gzclose(inp); + free(header); + itkGenericExceptionMacro( <<"Unknown EDF datatype \"" + << p + << "\""); + } + datalen = edf_datatype_table[k].sajzof; + switch(k) { + case U_CHAR_DATATYPE: + SetComponentType(itk::ImageIOBase::UCHAR); + break; + case CHAR_DATATYPE: + SetComponentType(itk::ImageIOBase::CHAR); + break; + case U_SHORT_DATATYPE: + SetComponentType(itk::ImageIOBase::USHORT); + break; + case SHORT_DATATYPE: + SetComponentType(itk::ImageIOBase::SHORT); + break; + case U_INT_DATATYPE: + SetComponentType(itk::ImageIOBase::UINT); + break; + case INT_DATATYPE: + SetComponentType(itk::ImageIOBase::INT); + break; + case U_L_INT_DATATYPE: + SetComponentType(itk::ImageIOBase::ULONG); + break; + case L_INT_DATATYPE: + SetComponentType(itk::ImageIOBase::LONG); + break; + case FLOAT_DATATYPE: + SetComponentType(itk::ImageIOBase::FLOAT); + break; + case DOUBLE_DATATYPE: + SetComponentType(itk::ImageIOBase::DOUBLE); + break; + } + } + + static const struct table edf_byteorder_table[] = + { + { "LowByteFirst", LittleEndian }, /* little endian */ + { "HighByteFirst", BigEndian }, /* big endian */ + { NULL, -1 } + }; + + int byteorder = LittleEndian; + if ( (p = edf_findInHeader(header, "ByteOrder") ) ) { + k = lookup_table_nth(edf_byteorder_table, p); + if (k >= 0) { + + byteorder = edf_byteorder_table[k].value; + if(byteorder==LittleEndian) + this->SetByteOrder(LittleEndian); + else + this->SetByteOrder(BigEndian); + } + } else + itkWarningMacro(<<"ByteOrder not specified in the header! Not swapping bytes (figure may not be correct)."); + // Get and verify size of the data: + int datasize = dim1 * dim2 * datalen; + if ( (p = edf_findInHeader(header, "Size") ) ) { + int d = atoi(p); + if (d != datasize) { + itkWarningMacro(<< "Size " << datasize << " is not " + << dim1 << 'x' << dim2 << "x" << datalen + << " = " << d << ". Supposing the latter."); + } + } + + // EHF files: binary data are in another file than the header file + m_BinaryFileName = m_FileName; + m_BinaryFileSkip = header_size; + if (otherfile_name) { + m_BinaryFileName = std::string(otherfile_name); + m_BinaryFileSkip = otherfile_skip; + } + + double spacing = 1.; + if ( (p = edf_findInHeader(header, "optic_used") ) ) + spacing = atof(p); + + free(header); + gzclose(inp); + + SetNumberOfDimensions(2); + SetDimensions(0, dim1); + SetDimensions(1, dim2); + SetSpacing(0, spacing); + SetSpacing(1, spacing); + SetOrigin(0, 0.); + SetOrigin(1, 0.); +} //// + +//-------------------------------------------------------------------- +// Read Image Information +bool rtk::EdfImageIO::CanReadFile(const char* FileNameToRead) +{ + std::string filename(FileNameToRead); + const std::string::size_type it = filename.find_last_of( "." ); + std::string fileExt( filename, it+1, filename.length() ); + + if (fileExt != std::string("edf") ) return false; + return true; +} //// + +//-------------------------------------------------------------------- +// Read Image Content +void rtk::EdfImageIO::Read(void * buffer) +{ + gzFile inp; + + inp = gzopen(m_BinaryFileName.c_str(), "rb"); + if (!inp) + itkGenericExceptionMacro(<< "Cannot open file \"" << m_FileName << "\""); + gzseek(inp, m_BinaryFileSkip, SEEK_SET); + + // read the data (image) + long numberOfBytesToBeRead = GetComponentSize(); + for(unsigned int i=0; i +#include +#include + +namespace rtk { + +/** \class EdfImageIO + * \brief Class for reading Edf image file format. Edf is the format of + * X-ray projection images at the ESRF. + * + * \author Simon Rit + * + * \ingroup IOFilters + */ +class EdfImageIO : public itk::ImageIOBase +{ +public: + /** Standard class typedefs. */ + typedef EdfImageIO Self; + typedef itk::ImageIOBase Superclass; + typedef itk::SmartPointer Pointer; + + EdfImageIO() : Superclass() { + } + + /** Method for creation through the object factory. */ + itkNewMacro(Self); + + /** Run-time type information (and related methods). */ + itkTypeMacro(EdfImageIO, ImageIOBase); + + /*-------- This part of the interface deals with reading data. ------ */ + virtual void ReadImageInformation(); + + virtual bool CanReadFile( const char* FileNameToRead ); + + virtual void Read(void * buffer); + + /*-------- This part of the interfaces deals with writing data. ----- */ + virtual void WriteImageInformation(bool keepOfStream); + + virtual void WriteImageInformation() { + WriteImageInformation(false); + } + + virtual bool CanWriteFile(const char* filename); + + virtual void Write(const void* buffer); + +protected: + std::string m_BinaryFileName; + int m_BinaryFileSkip; + + static char* edf_findInHeader( char* header, const char* key ); + + /* List of EDF supported datatypes + */ + enum DataType { + U_CHAR_DATATYPE = 0, CHAR_DATATYPE, // 8 bits = 1 B + U_SHORT_DATATYPE, SHORT_DATATYPE, // 16 bits = 2 B + U_INT_DATATYPE, INT_DATATYPE, // 32 bits = 4 B + U_L_INT_DATATYPE, L_INT_DATATYPE, // 32 bits = 4 B + FLOAT_DATATYPE, DOUBLE_DATATYPE, // 4 B, 8 B + UNKNOWN_DATATYPE = -1 + }; + + /* Note - compatibility: + Unsigned8 = 1, Signed8, Unsigned16, Signed16, + Unsigned32, Signed32, Unsigned64, Signed64, + FloatIEEE32, DoubleIEEE64 + */ + + /*************************************************************************** + * Tables + ***************************************************************************/ + + // table key-value structure + struct table { + const char *key; + int value; + }; + + struct table3 { + const char *key; + int value; + short sajzof; + }; + + /* Returns index of the table tbl whose key matches the beginning of the + * search string search_str. + * It returns index into the table or -1 if there is no match. + */ + static int + lookup_table_nth( const struct table *tbl, const char *search_str ) + { + int k = -1; + + while (tbl[++k].key) + if (tbl[k].key && !strncmp(search_str, tbl[k].key, strlen(tbl[k].key) ) ) + return k; + return -1; // not found + } + + static int + lookup_table3_nth( const struct table3 *tbl, const char *search_str ) + { + int k = -1; + + while (tbl[++k].key) + if (tbl[k].key && !strncmp(search_str, tbl[k].key, strlen(tbl[k].key) ) ) + return k; + return -1; // not found + } + + ///* Orientation of axes of the raster, as the binary matrix is saved in + // * the file. (Determines the scanning direction, or the "fastest" index + // * of the matrix in the data file.) + // */ + //enum EdfRasterAxes { + // RASTER_AXES_XrightYdown, // matricial format: rows, columns + // RASTER_AXES_XrightYup // cartesian coordinate system + // // other 6 combinations not available (not needed until now) + //}; + + //static const struct table rasteraxes_table[] = + //{ + // { "XrightYdown", RASTER_AXES_XrightYdown }, + // { "XrightYup", RASTER_AXES_XrightYup }, + // { NULL, -1 } + //}; + +}; // end class EdfImageIO + +} // end namespace + +#endif diff --git a/common/rtkEdfImageIOFactory.cxx b/common/rtkEdfImageIOFactory.cxx new file mode 100644 index 0000000..887c763 --- /dev/null +++ b/common/rtkEdfImageIOFactory.cxx @@ -0,0 +1,29 @@ +/*========================================================================= + * + * Copyright RTK Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ + +#include "rtkEdfImageIOFactory.h" + +//==================================================================== +rtk::EdfImageIOFactory::EdfImageIOFactory() +{ + this->RegisterOverride("itkImageIOBase", + "EdfImageIO", + "Edf Image IO", + 1, + itk::CreateObjectFunction::New() ); +} diff --git a/common/rtkEdfImageIOFactory.h b/common/rtkEdfImageIOFactory.h new file mode 100644 index 0000000..33eb518 --- /dev/null +++ b/common/rtkEdfImageIOFactory.h @@ -0,0 +1,77 @@ +/*========================================================================= + * + * Copyright RTK Consortium + * + * Licensed under the Apache License, Version 2.0 (the "License"); + * you may not use this file except in compliance with the License. + * You may obtain a copy of the License at + * + * http://www.apache.org/licenses/LICENSE-2.0.txt + * + * Unless required by applicable law or agreed to in writing, software + * distributed under the License is distributed on an "AS IS" BASIS, + * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. + * See the License for the specific language governing permissions and + * limitations under the License. + * + *=========================================================================*/ + +#ifndef __rtkEdfImageIOFactory_h +#define __rtkEdfImageIOFactory_h + +#include "rtkEdfImageIO.h" +#include +#include +#include + +namespace rtk { + +/** \class EdfImageIOFactory + * \brief ITK factory for Edf file I/O. + * + * \author Simon Rit + */ +class EdfImageIOFactory : public itk::ObjectFactoryBase +{ +public: + /** Standard class typedefs. */ + typedef EdfImageIOFactory Self; + typedef itk::ObjectFactoryBase Superclass; + typedef itk::SmartPointer Pointer; + typedef itk::SmartPointer ConstPointer; + + /** Class methods used to interface with the registered factories. */ + const char* GetITKSourceVersion(void) const { + return ITK_SOURCE_VERSION; + } + + const char* GetDescription(void) const { + return "Edf ImageIO Factory, allows the loading of Edf images into insight"; + } + + /** Method for class instantiation. */ + itkFactorylessNewMacro(Self); + + /** Run-time type information (and related methods). */ + itkTypeMacro(EdfImageIOFactory, ObjectFactoryBase); + + /** Register one factory of this type */ + static void RegisterOneFactory(void) { + ObjectFactoryBase::RegisterFactory( Self::New() ); + } + +protected: + EdfImageIOFactory(); + ~EdfImageIOFactory() {} + + typedef EdfImageIOFactory myProductType; + const myProductType* m_MyProduct; +private: + EdfImageIOFactory(const Self&); //purposely not implemented + void operator=(const Self&); //purposely not implemented + +}; + +} // end namespace + +#endif -- 2.47.1