Res: [vtkusers] Extracting Slices using vtkImageReslice

David Gobbi dgobbi at atamai.com
Sat Apr 21 21:40:05 EDT 2007


Hi Wagner,

If you look at the third column of the sagittal matrix, you'll notice
that it is negative, while the third column of the other matrices is
positive.  That does not mean that there is anything wrong with the
matrix.  If the third column was positive, then the determinant of
that matrix would be negative which is something you always want to
avoid.

When you calculate the "point" variable in your code, you multiply it
by the ResliceAxes matrix, and then set it as the new "translation"
portion of the ResliceMatrix transform.  For the sagittal case, try
writing out the matrix multiplication on paper, elements and all, and
you should be able to see what is going on.  The "spacing" is
multiplied by the negative element in that last column, and as a
result, you are stepping away from the volume rather than into it.

Try using

  point[2] = - m_ImageData->GetSpacing()[1];

for the sagittal case.

 - David


On 4/21/07, Wagner Sales <wsalles2003 at yahoo.com.br> wrote:
>
> Hi David,
>
> Thks, really a stupid mistake from me.
> Now I'm just setting using the image origin ( in my dataset are 0,0,0 ) to
> extract the slices.
> Spacing an extraction works very fine... Except by saggital orientation! In
> this case, no one slice are extracted.
> To start, I'm using the following matrices:
>
> ------------------------------------------------------------
> double SlicerExtractor::transversalElements[16] = {
>            1, 0, 0, 0,
>            0, 1, 0, 0,
>            0, 0, 1, 0,
>            0, 0, 0, 1 };
>
> double SlicerExtractor::coronalElements[16] = {
>            1, 0, 0, 0,
>            0, 0, 1, 0,
>            0,-1, 0, 0,
>            0, 0, 0, 1 };
>
> double SlicerExtractor::sagittalElements[16] ={
>            0, 0,-1, 0,
>            1, 0, 0, 0,
>            0,-1, 0, 0,
>            0, 0, 0, 1 };
> Then, when I'll choose the orientation:
> switch(m_Axy)
>     {
>         case VolumeRepresentation::Coronal:
>             m_Matrix->DeepCopy(coronalElements);
>               break;
>         case VolumeRepresentation::Saggital:
>             m_Matrix->DeepCopy(sagittalElements);
>             break;
>         case VolumeRepresentation::Transversal:
>             m_Matrix->DeepCopy(transversalElements);
>             break;
>     }
>     m_Matrix->SetElement(0, 3, origin[0]);
>       m_Matrix->SetElement(1, 3, origin[1]);
>       m_Matrix->SetElement(2, 3, origin[2]);
> ---------------------------------------------------------
>
>
> And, finally, extract my slices:
> ---------------------------------------------------------
> m_ImageData->UpdateInformation();
>     switch(m_Axy)
>     {
>         case VolumeRepresentation::Coronal:
>             point[2] = m_ImageData->GetSpacing()[0];
>               break;
>         case VolumeRepresentation::Saggital:
>             point[2] = m_ImageData->GetSpacing()[1];
>             break;
>         case VolumeRepresentation::Transversal:
>             point[2] = m_ImageData->GetSpacing()[2];
>             break;
>     }
>     vtkMatrix4x4 *matrix = m_Reslicer->GetResliceAxes();
>     point[3] = 1.0;
>         m_Matrix->MultiplyPoint(point, temp);
>         m_Matrix->SetElement(0, 3, temp[0]);
>         m_Matrix->SetElement(1, 3, temp[1]);
>         m_Matrix->SetElement(2, 3, temp[2]);
>     m_Reslicer->SetResliceAxes(m_Matrix);
>     m_Reslicer->Update();
>
> -----------------------------------------------------------
>
> If I use the center the slices at the sagittal axis, works well ( except, of
> course, by the only center-to-end slices are extracted ).
> I'm thinking thats my sagittal matrix are wrong.
> If you or someone can help, thks in advance.
>
> Regards,
>
> Wagner Sales
>
> ----- Mensagem original ----
> De: David Gobbi <dgobbi at atamai.com>
> Para: Wagner Sales <wsalles2003 at yahoo.com.br>
> Cc: vtkusers at vtk.org
> Enviadas: Sábado, 21 de Abril de 2007 17:13:54
> Assunto: Re: Res: [vtkusers] Extracting Slices using vtkImageReslice -
> Calculating extract spacing factor
>
> Hi Wagner,
>
> I think the reason this is happening is that you are extracting the
> first slice from the centre of the volume, and then by the time you have
> incremented 255 times you have reached the edge of the volume.  The rest
> of the slices will therefore be outside of the bounds of the image.  If
> you set the slice spacing to half, then you are still only extracting
> slices from half of the volume, but with half of the original slice spacing.
>
> In your code, you should always set the last column of the ResliceAxes
> matrix to the point that lies at the centre of the slice that you want
> to extract.  Right it looks like you are incrementing that point by one
> slice-spacing each time, but you are starting at the middle of the
> volume rather than at the first slice.
>
> The "factor" should never be set to anything except for 1, unless you
> want to interpolate new slices between the original slices.
>
>  - David
>
>
>
> Wagner Sales wrote:
> > Hi David,
> >
> > I was tested then now. But the results are empty slices yet.
> > The spacing values are (just to better explain):
> > Before SetOutputSpacing:
> > - 0.9765625
> > - 0.9765625
> > - 3.4000...
> >
> > When I set OutputSpacing to minSpacing, the spacing in all axis will
> > be 0.9765625, wich are the minimum spacing in all axis.
> > These two values are to coronal and sagittal axis, and the series are
> > in transversal axis, wich have the 3.4000.... spacing.
> > After the update information, the spacing will return to the original
> > value.
> > Just to make simple, I'll use only coronal axys:
> > - The result for an GetExtent(int[6]) are 0 to 511 to this axy;
> > - Using the value 0.97... ( spacing before SetOutputSpacing and after
> > UpdateInformation ), only 255 slices will be valid.
> > - After multiplication ( 0.97... * 0.5 ) all 511 slices are writed.
> > I was used the factor 0.5 just because only  half of the 511slices are ok.
> > I'm a little bit confused if the axis spacing value are not ok, or
> > information at my dataset ( I think a good dataset, are from OsiriX )
> > are wrong, or I'm completely wrong...
> > Just because that I was asked by "how factor I can safely use?".
> > If you have some idea, thks in advance.
> >
> >
> > Regards,
> > Wagner Sales
> >
> >
> >
> >
> > ----- Mensagem original ----
> > De: David Gobbi <dgobbi at atamai.com>
> > Para: Wagner Sales <wsalles2003 at yahoo.com.br>
> > Cc: vtkusers at vtk.org
> > Enviadas: Sábado, 21 de Abril de 2007 13:10:06
> > Assunto: Re: [vtkusers] Extracting Slices using vtkImageReslice -
> > Calculating extract spacing factor
> >
> > Hi Wagner,
> >
> > I think that I know what the error in your code is.
> >
> > Since you are setting the OutputSpacing() to minSpacing.  This means
> > the spacing of the output will always be minSpacing regardless of
> > whether you are looking at the x spacing, y spacing or z spacing.
> >
> > Maybe what you should do instead is get the input spacing:
> >
> >   m_ImageData->UpdateInformation();
> >   m_ImageData->GetSpacing()[i];
> >
> > This will return the original spacing as stored in the DICOM file,
> > which will allow you to properly go through the data set
> > slice-by-slice.
> >
> > - David
> >
> >
> > On 4/21/07, Wagner Sales <wsalles2003 at yahoo.com.br> wrote:
> > >
> > > Hi all,
> > >
> > > I'm extracting slices from a DICOM dataset and writing on disk. The
> > aproach
> > > are following:
> > > - Load using vtkDICOMImageReader;
> > > - Extract slices using vtkImageReslice.
> > > - Write using vtkPNGImageWriter
> > >
> > > -- CODE:
> > >  m_Reslicer->SetInput(m_ImageData);
> > >      m_Reslicer->SetInterpolationModeToLinear();
> > >      double minSpacing = fabs(spacing[0]);
> > >      if (fabs(spacing[1]) < minSpacing)
> > >      {
> > >      minSpacing = fabs(spacing[1]);
> > >      }
> > >      if (fabs(spacing[2]) < minSpacing)
> > >      {
> > >      minSpacing = fabs(spacing[2]);
> > >      }
> > >      m_Reslicer->SetOutputSpacing(minSpacing,
> minSpacing,
> > > minSpacing);
> > >      m_Reslicer->GetOutput()->UpdateInformation();
> > >      m_Reslicer->SetOutputDimensionality(2);
> > >        m_Reslicer->SetResliceAxes(m_Matrix);
> > >      m_Reslicer->GetOutput()->UpdateInformation();
> > >      double point[4];
> > >      double temp[4];
> > >      point[0] = 0.0;
> > >      point[1] = 0.0;
> > >      switch(m_Axy)
> > >      {
> > >          case VolumeRepresentation::Coronal:
> > >              point[2] =
> > > m_Reslicer->GetOutput()->GetSpacing()[0] * 0.5 //
> what's
> > > are the factor?;
> > >                break;
> > >          case VolumeRepresentation::Saggital:
> > >              point[2] =
> > > m_Reslicer->GetOutput()->GetSpacing()[1] * 0.5 //
> what's
> > > are the factor?;
> > >              break;
> > >          case VolumeRepresentation::Transversal:
> > >              point[2] =
> > > m_Reslicer->GetOutput()->GetSpacing()[2];
> > >              break;
> > >      }
> > >      vtkMatrix4x4 *matrix = m_Reslicer->GetResliceAxes();
> > >      point[3] = 1.0;
> > >          matrix->MultiplyPoint(point, temp);
> > >          matrix->SetElement(0, 3, temp[0]);
> > >          matrix->SetElement(1, 3, temp[1]);
> > >          matrix->SetElement(2, 3, temp[2]);
> > >      m_Reslicer->Update();
> > >      m_Slice->DeepCopy(m_Reslicer->GetOutput());
> > >      return m_Slice;
> > >
> > >
> > > When I extract slices from the transversal axis, all are extracted
> > ok, but
> > > when I extract from coronal and sagittal axis, a lot of slices are
> > empty.
> > > After some playing, I was discovered that's my error are these lines:
> > >
> > >     point[2] = m_Reslicer->GetOutput()->GetSpacing()[];
> > >
> > > Then I was multiplied, when extracting coronal or sagittal slices, like
> > > showed in the code.
> > > Works fine! No empty slices anymore! But, the factor are fixed, and
> > I don't
> > > know how to calculate the factor 0.5.
> > > I'm completely wrong? Someone can help?
> > >
> > > Thks in advance,
> > >
> > > Wagner Sales
> > >
> > >
> > >
> > >
> > >
> > > __________________________________________________
> > > Fale com seus amigos de graça com o novo Yahoo! Messenger
> > > http://br.messenger.yahoo.com/
> > > _______________________________________________
> > > This is the private VTK discussion list.
> > > Please keep messages on-topic. Check the FAQ at:
> > > http://www.vtk.org/Wiki/VTK_FAQ
> > > Follow this link to subscribe/unsubscribe:
> > > http://www.vtk.org/mailman/listinfo/vtkusers
> > >
> > >
> > _______________________________________________
> > This is the private VTK discussion list.
> > Please keep messages on-topic. Check the FAQ at:
> > http://www.vtk.org/Wiki/VTK_FAQ
> > Follow this link to subscribe/unsubscribe:
> > http://www.vtk.org/mailman/listinfo/vtkusers
> >
> >
> > __________________________________________________
> > Fale com seus amigos de graça com o novo Yahoo! Messenger
> > http://br.messenger.yahoo.com/
>
> _______________________________________________
> This is the private VTK discussion list.
> Please keep messages on-topic. Check the FAQ at:
> http://www.vtk.org/Wiki/VTK_FAQ
> Follow this link to subscribe/unsubscribe:
> http://www.vtk.org/mailman/listinfo/vtkusers
>
>
> __________________________________________________
> Fale com seus amigos de graça com o novo Yahoo! Messenger
> http://br.messenger.yahoo.com/



More information about the vtkusers mailing list