[Insight-users] Vital DICOM information/tags when registering images?
Andriy Fedorov
fedorov at bwh.harvard.edu
Wed Sep 16 09:14:41 EDT 2009
Patrick,
Registration has the metric sampling component, which is initialized
with a random seed by default. If you want to have consistent results
you might need to initialize registration with the same seed.
Another guess is to check that there is no type cast involved when you
read the image in Matlab. This might result in loss of precision.
Not sure if any of these is the source of the problem for you, just
some ideas to investigate.
Andriy Fedorov
> Date: Wed, 16 Sep 2009 09:02:59 +0200
> From: Patrik Brynolfsson <patrik.brynolfsson at gmail.com>
> Subject: Re: [Insight-users] Vital DICOM information/tags when
> registering images?
> To: Bill Lorensen <bill.lorensen at gmail.com>, ITK
> <insight-users at itk.org>
> Message-ID:
> <8d0068770909160002r1e0aa86el32f78b0d027c9368 at mail.gmail.com>
> Content-Type: text/plain; charset="iso-8859-1"
>
> Hello Bill, thanx for answering!
>
> The filter I use are:
> Registration part (same for both cases):
>
> itkImageRegistrationMethod
> itkMattesMutualInformationImageToImageMetric
> itkLinearInterpolateImageFunction
> itkVersorRigid3DTransform
> itkVersorRigid3DTransformOptimizer
> itkCenteredTransformInitializer
>
> Reading images from memory:
>
> itkImportImageFilter
> itkChangeInformationImageFilter
>
> Reading images from DICOM:
>
> itkImageSeriesReader
> itkGDCMImageIO
> itkGDCMSeriesFileNames
>
> *The code to read images from memory is:*
>
> typedef itk::ImportImageFilter<PixelType,Dimension> ImportFilterType;
> ImportFilterType::Pointer fixedImageReader = ImportFilterType::New();
> ImportFilterType::Pointer movingImageReader = ImportFilterType::New();
>
> const mwSize* dimsFix = mxGetDimensions(prhs[0]);
> const mwSize* dimsMov = mxGetDimensions(prhs[1]);
>
> ImportFilterType::SizeType sizeFix;
> ImportFilterType::SizeType sizeMov;
> sizeFix[0] = dimsFix[0];
> sizeFix[1] = dimsFix[1];
> sizeFix[2] = dimsFix[2];
>
> sizeMov[0] = dimsMov[0];
> sizeMov[1] = dimsMov[1];
> sizeMov[2] = dimsMov[2];
>
> ImportFilterType::IndexType start;
> start.Fill( 0 );
> ImportFilterType::RegionType regionFix;
> regionFix.SetIndex( start );
> regionFix.SetSize( sizeFix );
> ImportFilterType::RegionType regionMov;
> regionMov.SetIndex( start );
> regionMov.SetSize( sizeMov );
> fixedImageReader->SetRegion( regionFix );
> movingImageReader->SetRegion(regionMov );
>
> //Import image from memory, passed from MATLAB
> fixedImageReader->SetImportPointer((PixelType*)mxGetData(prhs[0]),mxGetNumberOfElements(prhs[0]),false);
> movingImageReader->SetImportPointer((PixelType*)mxGetData(prhs[1]),mxGetNumberOfElements(prhs[1]),false);
> fixedImageReader->Update();
> movingImageReader->Update();
>
> *...and the code to set the relevant info is:*
>
> double *fixedSpacingInput = mxGetPr(prhs[2]);
> double *movingSpacingInput = mxGetPr(prhs[3]);
>
> double *fixedOrientationInput = mxGetPr(prhs[4]);
> double *movingOrientationInput = mxGetPr(prhs[5]);
>
> double *fixedOriginInput = mxGetPr(prhs[6]);
> double *movingOriginInput = mxGetPr(prhs[7]);
>
> **Code for casting input to correct types are skipped, not important for
> this question**
>
> typedef itk::ChangeInformationImageFilter< FixedImageType >
> ChangeInfoType;
> ChangeInfoType::Pointer changeFixedInfo = ChangeInfoType::New();
> ChangeInfoType::Pointer changeMovingInfo = ChangeInfoType::New();
>
> changeFixedInfo->SetInput(fixedImageReader->GetOutput());
> changeMovingInfo->SetInput(movingImageReader->GetOutput());
>
> changeFixedInfo->SetOutputSpacing(fixedSpacing);
> changeMovingInfo->SetOutputSpacing(movingSpacing);
> changeFixedInfo->SetChangeSpacing(true);
> changeMovingInfo->SetChangeSpacing(true);
>
> changeFixedInfo->SetOutputDirection(fixedOrientation);
> changeMovingInfo->SetOutputDirection(movingOrientation);
> changeFixedInfo->SetChangeDirection(true);
> changeMovingInfo->SetChangeDirection(true);
>
> changeFixedInfo->SetOutputOrigin(fixedOrigin);
> changeMovingInfo->SetOutputOrigin(movingOrigin);
> changeFixedInfo->SetChangeOrigin(true);
> changeMovingInfo->SetChangeOrigin(true);
>
> changeFixedInfo->Update();
> changeMovingInfo->Update();
>
> registration->SetFixedImage( changeFixedInfo->GetOutput() );
> registration->SetMovingImage( changeMovingInfo->GetOutput() );
> registration->SetFixedImageRegion(
> changeFixedInfo->GetOutput()->GetBufferedRegion()
> );
>
> As a check I print the origin, spacing and direction for the images and they
> are identical to when reading from DICOM files. I read the DICOM series just
> as described in chapter 7.12 in the ITK manual.
>
> Any ideas to what might be wrong here?
>
> Thanks in advance!
>
> //Patrik Brynolfsson
>
>
>
> 2009/9/16 Bill Lorensen <bill.lorensen at gmail.com>
>
>> You should get the same answers. Exactly what filters are you using in
>> each of the cases. I can't tell from your description.
>>
>> Bill
>>
>> On Tue, Sep 15, 2009 at 4:38 PM, Patrik Brynolfsson
>> <patrik.brynolfsson at gmail.com> wrote:
>> > Hello,
>> > I'm working on a Matlab mex-implementation of an ITK registration. This
>> > means that instead of loading the images from files I will load them from
>> > memory as 3D matrices. All the DICOM information needs to be passed to
>> ITK
>> > from Matlab as well, and this is where my problems begin. The information
>> > that I thought was essential for a correct registration is:
>> > * Origin
>> > * Spacing
>> > * Direction
>> > I set these using the "ChangeInformationImageFilter" and it is working;
>> the
>> > images have the origin, spacing and directional cosines that are assigned
>> to
>> > them, and they are the same as if the dicom images were loaded from
>> files.
>> > However, the result of the registration is different from that when
>> loading
>> > the images from files, with all other components set the same, (of course
>> I
>> > reinitialize the seed). Not by much, but they should be exactly the same!
>> I
>> > use Mattes mutual information metric, linear interpolation, Versor rigid
>> 3D
>> > transform and the accompanying Versor rigid 3D optimizer together with
>> the
>> > "CenteredTransformInitializer".
>> > I know that the difference in results are due to the way I read the data.
>> If
>> > I change the input from memory to reading the dicom files I get the same
>> > result as with a stand-alone .exe-application, so this is where the error
>> > is.
>> > Does anyone have any idea as to what other information I need to pass
>> from
>> > matlab to get identical results?
>> > Thanks in advance
>> > --
>> > Patrik Brynolfsson
>> >
>> >
>> >
>> > _____________________________________
>> > Powered by www.kitware.com
>> >
>> > Visit other Kitware open-source projects at
>> > http://www.kitware.com/opensource/opensource.html
>> >
>> > Please keep messages on-topic and check the ITK FAQ at:
>> > http://www.itk.org/Wiki/ITK_FAQ
>> >
>> > Follow this link to subscribe/unsubscribe:
>> > http://www.itk.org/mailman/listinfo/insight-users
>> >
>> >
>>
>
>
>
> --
> Patrik Brynolfsson
> -------------- next part --------------
> An HTML attachment was scrubbed...
> URL: <http://www.itk.org/pipermail/insight-users/attachments/20090916/f0f7f731/attachment.htm>
>
> ------------------------------
>
> _______________________________________________
> Insight-users mailing list
> Insight-users at itk.org
> http://www.itk.org/mailman/listinfo/insight-users
>
>
> End of Insight-users Digest, Vol 65, Issue 56
> *********************************************
>
More information about the Insight-users
mailing list