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