Go to the documentation of this file.00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033 #include <Pipeline/cmtkSlicer.h>
00034
00035 #include <Base/cmtkAffineXform.h>
00036 #include <Base/cmtkUniformVolumeInterpolator.h>
00037 #include <Base/cmtkLinearInterpolator.h>
00038 #include <Base/cmtkCubicInterpolator.h>
00039 #include <Base/cmtkSincInterpolator.h>
00040 #include <Base/cmtkNearestNeighborInterpolator.h>
00041 #include <Base/cmtkUniformVolumeInterpolatorPartialVolume.h>
00042
00043 namespace
00044 cmtk
00045 {
00046
00049
00050 Slicer::Slicer()
00051 {
00052 ApplyWarp = false;
00053 this->m_Plane = NULL;
00054 InterpolationMode = cmtk::Interpolators::LINEAR;
00055 ZoomFactor = 1;
00056 TempImage = NULL;
00057 }
00058
00059 Slicer::~Slicer ()
00060 {
00061 if ( this->m_Plane )
00062 this->m_Plane->Delete();
00063 if ( TempImage )
00064 TempImage->Delete();
00065 }
00066
00067 void
00068 Slicer::SetPlane( Plane *const plane )
00069 {
00070 this->ReplaceObject( this->m_Plane, plane );
00071 }
00072
00073 long
00074 Slicer::Update()
00075 {
00076 this->CheckInputForUpdate( this->m_Plane );
00077 this->CheckInputForUpdate( Input );
00078 return this->ExecuteIfNecessary();
00079 }
00080
00081 void
00082 Slicer::Execute()
00083 {
00084 if ( Input == NULL ) return;
00085 if ( this->m_Plane == NULL ) return;
00086 UniformVolume::SmartPtr volume = this->Input->GetVolume();
00087 if ( ! volume ) return;
00088
00089 WarpXform::SmartPtr warpXform = Input->GetWarpXform();
00090
00091 Image *image = NULL;
00092
00093 if ( ZoomFactor == 1 )
00094 image = this->GetOutput();
00095 else {
00096 if ( ! TempImage ) {
00097 TempImage = Image::New();
00098 }
00099 image = TempImage;
00100 }
00101 image->CopyStructure( this->m_Plane );
00102
00103 image->SetDataType( volume->GetData()->GetType() );
00104 TypedArray::SmartPtr data = image->GetData();
00105
00106 FixedVector<3,Types::Coordinate> p( this->m_Plane->GetOrigin() );
00107 FixedVector<3,Types::Coordinate> dirX( this->m_Plane->GetDirectionX() );
00108 FixedVector<3,Types::Coordinate> dirY( this->m_Plane->GetDirectionY() );
00109
00110 dirX *= (this->m_Plane->GetSpacing()[0] / dirX.RootSumOfSquares());
00111 dirY *= (this->m_Plane->GetSpacing()[1] / dirY.RootSumOfSquares());
00112
00113 const unsigned int *dims = this->m_Plane->GetDims();
00114
00115 cmtk::UniformVolumeInterpolatorBase::SmartPtr interpolator;
00116 switch ( this->InterpolationMode )
00117 {
00118 default:
00119 case cmtk::Interpolators::LINEAR:
00120 {
00121 typedef cmtk::UniformVolumeInterpolator<cmtk::Interpolators::Linear> TInterpolator;
00122 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00123 break;
00124 }
00125 case cmtk::Interpolators::NEAREST_NEIGHBOR:
00126 {
00127 typedef cmtk::UniformVolumeInterpolator<cmtk::Interpolators::NearestNeighbor> TInterpolator;
00128 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00129 break;
00130 }
00131 case cmtk::Interpolators::CUBIC:
00132 {
00133 typedef cmtk::UniformVolumeInterpolator<cmtk::Interpolators::Cubic> TInterpolator;
00134 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00135 break;
00136 }
00137 case cmtk::Interpolators::COSINE_SINC:
00138 {
00139 typedef cmtk::UniformVolumeInterpolator< cmtk::Interpolators::CosineSinc<> > TInterpolator;
00140 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00141 break;
00142 }
00143 case cmtk::Interpolators::PARTIALVOLUME:
00144 {
00145 typedef cmtk::UniformVolumeInterpolatorPartialVolume TInterpolator;
00146 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00147 break;
00148 }
00149 }
00150
00151 if ( (!ApplyWarp) || ( !warpXform) )
00152 {
00153 AffineXform::SmartPtr affineXform( NULL );
00154 AffineXform::SmartPtr originalXform = Input->GetAffineXform();
00155 if ( originalXform ) affineXform = originalXform;
00156
00157 if ( affineXform )
00158 {
00159 affineXform->ApplyInPlace( dirX += p );
00160 affineXform->ApplyInPlace( dirY += p );
00161 affineXform->ApplyInPlace( p );
00162 dirX -= p;
00163 dirY -= p;
00164 }
00165
00166 size_t index = 0;
00167 Types::DataItem value;
00168 Vector3D rowOffset;
00169 for ( unsigned int y = 0; y<dims[1]; ++y, p += dirY )
00170 {
00171 rowOffset = p;
00172 for ( unsigned int x = 0; x<dims[0]; ++x, p += dirX, ++index )
00173 {
00174 if ( interpolator->GetDataAt( p, value ) )
00175 data->Set( value, index );
00176 else
00177 data->SetPaddingAt( index );
00178 }
00179 p = rowOffset;
00180 }
00181 }
00182 else
00183 {
00184 SplineWarpXform::SmartPtr splineWarpXform = SplineWarpXform::SmartPtr::DynamicCastFrom( warpXform );
00185 if ( splineWarpXform )
00186 {
00187 this->ExecuteSplineWarp( data, splineWarpXform, dims, p, dirX, dirY );
00188 }
00189 else
00190 {
00191 fputs( "Coding error: Unsupported WarpXform descendant encountered or dynamic_cast() failed.\n", stderr );
00192 }
00193 }
00194
00196 if ( ZoomFactor != 1 )
00197 {
00198 Image *output = this->GetOutput();
00199
00200 unsigned int dims[2];
00201 image->GetDims( dims );
00202
00203 dims[0] = static_cast<unsigned int>( dims[0] * ZoomFactor );
00204 dims[1] = static_cast<unsigned int>( dims[1] * ZoomFactor );
00205
00206 Types::Coordinate spacing[2];
00207 image->GetSpacing( spacing );
00208 spacing[0] /= ZoomFactor;
00209 spacing[1] /= ZoomFactor;
00210
00211 output->SetDims( dims );
00212 output->SetSpacing( spacing );
00213 output->SetDataType( image->GetDataType() );
00214
00215 TypedArray::SmartPtr outData = output->GetData();
00216
00217 unsigned int offset = 0;
00218 for ( unsigned int y = 0; y < dims[1]; ++y )
00219 for ( unsigned int x = 0; x < dims[0]; ++x, ++offset )
00220 {
00221 outData->Set( image->GetDataAt( x * spacing[0], y * spacing[1] ), offset );
00222 }
00223 }
00224
00225 this->UpdateExecuteTime();
00226 }
00227
00228 void
00229 Slicer::ExecuteSplineWarp
00230 ( TypedArray::SmartPtr& data, const SplineWarpXform* warpXform,
00231 const unsigned int* dims, const Vector3D& offset,
00232 const Vector3D& dX, const Vector3D& dY )
00233 {
00234 UniformVolume::SmartPtr volume = this->Input->GetVolume();
00235 if ( ! volume) return;
00236
00237 #ifdef DEBUG
00238 puts( "Slicer: Applying spline warp." );
00239 #endif
00240 unsigned int index = 0;
00241 Types::DataItem value;
00242 Vector3D pDeformed, rowOffset;
00243 Vector3D p( offset );
00244
00245 cmtk::UniformVolumeInterpolatorBase::SmartPtr interpolator;
00246 switch ( this->InterpolationMode )
00247 {
00248 default:
00249 case cmtk::Interpolators::LINEAR:
00250 {
00251 typedef cmtk::UniformVolumeInterpolator<cmtk::Interpolators::Linear> TInterpolator;
00252 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00253 break;
00254 }
00255 case cmtk::Interpolators::NEAREST_NEIGHBOR:
00256 {
00257 typedef cmtk::UniformVolumeInterpolator<cmtk::Interpolators::NearestNeighbor> TInterpolator;
00258 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00259 break;
00260 }
00261 case cmtk::Interpolators::CUBIC:
00262 {
00263 typedef cmtk::UniformVolumeInterpolator<cmtk::Interpolators::Cubic> TInterpolator;
00264 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00265 break;
00266 }
00267 case cmtk::Interpolators::COSINE_SINC:
00268 {
00269 typedef cmtk::UniformVolumeInterpolator< cmtk::Interpolators::CosineSinc<> > TInterpolator;
00270 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00271 break;
00272 }
00273 case cmtk::Interpolators::PARTIALVOLUME:
00274 {
00275 typedef cmtk::UniformVolumeInterpolatorPartialVolume TInterpolator;
00276 interpolator = cmtk::UniformVolumeInterpolatorBase::SmartPtr( new TInterpolator( *volume ) );
00277 break;
00278 }
00279 }
00280
00281 for ( unsigned int y = 0; y<dims[1]; ++y, p += dY )
00282 {
00283 rowOffset = p;
00284 for ( unsigned int x = 0; x<dims[0]; ++x, p += dX, ++index )
00285 {
00286 pDeformed = p;
00287 warpXform->ApplyInPlace( pDeformed );
00288 if ( interpolator->GetDataAt( pDeformed, value ) )
00289 data->Set( value, index );
00290 else
00291 data->SetPaddingAt( index );
00292 }
00293 p = rowOffset;
00294 }
00295 }
00296
00297 }