From a8d55192c4bf80942ee940919b6c22a299493557 Mon Sep 17 00:00:00 2001
From: Mickael Savinaud <mickael.savinaud@c-s.fr>
Date: Mon, 7 Mar 2011 16:43:23 +0100
Subject: [PATCH] ENH: clean code to read complex data
---
Code/IO/otbGDALImageIO.cxx | 222 ++++++++----------
Code/IO/otbImageFileReader.txx | 16 +-
.../ITK/Code/IO/itkConvertPixelBuffer.txx | 1 -
3 files changed, 98 insertions(+), 141 deletions(-)
diff --git a/Code/IO/otbGDALImageIO.cxx b/Code/IO/otbGDALImageIO.cxx
index b09b44a31b..5440596469 100644
--- a/Code/IO/otbGDALImageIO.cxx
+++ b/Code/IO/otbGDALImageIO.cxx
@@ -38,6 +38,55 @@
namespace otb
{
+template<class InputType>
+void printOutputData(InputType *pData, int nbBands, int nbPixelToRead)
+{
+ for (unsigned int itPxl = 0; itPxl < (unsigned int) (nbPixelToRead * nbBands); itPxl++)
+ {
+ std::cout << "Buffer["<< itPxl << "] = " << *(pData + itPxl) << std::endl;
+ }
+};
+
+void printDataBuffer(unsigned char *pData, GDALDataType pxlType, int nbBands, int nbPixelToRead)
+{
+ if (pxlType == GDT_Int16)
+ {
+ printOutputData( static_cast<short*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else if (pxlType == GDT_Int32)
+ {
+ printOutputData( static_cast<int*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else if (pxlType == GDT_Float32)
+ {
+ printOutputData( static_cast<float*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else if (pxlType == GDT_Float64)
+ {
+ printOutputData( static_cast<double*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else if (pxlType == GDT_CInt16)
+ {
+ printOutputData( static_cast<std::complex<short>*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else if (pxlType == GDT_CInt32)
+ {
+ printOutputData( static_cast<std::complex<int>*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else if (pxlType == GDT_CFloat32)
+ {
+ printOutputData( static_cast<std::complex<float>*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else if (pxlType == GDT_CFloat64)
+ {
+ printOutputData( static_cast<std::complex<double>*>( static_cast<void*>(pData) ), nbBands, nbPixelToRead);
+ }
+ else
+ {
+ std::cerr << "Pixel type unknown" << std::endl;
+ }
+};
+
// only two states : the Pointer is Null or GetDataSet() returns a valid dataset
class GDALDatasetWrapper : public itk::LightObject
@@ -247,7 +296,6 @@ void GDALImageIO::Read(void* buffer)
&& (m_PxType != GDT_CFloat32)
&& (m_PxType != GDT_CFloat64))
{
- std::cout << "*** GDT_INT32 AND GDT_INT16 CASE ***" << std::endl;
int pixelOffset = m_BytePerPixel * m_NbBands;
int lineOffset = m_BytePerPixel * m_NbBands * lNbColumns;
int bandOffset = m_BytePerPixel;
@@ -257,19 +305,16 @@ void GDALImageIO::Read(void* buffer)
std::streamoff nbBytes = static_cast<std::streamoff>(m_NbBands) * static_cast<std::streamoff>(nbPixelToRead) * static_cast<std::streamoff>(m_BytePerPixel);
unsigned char *pBufferTemp = new unsigned char[static_cast<unsigned int>(nbBytes)];
-
// keep it for the moment
- std::cout << "*** GDALimageIO::Read: nominal case ***"<<std::endl;
- std::cout << "Paremeters RasterIO :" \
- << ", indX = " << lFirstColumn \
- << ", indY = " << lFirstLine \
- << ", sizeX = " << lNbColumns \
- << ", sizeY = " << lNbLines \
- << ", GDAL Data Type = " << GDALGetDataTypeName(m_PxType) \
- << ", pixelOffset = " << pixelOffset \
- << ", lineOffset = " << lineOffset
- << ", bandOffset = " << bandOffset << std::endl;
-
+ otbMsgDevMacro(<< "Parameters RasterIO (case CInt and CShort):"
+ << "\n, indX = " << lFirstColumn
+ << "\n, indY = " << lFirstLine
+ << "\n, sizeX = " << lNbColumns
+ << "\n, sizeY = " << lNbLines
+ << "\n, GDAL Data Type = " << GDALGetDataTypeName(m_PxType)
+ << "\n, pixelOffset = " << pixelOffset
+ << "\n, lineOffset = " << lineOffset
+ << "\n, bandOffset = " << bandOffset);
CPLErr lCrGdal = m_Dataset->GetDataSet()->RasterIO(GF_Read,
lFirstColumn,
@@ -290,58 +335,45 @@ void GDALImageIO::Read(void* buffer)
if (lCrGdal == CE_Failure)
{
itkExceptionMacro(<< "Error while reading image (GDAL format) " << m_FileName );
+ return;
}
+ //std::cout << "RAW BUFFER:" <<std::endl;
+ //printDataBuffer(pBufferTemp, m_PxType, m_NbBands, lNbColumns*lNbLines);
- // Reorganiser le buffer pour le faire fitter avec p qui doit contenir des GDT_Float64
+ // Convert the buffer to GDT_Float64 type
typedef std::complex<double> RealType;
typedef double ScalarRealType;
- RealType *realPxlValue = new RealType[nbPixelToRead * m_NbBands];
if (m_PxType == GDT_CInt32)
{
- std::cout << "Convert from input File from GDT_CInt32 to GDT_CFloat64" << std::endl;
+ //std::cout << "Convert input File from GDT_CInt32 to GDT_CFloat64" << std::endl;
typedef std::complex<int> ComplexIntType;
- ComplexIntType *pxlValue = new ComplexIntType[nbPixelToRead * m_NbBands];
for (unsigned int itPxl = 0; itPxl < (unsigned int) (nbPixelToRead * m_NbBands); itPxl++)
{
- pxlValue[itPxl] = *(static_cast<ComplexIntType*>( static_cast<void*>(pBufferTemp)) + itPxl );
- std::cout << "originalBuffer["<< itPxl << "] = " << pxlValue[itPxl] << std::endl;
+ ComplexIntType pxlValue = *(static_cast<ComplexIntType*>( static_cast<void*>(pBufferTemp)) + itPxl );
- RealType pxlValueReal( static_cast<ScalarRealType>(pxlValue[itPxl].real()), static_cast<ScalarRealType>(pxlValue[itPxl].imag()) );
- realPxlValue[itPxl] = pxlValueReal;
- std::cout << "newBuffer["<< itPxl << "] = " << realPxlValue[itPxl] << std::endl;
+ RealType pxlValueReal( static_cast<ScalarRealType>(pxlValue.real()), static_cast<ScalarRealType>(pxlValue.imag()) );
- memcpy((void*) (&(p[itPxl*sizeof(RealType)])), (const void*) (&(realPxlValue[itPxl])), (size_t) (sizeof(RealType)));
-
- RealType pxlValueRealOut = *(static_cast<RealType*>( static_cast<void*>(p)) + itPxl );
- std::cout << "p["<< itPxl << "] = " << pxlValueRealOut << std::endl;
+ memcpy((void*) (&(p[itPxl*sizeof(RealType)])), (const void*) (&(pxlValueReal)), (size_t) (sizeof(RealType)));
}
- delete[] pxlValue;
}
else if (m_PxType == GDT_CInt16)
{
- std::cout << "Convert from input File from GDT_CInt16 to GDT_CFloat64" << std::endl;
+ //std::cout << "Convert input File from GDT_CInt16 to GDT_CFloat64" << std::endl;
typedef std::complex<short> ComplexShortType;
- ComplexShortType *pxlValue = new ComplexShortType[nbPixelToRead * m_NbBands];
+
for (unsigned int itPxl = 0; itPxl < (unsigned int) (nbPixelToRead * m_NbBands); itPxl++)
{
- pxlValue[itPxl] = *(static_cast<ComplexShortType*>( static_cast<void*>(pBufferTemp)) + itPxl );
- std::cout << "originalBuffer["<< itPxl << "] = " << pxlValue[itPxl] << std::endl;
+ ComplexShortType pxlValue = *(static_cast<ComplexShortType*>( static_cast<void*>(pBufferTemp)) + itPxl );
- RealType pxlValueReal( static_cast<ScalarRealType>(pxlValue[itPxl].real()), static_cast<ScalarRealType>(pxlValue[itPxl].imag()) );
- realPxlValue[itPxl] = pxlValueReal;
- std::cout << "newBuffer["<< itPxl << "] = " << realPxlValue[itPxl] << std::endl;
+ RealType pxlValueReal( static_cast<ScalarRealType>(pxlValue.real()), static_cast<ScalarRealType>(pxlValue.imag()) );
- memcpy((void*) (&(p[itPxl*sizeof(RealType)])), (const void*) (&(realPxlValue[itPxl])), (size_t) (sizeof(RealType)));
-
- RealType pxlValueRealOut = *(static_cast<RealType*>( static_cast<void*>(p)) + itPxl );
- std::cout << "p["<< itPxl << "] = " << pxlValueRealOut << std::endl;
+ memcpy((void*) (&(p[itPxl*sizeof(RealType)])), (const void*) (&(pxlValueReal)), (size_t) (sizeof(RealType)));
}
- delete[] pxlValue;
}
- delete[] realPxlValue;
-
+ //std::cout << "CONVERTED BUFFER:" <<std::endl;
+ //printDataBuffer(p, GDT_CFloat64, m_NbBands, lNbColumns*lNbLines);
}
// In the indexed case, one has to retrieve the index image and the
@@ -390,61 +422,38 @@ void GDALImageIO::Read(void* buffer)
}
else
{
- // Nominal case
- std::cout << "*** NOMINAL CASE ***" << std::endl;
+ /******** Nominal case ***********/
int pixelOffset = m_BytePerPixel * m_NbBands;
int lineOffset = m_BytePerPixel * m_NbBands * lNbColumns;
int bandOffset = m_BytePerPixel;
int nbBands = m_NbBands;
// In some cases, we need to change some parameters for RasterIO
-
- // if the file is complex and the reader is based on a vector of scalar,
- // output 2 times the number of bands, with real and imaginary parts interleaved
- /*if(GDALDataTypeIsComplex(m_PxType) && !m_IsComplex && m_IsVectorImage)
- {
- // ImageIO NbComponents is set to 2 * m_NbBands
- // m_BytePerPixel is already sizeof(std::complex<m_PxType>) / 2
- pixelOffset = m_BytePerPixel * this->GetNumberOfComponents();
- lineOffset = pixelOffset * lNbColumns;
- bandOffset = 2 * m_BytePerPixel;
- nbBands = this->GetNumberOfComponents() / 2;
- }*/
-
- // if the file is scalar with only one band and the reader is based on a vector of complex
- /*if(!GDALDataTypeIsComplex(m_PxType) && (m_NbBands == 1) && m_IsComplex && m_IsVectorImage )
- {
- pixelOffset = m_BytePerPixel / 2;
- lineOffset =pixelOffset * lNbColumns;
- bandOffset = m_BytePerPixel / 2;
- }*/
- if(!GDALDataTypeIsComplex(m_PxType) /*&& ((unsigned int)m_NbBands == this->GetNumberOfComponents())*/ && m_IsComplex && m_IsVectorImage && (m_NbBands > 1))
+ if(!GDALDataTypeIsComplex(m_PxType) && m_IsComplex && m_IsVectorImage && (m_NbBands > 1))
{
pixelOffset = m_BytePerPixel * 2 ;
lineOffset = pixelOffset * lNbColumns;
- bandOffset = m_BytePerPixel /*/ 2*/;
+ bandOffset = m_BytePerPixel;
}
// keep it for the moment
-
- std::cout << "*** GDALimageIO::Read: nominal case ***"<<std::endl;
- std::cout << "Paremeters RasterIO :" \
- << ", indX = " << lFirstColumn \
- << ", indY = " << lFirstLine \
- << ", sizeX = " << lNbColumns \
- << ", sizeY = " << lNbLines \
- << ", GDAL Data Type = " << GDALGetDataTypeName(m_PxType) \
- << ", pixelOffset = " << pixelOffset \
- << ", lineOffset = " << lineOffset
- << ", bandOffset = " << bandOffset << std::endl;
-
+ //otbMsgDevMacro(<< "Number of bands inside input file: " << m_NbBands);
+ otbMsgDevMacro(<< "Parameters RasterIO : \n"
+ << ", indX = " << lFirstColumn << "\n"
+ << ", indY = " << lFirstLine << "\n"
+ << ", sizeX = " << lNbColumns << "\n"
+ << ", sizeY = " << lNbLines << "\n"
+ << ", GDAL Data Type = " << GDALGetDataTypeName(m_PxType) << "\n"
+ << ", pixelOffset = " << pixelOffset << "\n"
+ << ", lineOffset = " << lineOffset << "\n"
+ << ", bandOffset = " << bandOffset);
CPLErr lCrGdal = m_Dataset->GetDataSet()->RasterIO(GF_Read,
lFirstColumn,
lFirstLine,
lNbColumns,
lNbLines,
- p, // pData
+ p,
lNbColumns,
lNbLines,
m_PxType,
@@ -458,21 +467,12 @@ void GDALImageIO::Read(void* buffer)
if (lCrGdal == CE_Failure)
{
itkExceptionMacro(<< "Error while reading image (GDAL format) " << m_FileName.c_str() << ".");
+ return;
}
- typedef std::complex<int> ComplexType;
- int nbPixelToRead = lNbColumns * lNbLines;
- ComplexType *pxlValue = new ComplexType[nbPixelToRead * m_NbBands];
- unsigned int count = 0;
- ComplexType expectedValue;
- for (unsigned int itPxl = 0; itPxl < (unsigned int) (nbPixelToRead * m_NbBands); itPxl++)
- {
- pxlValue[itPxl] = *(static_cast<ComplexType*>( static_cast<void*>(p)) + itPxl );
- std::cout << "loadBuffer["<< itPxl << "] = " << pxlValue[itPxl] << std::endl;
- }
+ //printDataBuffer(p, m_PxType, m_NbBands, lNbColumns*lNbLines);
}
}
-
bool GDALImageIO::GetSubDatasetInfo(std::vector<std::string> &names, std::vector<std::string> &desc)
{
// Note: we assume that the subdatasets are in order : SUBDATASET_ID_NAME, SUBDATASET_ID_DESC, SUBDATASET_ID+1_NAME, SUBDATASET_ID+1_DESC
@@ -682,13 +682,8 @@ void GDALImageIO::InternalReadImageInformation()
}
/******************************************************************/
- // Pixel Type always set to Scalar for GDAL ? maybe also to vector ?
-
- //Once all sorts of gdal complex image are handle, this won't be
- //necessary any more
-
// Set the pixel type with some special cases linked to the fact
- // we read some data with complex type
+ // we read some data with complex type.
if ( GDALDataTypeIsComplex(m_PxType) ) // Try to read data with complex type with GDAL
{
if ( !m_IsComplex && m_IsVectorImage )
@@ -706,7 +701,6 @@ void GDALImageIO::InternalReadImageInformation()
}
}
else // Try to read data with scalar type with GDAL
- //if ( !GDALDataTypeIsComplex(m_PxType) )
{
this->SetNumberOfComponents(m_NbBands);
if (this->GetNumberOfComponents() == 1)
@@ -719,39 +713,11 @@ void GDALImageIO::InternalReadImageInformation()
}
}
- /*if (GDALDataTypeIsComplex(m_PxType) && )
- {
- // we are reading a complex data set into an image where the pixel
- // type is Vector<real>: we have to double the number of component
- // for that to work
- std::cout << "GDALtypeIO = Complex and IFReader::InternalPixelType = Scalar and IFReader::PixelType = Vector" << std::endl;
- //if ( (m_PxType != GDT_CInt32) && (m_PxType != GDT_CInt16) )
- // m_BytePerPixel = m_BytePerPixel / 2;
- this->SetNumberOfComponents(m_NbBands*2);
- this->SetPixelType(VECTOR);
- }*/
-
- /*if (!GDALDataTypeIsComplex(m_PxType) && m_IsComplex && m_IsVectorImage && (m_NbBands > 1))
- {
- // we are reading a non-complex data set into an image where the pixel
- // type is Vector<complex>: we have to double the number of byte per pixel
- // for that to work
- std::cout << "GDALtypeIO = Scalar and IFReader::InternalPixelType = Complex and IFReader::PixelType = Vector" << std::endl;
- //if ( (m_PxType != GDT_CInt32) && (m_PxType != GDT_CInt16) )
- //m_BytePerPixel = m_BytePerPixel * 2;
- this->SetPixelType(VECTOR);
- }*/
-
- /* if (!GDALDataTypeIsComplex(m_PxType) && m_IsComplex && !m_IsVectorImage && ( (m_PxType == GDT_CInt32) || (m_PxType == GDT_CInt16) ) )
- {
- m_BytePerPixel = m_BytePerPixel / 2;
- }*/
-
/*** Parameters set by Internal Read function ***/
- otbMsgDevMacro( << " Pixel Type IFReader = " << GetPixelTypeAsString(this->GetPixelType()) )
- otbMsgDevMacro( << " Number of component IFReader = " << this->GetNumberOfComponents() )
- otbMsgDevMacro( << " Byte per pixel set = " << m_BytePerPixel )
- otbMsgDevMacro( << " Component Type set = " << GetComponentTypeAsString(this->GetComponentType()) );
+ otbMsgDevMacro( << "Pixel Type IFReader = " << GetPixelTypeAsString(this->GetPixelType()) )
+ otbMsgDevMacro( << "Number of component IFReader = " << this->GetNumberOfComponents() )
+ otbMsgDevMacro( << "Byte per pixel set = " << m_BytePerPixel )
+ otbMsgDevMacro( << "Component Type set = " << GetComponentTypeAsString(this->GetComponentType()) );
/*----------------------------------------------------------------------*/
/*-------------------------- METADATA ----------------------------------*/
diff --git a/Code/IO/otbImageFileReader.txx b/Code/IO/otbImageFileReader.txx
index c94d41d679..a72fe18dd1 100644
--- a/Code/IO/otbImageFileReader.txx
+++ b/Code/IO/otbImageFileReader.txx
@@ -92,23 +92,16 @@ ImageFileReader<TOutputImage>
output->SetBufferedRegion(output->GetRequestedRegion());
output->Allocate();
-//otbMsgDebugMacro( <<"ImageFileReader<TOutputImage>::GenerateData : ");
-//otbMsgDebugMacro( <<" output->GetRequestedRegion() : "<<output->GetRequestedRegion());
-
// Test if the file exist and if it can be open.
// An exception will be thrown otherwise.
this->TestFileExistanceAndReadability();
// Tell the ImageIO to read the file
- //
OutputImagePixelType *buffer =
output->GetPixelContainer()->GetBufferPointer();
this->m_ImageIO->SetFileName(this->m_FileName.c_str());
itk::ImageIORegion ioRegion(TOutputImage::ImageDimension);
-//otbMsgDebugMacro( <<" Avant ioRegion : "<<ioRegion);
-
-// itk::ImageIORegion ioRegionStreaming = output->GetRequestedRegion();
itk::ImageIORegion::SizeType ioSize = ioRegion.GetSize();
itk::ImageIORegion::IndexType ioStart = ioRegion.GetIndex();
@@ -148,8 +141,6 @@ ImageFileReader<TOutputImage>
ioRegion.SetSize(ioSize);
ioRegion.SetIndex(ioStart);
-//otbMsgDebugMacro( <<" Apres ioRegion : "<<ioRegion);
-
this->m_ImageIO->SetIORegion(ioRegion);
typedef itk::DefaultConvertPixelTraits<ITK_TYPENAME TOutputImage::IOPixelType> ConvertIOPixelTraits;
@@ -175,11 +166,12 @@ ImageFileReader<TOutputImage>
* static_cast<std::streamoff>(region.GetNumberOfPixels());
char * loadBuffer = new char[nbBytes];
- std::cout<< "*** IFReader : read with conversion ***" << std::endl;
- std::cout<< "size of Buffer to GDALImageIO::read = " << nbBytes << " = " \
+
+ otbMsgDevMacro(<< "size of Buffer to GDALImageIO::read = " << nbBytes << " = \n"
<< "ComponentSize ("<< this->m_ImageIO->GetComponentSize() << ") x " \
<< "Nb of Component (" << this->m_ImageIO->GetNumberOfComponents() << ") x " \
- << "Nb of Pixel to read (" << region.GetNumberOfPixels() << ")" << std::endl;
+ << "Nb of Pixel to read (" << region.GetNumberOfPixels() << ")" );
+
this->m_ImageIO->Read(loadBuffer);
this->DoConvertBuffer(loadBuffer, region.GetNumberOfPixels());
diff --git a/Utilities/ITK/Code/IO/itkConvertPixelBuffer.txx b/Utilities/ITK/Code/IO/itkConvertPixelBuffer.txx
index fe982f7f31..5976db4b70 100644
--- a/Utilities/ITK/Code/IO/itkConvertPixelBuffer.txx
+++ b/Utilities/ITK/Code/IO/itkConvertPixelBuffer.txx
@@ -726,7 +726,6 @@ SpecialCast(const std::complex<InputType>& in, const OutputType& itkNotUsed(dumm
typedef typename itk::NumericTraits<std::complex<InputType> >::RealType RealType;
typedef typename itk::NumericTraits<std::complex<InputType> >::ScalarRealType ScalarRealType;
- //RealType inReal = static_cast<RealType>(in);
RealType inReal( static_cast<ScalarRealType>(in.real()), static_cast<ScalarRealType>(in.imag()) );
return static_cast < OutputType >( vcl_abs(inReal) );
--
GitLab