[Insight-users] index problem
Lydia Ng
lng at insightful.com
Mon, 12 Jan 2004 10:29:14 -0800
Hi Corinne,
Your code looks a lot like the functionality of the
MultiResolutionPDEDeformableRegistration. Over the few days I made some
modifications to this class and VectorResampleFilter so that it handles
arbitrary image sizes and starting indices (i.e. does not have to be =
multiple
of the shrink factors).
Thought the revised code might be of interest to you.
Cheers,
Lydia
> -----Original Message-----
> From: Corinne Mattmann [mailto:mattmaco at ee.ethz.ch]
> Sent: Thursday, January 08, 2004 2:02 PM
> To: Lydia Ng; 'Luis Ibanez'
> Cc: insight-users at itk.org
> Subject: RE: [Insight-users] index problem
>=20
> I just got an email from Luis, so I don't know if you still need the
> code. But as I have already prepared the email, I send it anyway.
> Concerning the truncation: You have the same problem with the size of
> the image, but you can solve this problem by setting the expand =
factors
> correctly (a floating point number). Shouldn't this also be possible
> with the starting index?
>=20
>=20
> Hi Lydia,
>=20
> Here is the relevant code of how I use the two filters.
> Below the code you can see part of the couts representing what's
> happening if I input an image with starting index [3,1,1].
>=20
> Hope that nothing is missing,
> Corinne
>=20
> =
------------------------------------------------------------------------
> --------
>=20
> typedef
> =
itk::RecursiveMultiResolutionPyramidImageFilter<InputImageType,InternalI
> mageType> PyramidType;
> typedef
> =
itk::DemonsRegistrationFilter<InternalImageType,InternalImageType,Output
> ImageType> RegistrationFilterType;
> typedef
> itk::VectorExpandImageFilter<OutputImageType,OutputImageType>
> ExpanderType;
> typedef typename ExpanderType::ExpandFactorsType =
ExpandFactorsType;
>=20
> std::vector<RegistrationFilterType::Pointer> filter;
> filter.resize(m_NumberOfLevels);
> std::vector<ExpanderType::Pointer> FieldExpander;
> FieldExpander.resize(m_NumberOfLevels);
>=20
> for (int CurrentLevel=3D0; CurrentLevel<m_NumberOfLevels;
> CurrentLevel++)
> {
> std::cout << " beginning level " << CurrentLevel << =
std::endl;
>=20
> MovingPyramid->GetOutput(CurrentLevel)->Update();
> FixedPyramid->GetOutput(CurrentLevel)->Update();
>=20
> ... (some couts)
>=20
> filter[CurrentLevel] =3D RegistrationFilterType::New();
>=20
> =
filter[CurrentLevel]->SetFixedImage(FixedPyramid->GetOutput(CurrentLevel
> ));
>=20
> =
filter[CurrentLevel]->SetMovingImage(MovingPyramid->GetOutput(CurrentLev
> el));
>=20
> =
filter[CurrentLevel]->SetNumberOfIterations(m_NumberOfIterations[Current
> Level]);
>=20
> =
filter[CurrentLevel]->SetStandardDeviations(m_StandardDeviations[Current
> Level]);
> filter[CurrentLevel]->SetReleaseDataFlag(true);
> if(CurrentLevel>0){
>=20
> =
filter[CurrentLevel]->SetInitialDeformationField(FieldExpander[CurrentLe
> vel-1]->GetOutput());
> }
> filter[CurrentLevel]->Update();
>=20
> .... (calculate the expand factors)
>=20
> // expand the deformation field for the next level
> FieldExpander[CurrentLevel] =3D ExpanderType::New();
>=20
> =
FieldExpander[CurrentLevel]->SetInput(filter[CurrentLevel]->GetOutput())
> ;
> FieldExpander[CurrentLevel]->SetExpandFactors(expandFactors);
> FieldExpander[CurrentLevel]->UpdateLargestPossibleRegion();
> FieldExpander[CurrentLevel]->SetReleaseDataFlag(true);
>=20
>=20
> // at the last level, the pipeline has to be connected
> if(CurrentLevel =3D=3D (m_NumberOfLevels-1)){
>=20
> =
this->GetOutput()->SetRequestedRegion(FieldExpander[CurrentLevel]->GetOu
> tput()->GetRequestedRegion());
> FieldExpander[CurrentLevel]->GraftOutput(outputPtr);
> FieldExpander[CurrentLevel]->Update();
> this->GraftOutput(FieldExpander[CurrentLevel]->GetOutput());
> }
>=20
> std::cout << " end level " << CurrentLevel << std::endl;
>=20
> }
>=20
> =
------------------------------------------------------------------------
> --------
>=20
> Image region of the input images:
> Dimension: 3
> Index: [3, 1, 1]
> Size: [47, 49, 49]
>=20
> beginning level 0
> Image region of the pyramid output:
> Dimension: 3
> Index: [2, 1, 1]
> Size: [23, 24, 24]
>=20
> end level 0
> beginning level 1
> Image region of the pyramid output:
> Dimension: 3
> Index: [3, 1, 1]
> Size: [47, 49, 49]
>=20
> Image region of the expanded deformation field:
> Dimension: 3
> Index: [4, 2, 2]
> Size: [47, 49, 49]
> =3D> Program crashes
>=20
>=20
>=20
>=20
> -----Original Message-----
> From: Lydia Ng [mailto:lng at insightful.com]
> Sent: Thursday, January 08, 2004 12:23 PM
> To: Corinne Mattmann; Luis Ibanez
> Cc: insight-users at itk.org
> Subject: RE: [Insight-users] index problem
>=20
>=20
> Hi Corinne,
>=20
> Would you mind posting some code that illustrates this problem? I.e. =
how
> the various filter is being used. I am curious if the problem lies =
with
> the indexing or the interaction between the filters or if some
> fundamental assumption were made in the algorithm that assume zero
> starting index.
>=20
> Thanks,
> Lydia
>=20
>=20
> > -----Original Message-----
> > From: Corinne Mattmann [mailto:mattmaco at ee.ethz.ch]
> > Sent: Tuesday, January 06, 2004 4:26 PM
> > To: 'Luis Ibanez'
> > Cc: insight-users at itk.org
> > Subject: RE: [Insight-users] index problem
> >
> > Hi Luis,
> >
> > Thanks for your answer. I would like to make my image smaller with =
the
>=20
> > same content (cut some pixels away at the border and set another
> > index), so I think A) should work. At the moment I have another
> > questions: I am programming a filter which registers two images with
> > the demons algorithm using a multi-resolution pyramid. To set up the
> > pyramid I'm using itkRecursiveMultiResolutionPyramidImageFilter.
> > If I have an image with index [3,1,1] and if I reduce it by a factor
> of
> > 2, the output of the pyramid would be an image with index [2,1,1].
> > This images I input into the demons-filter and get as output a
> > deformable field with the same index. So far so good.
> > But if I now expand this deformable field with the
> > itkVectorExpandImageFilter, I get an image with index [4,2,2], which
> > results in an error when I apply it to the images with index =
[3,1,1].
> > Took me quite long to find this problem ;-)
> > I wanted to set the "right" index at the expand filter but I =
couldn't
> > find an appropriate function. What solution do you suggest?
> >
> > Thanks very much for your answer and have a nice evening, Corinne
> >
> >
> > -----Original Message-----
> > From: Luis Ibanez [mailto:luis.ibanez at kitware.com]
> > Sent: Tuesday, January 06, 2004 4:48 PM
> > To: Corinne Mattmann
> > Cc: insight-users at itk.org
> > Subject: Re: [Insight-users] resize image
> >
> >
> >
> > Hi Corinne
> >
> > What is the goal that you pursue by resizing the image ?
> >
> >
> > A) Do you want to have as a result the same image
> > with an additional band of pixels around ?
> >
> > or
> >
> > B) Do you simply want to get a larger size and
> > don't care about discarding the original
> > content of the image ?
> >
> > If you are in case (A) you may want to use the PasteImageFilter
> > =
http://www.itk.org/Insight/Doxygen/html/classitk_1_1PasteImageFilter.h
> > tm
> > l
> >
> >
> > If you are in case (B) the best way to go is simply
> > to create a new image. You could recover the origin
> > and spacing information from the original image...
> > but you have to think twice if the origin should be
> > conserved or should be shifted...
> >
> >
> > The code will look like:
> >
> > // get the original image
> > ImageType::ConstPointer image1 =3D GetImageSomehow();
> >
> > // prepare new size
> > ImageType::SizeType size;
> > ImageType::IndexType start;
> > ImageType::RegionType region;
> >
> > size[0] =3D 512; // new image size
> > size[1] =3D 512;
> > size[2] =3D 512;
> >
> > start[0] =3D 0; // or maybe start somewhere else... ?
> > start[1] =3D 0;
> > start[2] =3D 0;
> >
> > region.SetSize( size );
> > region.SetIndex( start );
> >
> > ImageType::Pointer image2 =3D ImageType::New();
> >
> > // copy origin and spacing
> > image2->CopyInformation( image1 );
> >
> > image2->SetRegions( region );
> > image2->Allocate();
> >
> >
> >
> > The error message that you are getting seem
> > to indicate that you are trying to call SetRegions()
> > in a ConstPointer. Const-ness will not allow you
> > to modify the regions of an image. Note also that
> > you don't want to do this on the output of a filter
> > since the image is owned and controled by the filter.
> >
> > Any subsequent call to Update() on the filter will
> > reset the content of the image and override your modifications.
> >
> >
> >
> > Please let us know if you have further questions,
> >
> >
> > Thanks
> >
> >
> > Luis
> >
> >
> >
> > --------------------------
> > Corinne Mattmann wrote:
> >
> > > Hi,
> > >
> > > I would like to resize my image (index and size). What's the =
easiest
>=20
> > > way to do this? I tried something like this:
> > >
> > > typedef ImageRegion<ImageDimension> ImageRegionType;
> > > ImageRegionType region;
> > > typedef Index<ImageDimension> IndexType;
> > > IndexType index;
> > > typedef Size<ImageDimension> SizeType;
> > > SizeType size;
> > > const long in[3] =3D {0,0,0};
> > > index.SetIndex(in);
> > > const unsigned long si[3] =3D {50,50,50};
> > > size.SetSize(si);
> > > region.SetIndex(index);
> > > region.SetSize(size);
> > > image->SetRegions(region);
> > >
> > > But I got the following error:
> > >
> > > error C2663: 'SetRegions' : 2 overloads have no legal conversion
> > for
> > > 'this' pointer
> > >
> > >
> > > Thanks very much,
> > > Corinne Mattmann
> > >
> > > _______________________________________________
> > > Insight-users mailing list
> > > Insight-users at itk.org
> > > http://www.itk.org/mailman/listinfo/insight-users
> > >
> >
> >
> >
> > _______________________________________________
> > Insight-users mailing list
> > Insight-users at itk.org
> > http://www.itk.org/mailman/listinfo/insight-users