Commit 162a8787 authored by Mathias Goldau's avatar Mathias Goldau
Browse files

[CHANGE] Now the Nifitreader also supports reading the bvalues. Additionally...

[CHANGE] Now the Nifitreader also supports reading the bvalues. Additionally the source for reading the bvec and bval is refactored into two functions where the expected file format is given in the doc-strings.
parent 819abadb
......@@ -94,6 +94,121 @@ WMatrix< double > WReaderNIfTI::convertMatrix( const mat44& in )
return out;
}
WReaderNIfTI::GradVec WReaderNIfTI::readGradientsIfAvailable( unsigned int vDim )
{
std::string gradientFileName = m_fname;
std::string suffix = getSuffix( m_fname );
if( suffix == ".gz" )
{
WAssert( gradientFileName.length() > 3, "" );
gradientFileName.resize( gradientFileName.length() - 3 );
suffix = getSuffix( gradientFileName );
}
WAssert( suffix == ".nii", "Input file is not a nifti file." );
WAssert( gradientFileName.length() > 4, "" );
gradientFileName.resize( gradientFileName.length() - 4 );
gradientFileName += ".bvec";
GradVec result; // incase of error return NULL_ptr
// check if the file exists
std::ifstream i( gradientFileName.c_str() );
if( i.bad() || !i.is_open() )
{
if( i.is_open() )
{
i.close();
}
wlog::debug( "WReaderNIfTI" ) << "Could not find gradient file expected at: \"" << gradientFileName << "\", skipping this.";
}
else
{
wlog::debug( "WReaderNIfTI" ) << "Found b-vectors file: " << gradientFileName << " will try reading...";
result = GradVec( new std::vector< WVector3d >( vDim ) );
// the file should contain the x-coordinates in line 0, y-coordinates in line 1,
// z-coordinates in line 2
for( unsigned int j = 0; j < 3; ++j )
{
for( unsigned int k = 0; k < vDim; ++k )
{
i >> result->operator[] ( k )[ j ];
}
}
bool success = !i.eof();
i.close();
if( !success )
{
wlog::error( "WReaderNIfTI" ) << "Error while reading gradient file: did not contain enough gradients: " << result->size();
return GradVec(); // return Null_ptr
}
else
{
wlog::debug( "WReaderNIfTI" ) << "Successfully loaded " << result->size() << " gradients";
}
}
return result;
}
WReaderNIfTI::BValues WReaderNIfTI::readBValuesIfAvailable( unsigned int vDim )
{
std::string bvaluesFileName = m_fname;
std::string suffix = getSuffix( m_fname );
if( suffix == ".gz" )
{
WAssert( bvaluesFileName.length() > 3, "" );
bvaluesFileName.resize( bvaluesFileName.length() - 3 );
suffix = getSuffix( bvaluesFileName );
}
WAssert( suffix == ".nii", "Input file is not a nifti file." );
WAssert( bvaluesFileName.length() > 4, "" );
bvaluesFileName.resize( bvaluesFileName.length() - 4 );
bvaluesFileName += ".bval";
BValues result; // return NULL_ptr in case of error
// check if the file exists
std::ifstream i( bvaluesFileName.c_str() );
if( i.bad() || !i.is_open() )
{
if( i.is_open() )
{
i.close();
}
wlog::debug( "WReaderNIfTI" ) << "Could not find b-values file expected at: \"" << bvaluesFileName << "\", skipping this.";
}
else
{
//read b-values
char value[ 8 ];
// there should be 3 * vDim values in the file
result = BValues( new std::vector< float >( vDim * 3, 0 ) );
size_t numValues = 0;
while( i.good() && !i.eof() )
{
i.getline( value, 8 );
float fVal;
std::istringstream( value ) >> fVal;
( *result )[ numValues ] = fVal;
if( numValues > vDim * 3 )
{
wlog::error( "WReaderNIfTI" ) << "Too many b-Values: " << numValues << " but expected: " << vDim * 3;
return BValues(); // return Null_ptr
}
numValues++;
}
i.close();
wlog::debug( "WReaderNIfTI" ) << "Found b-values file and loaded " << result->size() << " values.";
}
return result;
}
boost::shared_ptr< WDataSet > WReaderNIfTI::load( DataSetType dataSetType )
{
boost::shared_ptr< nifti_image > filedata( nifti_image_read( m_fname.c_str(), 1 ), &nifti_image_free );
......@@ -358,55 +473,31 @@ boost::shared_ptr< WDataSet > WReaderNIfTI::load( DataSetType dataSetType )
{
wlog::debug( "WReaderNIfTI" ) << "Load as WDataSetRawHARDI";
std::string gradientFileName = m_fname;
std::string suffix = getSuffix( m_fname );
if( suffix == ".gz" )
{
WAssert( gradientFileName.length() > 3, "" );
gradientFileName.resize( gradientFileName.length() - 3 );
suffix = getSuffix( gradientFileName );
}
WAssert( suffix == ".nii", "Input file is not a nifti file." );
GradVec newGradients = readGradientsIfAvailable( vDim );
BValues newBValues = readBValuesIfAvailable( vDim );
WAssert( gradientFileName.length() > 4, "" );
gradientFileName.resize( gradientFileName.length() - 4 );
gradientFileName += ".bvec";
// check if the file exists
std::ifstream i( gradientFileName.c_str() );
if( i.bad() || !i.is_open() )
if( !newGradients )
{
if( i.is_open() )
{
i.close();
}
// cannot find the appropriate gradient vectors, build a dataSetSingle instead of hardi
// cannot find the appropriate gradient vectors and/or bvalues, build a dataSetSingle instead of hardi
newDataSet = boost::shared_ptr< WDataSet >( new WDataSetSingle( newValueSet, newGrid ) );
wlog::debug( "WReaderNIfTI" ) << "Could not find corresponding gradients file \"" << gradientFileName.c_str() <<
"\", loading as WDataSetSingle instead.";
wlog::debug( "WReaderNIfTI" ) << "No gradients given. See above for expected filename. Loading as WDataSetSingle instead.";
}
else
{
// read gradients, there should be 3 * vDim values in the file
typedef std::vector< WVector3d > GradVec;
boost::shared_ptr< GradVec > newGradients( new GradVec( vDim ) );
// the file should contain the x-coordinates in line 0, y-coordinates in line 1,
// z-coordinates in line 2
for( unsigned int j = 0; j < 3; ++j )
if( !newBValues )
{
for( unsigned int k = 0; k < vDim; ++k )
{
i >> newGradients->operator[] ( k )[ j ];
}
newDataSet = boost::shared_ptr< WDataSet >( new WDataSetRawHARDI( newValueSet, newGrid, newGradients ) );
}
bool success = !i.eof();
i.close();
WAssert( success, "Gradient file did not contain enough gradients." );
else
{
newDataSet = boost::shared_ptr< WDataSet >( new WDataSetRawHARDI( newValueSet, newGrid, newGradients, newBValues ) );
}
}
newDataSet = boost::shared_ptr< WDataSet >( new WDataSetRawHARDI( newValueSet, newGrid, newGradients ) );
if( newBValues && newGradients && ( newBValues->size() > 1 && newBValues->size() != newGradients->size() * 3 ) )
{
wlog::warn( "WReaderNIfTI" ) << "Be careful: there are " << newBValues->size() / 3 << " b-values but only "
<< newGradients->size() << " gradients.";
}
}
else if( filedata->intent_code == NIFTI_INTENT_SYMMATRIX )
......
......@@ -86,6 +86,75 @@ public:
WMatrix< double > getQFormTransform() const;
protected:
/**
* Shorthand type for a vector of gradients.
*/
typedef boost::shared_ptr< std::vector< WVector3d > > GradVec;
/** Shorthand type for a vector of bvalues. The comopents are typically given as floats, as partical b-values start
* at 0 and are at maximum around 17000 (see: https://dx.doi.org/10.1002/mrm.20642).
*/
typedef boost::shared_ptr< std::vector< float > > BValues;
/**
* Reads the additional bval file if available. The file format should be (ASCII file):
* - for each (x,y,z) component the bvalue must be repeated in such a form:
*
* bval for_x-component of the first gradient followed by a newline
* bval for x-component of the second gradient followed by a newline
* ...
* bval for x-compoennt of the last gradient followed by a newline
* bval for_y-component of the first gradient followed by a newline
* ...
* bval for y-compoennt of the last gradient followed by a newline
* bval for_z-component of the first gradient followed by a newline
* ...
* bval for z-compoennt of the last gradient followed by a newline
*
* For example: three gradients, with bvalues: 0, 1000, and 1000.
*
* 0
* 1000
* 1000
* 0
* 1000
* 1000
* 0
* 1000
* 1000
*
* for three gradients, where first all bvalues for all x-components must be given, each in a separate line, then all bavlues
* for all y-components, and finally the same for the z-components. Thus there should be three times of the number of
* gradients be present in such a file.
*
* \param vDim Nifti dimension
*
* \return Null ptr if not available or in case of error during read, otherwise the pointer of the vector of BValues.
*/
BValues readBValuesIfAvailable( unsigned int vDim );
/**
* Reads the additional gradient file if available. The file format should be (ASCII file):
* There are three lines, each containing the x-, y-, and z-components of the gradients.
* The first line contains all x-components of all gradients. Thus each lines should have the same number of values, each
* separated by whitespace.
*
* x-comp._of_the_first_gradient x-comp._of_the_second_gradient ... x-comp._of_the_last_gradient followed by a newline
* y-comp._of_the_first_gradient y-comp._of_the_second_gradient ... y-comp._of_the_last_gradient followed by a newline
* z-comp._of_the_first_gradient z-comp._of_the_second_gradient ... z-comp._of_the_last_gradient followed by a newline
*
* For example: three gradients:
*
* 0 -0.530466 -0.294864
* 0 -0.558441 -0.285613
* 0 0.637769 -0.911856
*
* \param vDim Nifti dimension
*
* \return Null ptr if not available or in case of error during read, otherwise the pointer of the vector of gradients.
*/
GradVec readGradientsIfAvailable( unsigned int vDim );
private:
/**
* This function allows one to copy the data given as a T*
......
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment