[Insight-users] Inhomogenity correction using MRIBiasField Correction

cspl affable at hd2 . dot . net . in
Sat, 3 Aug 2002 16:49:25 +0530


This is a multi-part message in MIME format.

------=_NextPart_000_002E_01C23B0D.B9CC1AD0
Content-Type: text/plain;
	charset="iso-8859-1"
Content-Transfer-Encoding: quoted-printable

Dear Mr.bjorn & Friends,

I am working on itkMRIBiasFieldCorrectionFilter.I got hte problem with =
array to assign mean and sigma values. I got exception while running at =
setting parameters like as follows.

Volume* CIPServerDoc::Inhomogenity(Volume *vol,Volume *mask)
{
 if (vol =3D=3D NULL) return vol;
 short *Buffer =3D new short[vol->width * vol->height *  vol->depth];
 short *Buffer1 =3D new short[mask->width * mask->height *  =
mask->depth];
=20
 Buffer =3D (short *) vol->Mem;
 Buffer1=3D (short *) mask->Mem;=20
=20
 typedef itk::Image<float,3> ImageType;
 typedef itk::Image<float,3> ImageType1;
 typedef itk::Image<float,3> ImageType2;

 itk::RawImageIO<float, 3>::Pointer output_io
      =3D itk::RawImageIO<float, 3>::New();
=20
     output_io->SetByteOrderToLittleEndian();
     output_io->SetFileTypeToBinary();
     output_io->SetFileDimensionality(3);

     // Initialize the initial and target volumes.
        ImageType::Pointer vol_init =3D ImageType::New();
 ImageType1::Pointer vol_init1 =3D ImageType1::New();
 ImageType2::Pointer vol_init2 =3D ImageType2::New();
=20
   =20
    ImageType::SizeType sz;
    sz[0] =3D vol->width;
    sz[1] =3D vol->height;
 sz[2]=3Dvol->depth;=20

 ImageType1::SizeType sz1;
 sz1[0] =3D mask->width;
 sz1[1] =3D mask->height;
 sz1[2]=3D  mask->depth;=20

   =20
 ImageType::IndexType idx;
 idx[0] =3D 0;
 idx[1] =3D 0;
 idx[2] =3D 0 ;

 ImageType1::IndexType idx1;
        idx1[0] =3D 0;
        idx1[1] =3D 0;
 idx1[2] =3D 0 ;
   =20
        ImageType::RegionType reg; =20
        reg.SetSize(sz);
        reg.SetIndex(idx);

 ImageType1::RegionType reg1; =20
        reg.SetSize(sz1);
        reg.SetIndex(idx1);
   =20
        vol_init->SetRegions(reg);
 vol_init1->SetRegions(reg1);
 vol_init2->SetRegions(reg);
        vol_init->Allocate();
 vol_init1->Allocate();
 vol_init2->Allocate();


=20
 for(int z=3D0; z<vol->depth; ++z )
 {=20
  idx[2]=3Dz;
=20
         for (int y =3D 0; y < vol->height; ++y)
          {
           idx[1] =3D y;
           for (int x =3D 0; x < vol->width; ++x)
            {
              idx[0] =3D x;
      vol_init->SetPixel(idx, *(Buffer   +(z*vol->height*vol->width)+y =
*vol->width +x));
         }
     }
 =20
 }
=20
=20
 ////mask
 for(int z1=3D0; z1<mask->depth; ++z1 )
 {=20
     idx1[2]=3Dz1;
    for (int y1 =3D 0; y1<mask->height; ++y1)
        {
          idx1[1] =3D y1;
    //AfxMessageBox("in mask height loop");
          for (int x1=3D0;  x1 < mask->width; ++x1)
            {
              idx1[0] =3D x1;
       vol_init2->SetPixel(idx1,*(Buffer1   =
+(z1*mask->height*mask->width)+y1 *mask->width +x1));
      }
 =20
         }
  =20
 }
   =20
typedef BiasFieldCorrectionFilter<ImageType,ImageType,ImageType1>=20
 Corrector;
 Corrector::Pointer filter =3D Corrector::New() ;
//parameters initiliazation
  bool useLog =3D true;
  int degree =3D 3;
  int sliceDirection =3D 2;
 =20

  vnl_vector<double> coefficientVector ;
  itk::Array<double> *classMeans =3D new itk::Array<double>;
 =20
 classMeans[0]=3D1500;
 classMeans[1]=3D570;

 itk::Array<double> *classSigmas=3D new itk::Array<double> ;
 classSigmas[0]=3D100;
 classSigmas[1]=3D70;

 int volumeMaximumIteration =3D 2000;=20
 int interSliceMaximumIteration =3D 20;=20
// double initialRadius ;
 double growth =3D 1.05;
 double shrink =3D 0.0;
 AfxMessageBox("before  set input");
 filter->SetInput(vol_init);
 filter->SetInputMask(vol_init2);


 filter->IsBiasFieldMultiplicative(useLog) ;
///////////////////////////////////////////////////////////////////////
 filter->SetTissueClassStatistics(*classMeans,*classSigmas) ;/// /* here =
Getting exception  */
//////////////////////////////////////////////////////////////////////

 filter->SetOptimizerGrowthFactor(growth) ;

 filter->SetOptimizerShrinkFactor(shrink) ;

 filter->SetVolumeCorrectionMaximumIteration(volumeMaximumIteration) ;
=20
 =
//filter->SetInterSliceCorrectionMaximumIteration(interSliceMaximumIterat=
ion) ;
 //filter->SetOptimizerInitialRadius(initialRadius) ;
 filter->SetBiasFieldDegree(degree) ;
=20
// filter->SetUsingSlabIdentification(usingSlabIdentification) ;
 filter->SetSlicingDirection(sliceDirection) ;
=20
 =20
 AfxMessageBox("before update");
 filter->Update();
=20
 vol_init2=3Dfilter->GetOutput();
=20
 for(unsigned int z=3D0; z<vol->depth; ++z )
 {
  idx[2]=3Dz;
=20
      for (unsigned int y =3D 0; y < vol->height; ++y)
        {
          idx[1] =3D y;
          for (unsigned int x =3D 0; x < vol->width; ++x)
            {
              idx[0] =3D x;
              *(Buffer +(z*vol->height*vol->width)+y *vol->width    =
+x)=3Dvol_init->GetPixel(idx);
            }
        }
 }

 vol->Mem =3D Buffer;
        return vol;

}


Why I am getting exception when call  =
filter->SetTissueClassStatistics(*classMeans,*classSigmas) .

Please give me suggestion.

Thanking you,
Regards,
Ramakrishna

------=_NextPart_000_002E_01C23B0D.B9CC1AD0
Content-Type: text/html;
	charset="iso-8859-1"
Content-Transfer-Encoding: quoted-printable

<!DOCTYPE HTML PUBLIC "-//W3C//DTD HTML 4.0 Transitional//EN">
<HTML><HEAD>
<META http-equiv=3DContent-Type content=3D"text/html; =
charset=3Diso-8859-1">
<META content=3D"MSHTML 5.50.4134.600" name=3DGENERATOR>
<STYLE></STYLE>
</HEAD>
<BODY bgColor=3D#ffffff>
<DIV><FONT face=3DArial size=3D2>Dear Mr.bjorn &amp; =
Friends,</FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>I am working on =
itkMRIBiasFieldCorrectionFilter.I=20
got hte problem with array to assign mean and sigma values. I got =
exception=20
while running at setting parameters like as follows.</FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>Volume* =
CIPServerDoc::Inhomogenity(Volume=20
*vol,Volume *mask)<BR>{<BR>&nbsp;if (vol =3D=3D NULL) return =
vol;<BR>&nbsp;short=20
*Buffer =3D new short[vol-&gt;width * vol-&gt;height *=20
&nbsp;vol-&gt;depth];<BR>&nbsp;short *Buffer1 =3D new =
short[mask-&gt;width *=20
mask-&gt;height * &nbsp;mask-&gt;depth];<BR>&nbsp;<BR>&nbsp;Buffer =3D =
(short *)=20
vol-&gt;Mem;<BR>&nbsp;Buffer1=3D (short *) mask-&gt;Mem;=20
<BR>&nbsp;<BR>&nbsp;typedef itk::Image&lt;float,3&gt;=20
ImageType;<BR>&nbsp;typedef itk::Image&lt;float,3&gt;=20
ImageType1;<BR>&nbsp;typedef itk::Image&lt;float,3&gt; =
ImageType2;</FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>&nbsp;itk::RawImageIO&lt;float, =
3&gt;::Pointer=20
output_io<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; =3D =
itk::RawImageIO&lt;float,=20
3&gt;::New();<BR>&nbsp;<BR>&nbsp;&nbsp;&nbsp;&nbsp;=20
output_io-&gt;SetByteOrderToLittleEndian();<BR>&nbsp;&nbsp;&nbsp;&nbsp;=20
output_io-&gt;SetFileTypeToBinary();<BR>&nbsp;&nbsp;&nbsp;&nbsp;=20
output_io-&gt;SetFileDimensionality(3);</FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>&nbsp;&nbsp;&nbsp;&nbsp; // Initialize =
the initial=20
and target volumes.<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
ImageType::Pointer vol_init =3D =
ImageType::New();<BR>&nbsp;ImageType1::Pointer=20
vol_init1 =3D ImageType1::New();<BR>&nbsp;ImageType2::Pointer vol_init2 =
=3D=20
ImageType2::New();<BR>&nbsp;<BR>&nbsp;&nbsp;&nbsp; =
<BR>&nbsp;&nbsp;&nbsp;=20
ImageType::SizeType sz;<BR>&nbsp;&nbsp;&nbsp; sz[0] =3D=20
vol-&gt;width;<BR>&nbsp;&nbsp;&nbsp; sz[1] =3D=20
vol-&gt;height;<BR>&nbsp;sz[2]=3Dvol-&gt;depth; </FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>&nbsp;ImageType1::SizeType =
sz1;<BR>&nbsp;sz1[0] =3D=20
mask-&gt;width;<BR>&nbsp;sz1[1] =3D =
mask-&gt;height;<BR>&nbsp;sz1[2]=3D&nbsp;=20
mask-&gt;depth; </FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>&nbsp;&nbsp;&nbsp; =
<BR>&nbsp;ImageType::IndexType=20
idx;<BR>&nbsp;idx[0] =3D 0;<BR>&nbsp;idx[1] =3D 0;<BR>&nbsp;idx[2] =3D 0 =

;</FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>&nbsp;ImageType1::IndexType=20
idx1;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; idx1[0] =3D=20
0;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; idx1[1] =3D =
0;<BR>&nbsp;idx1[2] =3D=20
0 ;<BR>&nbsp;&nbsp;&nbsp; <BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; =

ImageType::RegionType reg;&nbsp; =
<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
reg.SetSize(sz);<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
reg.SetIndex(idx);</FONT></DIV>
<DIV>&nbsp;</DIV>
<DIV><FONT face=3DArial size=3D2>&nbsp;ImageType1::RegionType =
reg1;&nbsp;=20
<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
reg.SetSize(sz1);<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
reg.SetIndex(idx1);<BR>&nbsp;&nbsp;&nbsp;=20
<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
vol_init-&gt;SetRegions(reg);<BR>&nbsp;vol_init1-&gt;SetRegions(reg1);<BR=
>&nbsp;vol_init2-&gt;SetRegions(reg);<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&n=
bsp;&nbsp;=20
vol_init-&gt;Allocate();<BR>&nbsp;vol_init1-&gt;Allocate();<BR>&nbsp;vol_=
init2-&gt;Allocate();</FONT></DIV>
<DIV>&nbsp;</DIV><FONT face=3DArial size=3D2>
<DIV><BR>&nbsp;<BR>&nbsp;for(int z=3D0; z&lt;vol-&gt;depth; ++z =
)<BR>&nbsp;{=20
<BR>&nbsp;&nbsp;idx[2]=3Dz;<BR>&nbsp;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&n=
bsp;&nbsp;&nbsp;=20
for (int y =3D 0; y &lt; vol-&gt;height;=20
++y)<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
{<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; idx[1] =
=3D=20
y;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; for =
(int x =3D=20
0; x &lt; vol-&gt;width;=20
++x)<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp=
;=20
{<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&n=
bsp;&nbsp;=20
idx[0] =3D x;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; =
vol_init-&gt;SetPixel(idx,=20
*(Buffer &nbsp;&nbsp;+(z*vol-&gt;height*vol-&gt;width)+y *vol-&gt;width=20
+x));<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
}<BR>&nbsp;&nbsp;&nbsp;&nbsp; }<BR>&nbsp;=20
<BR>&nbsp;}<BR>&nbsp;<BR>&nbsp;<BR>&nbsp;////mask<BR>&nbsp;for(int =
z1=3D0;=20
z1&lt;mask-&gt;depth; ++z1 )<BR>&nbsp;{ <BR>&nbsp;&nbsp;&nbsp;=20
&nbsp;idx1[2]=3Dz1;<BR>&nbsp;&nbsp;&nbsp; for (int y1 =3D 0; =
y1&lt;mask-&gt;height;=20
++y1)<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
{<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; idx1[1] =3D=20
y1;<BR>&nbsp;&nbsp;&nbsp; //AfxMessageBox("in mask height=20
loop");<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; for =
(int=20
x1=3D0;&nbsp; x1 &lt; mask-&gt;width;=20
++x1)<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbs=
p;=20
{<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&n=
bsp;&nbsp;=20
idx1[0] =3D x1;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
vol_init2-&gt;SetPixel(idx1,*(Buffer1=20
&nbsp;&nbsp;+(z1*mask-&gt;height*mask-&gt;width)+y1 *mask-&gt;width=20
+x1));<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<BR>&nbsp;=20
<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<BR>&nbsp;&nbsp;=20
<BR>&nbsp;}<BR>&nbsp;&nbsp;&nbsp; <BR>typedef=20
BiasFieldCorrectionFilter&lt;ImageType,ImageType,ImageType1&gt;=20
<BR>&nbsp;Corrector;<BR>&nbsp;Corrector::Pointer filter =3D =
Corrector::New()=20
;<BR>//parameters initiliazation<BR>&nbsp; bool useLog =3D =
true;<BR>&nbsp; int=20
degree =3D 3;<BR>&nbsp; int sliceDirection =3D 2;<BR>&nbsp; </DIV>
<DIV>&nbsp;</DIV>
<DIV>&nbsp; vnl_vector&lt;double&gt; coefficientVector ;<BR>&nbsp;=20
itk::Array&lt;double&gt; *classMeans =3D new =
itk::Array&lt;double&gt;;<BR>&nbsp;=20
<BR>&nbsp;classMeans[0]=3D1500;<BR>&nbsp;classMeans[1]=3D570;</DIV>
<DIV>&nbsp;</DIV>
<DIV>&nbsp;itk::Array&lt;double&gt; *classSigmas=3D new =
itk::Array&lt;double&gt;=20
;<BR>&nbsp;classSigmas[0]=3D100;<BR>&nbsp;classSigmas[1]=3D70;</DIV>
<DIV>&nbsp;</DIV>
<DIV>&nbsp;int volumeMaximumIteration =3D 2000; <BR>&nbsp;int=20
interSliceMaximumIteration =3D 20; <BR>//&nbsp;double initialRadius=20
;<BR>&nbsp;double growth =3D 1.05;<BR>&nbsp;double shrink =3D=20
0.0;<BR>&nbsp;AfxMessageBox("before&nbsp; set=20
input");<BR>&nbsp;filter-&gt;SetInput(vol_init);<BR>&nbsp;filter-&gt;SetI=
nputMask(vol_init2);</DIV>
<DIV>&nbsp;</DIV>
<DIV><BR>&nbsp;filter-&gt;IsBiasFieldMultiplicative(useLog)=20
;<BR>////////////////////////////////////////////////////////////////////=
///<BR>&nbsp;filter-&gt;SetTissueClassStatistics(*classMeans,*classSigmas=
)=20
;/// /*&nbsp;here Getting exception&nbsp;=20
*/<BR>///////////////////////////////////////////////////////////////////=
///</DIV>
<DIV>&nbsp;</DIV>
<DIV>&nbsp;filter-&gt;SetOptimizerGrowthFactor(growth) ;</DIV>
<DIV>&nbsp;</DIV>
<DIV>&nbsp;filter-&gt;SetOptimizerShrinkFactor(shrink) ;</DIV>
<DIV>&nbsp;</DIV>
<DIV>&nbsp;filter-&gt;SetVolumeCorrectionMaximumIteration(volumeMaximumIt=
eration)=20
;<BR>&nbsp;<BR>&nbsp;//filter-&gt;SetInterSliceCorrectionMaximumIteration=
(interSliceMaximumIteration)=20
;<BR>&nbsp;//filter-&gt;SetOptimizerInitialRadius(initialRadius)=20
;<BR>&nbsp;filter-&gt;SetBiasFieldDegree(degree)=20
;<BR>&nbsp;<BR>//&nbsp;filter-&gt;SetUsingSlabIdentification(usingSlabIde=
ntification)=20
;<BR>&nbsp;filter-&gt;SetSlicingDirection(sliceDirection) =
;<BR>&nbsp;<BR>&nbsp;=20
<BR>&nbsp;AfxMessageBox("before=20
update");<BR>&nbsp;filter-&gt;Update();<BR>&nbsp;<BR>&nbsp;vol_init2=3Dfi=
lter-&gt;GetOutput();<BR>&nbsp;<BR>&nbsp;for(unsigned=20
int z=3D0; z&lt;vol-&gt;depth; ++z=20
)<BR>&nbsp;{<BR>&nbsp;&nbsp;idx[2]=3Dz;<BR>&nbsp;<BR>&nbsp;&nbsp;&nbsp;&n=
bsp;&nbsp;=20
for (unsigned int y =3D 0; y &lt; vol-&gt;height;=20
++y)<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
{<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; idx[1] =3D=20
y;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; for =
(unsigned int x=20
=3D 0; x &lt; vol-&gt;width;=20
++x)<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp=
;=20
{<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&n=
bsp;&nbsp;=20
idx[0] =3D=20
x;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&=
nbsp;&nbsp;=20
*(Buffer +(z*vol-&gt;height*vol-&gt;width)+y *vol-&gt;width=20
&nbsp;&nbsp;&nbsp;+x)=3Dvol_init-&gt;GetPixel(idx);<BR>&nbsp;&nbsp;&nbsp;=
&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
}<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp; }<BR>&nbsp;}</DIV>
<DIV>&nbsp;</DIV>
<DIV>&nbsp;vol-&gt;Mem =3D =
Buffer;<BR>&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;&nbsp;=20
return vol;</DIV>
<DIV>&nbsp;</DIV>
<DIV>}</DIV>
<DIV>&nbsp;</DIV>
<DIV><BR>Why I am getting exception when call&nbsp;=20
filter-&gt;SetTissueClassStatistics(*classMeans,*classSigmas) .</DIV>
<DIV>&nbsp;</DIV>
<DIV>Please give me suggestion.</DIV>
<DIV>&nbsp;</DIV>
<DIV>Thanking you,<BR>Regards,<BR>Ramakrishna</FONT></DIV></BODY></HTML>

------=_NextPart_000_002E_01C23B0D.B9CC1AD0--