The output of the Joseph filter (of any forward projector) is a stack of projection images, i.e., a 3D image. If you slice this 3D image in the z direction, you get one projection images. If you slice this 3D image in the y direction, you get a sinogram, e.g., the 2D sinogram of the central slice (plane of the source trajectory) if y=0. The best way to have this sinogram is to set size[1]=1 and there is no better way to do it in RTK. For zero-padding, you need to set extendRegion[?]=1, with (for all dimensions, i.e., ?=0,1,2) and SetConstant(0). I'm surprised though that your image is a black square, something else must be going on. One thing to remember, the typical RTK rotation is around the y axis so your 3D slice must be 1 in the y direction. A simple way to change the orientation of one slice is to change its direction with ChangeImageFilter Simon
On Sat, May 14, 2016 at 11:48 AM, Tobias Stein <[email protected]> wrote: > Hello Simon, > > > > Thanks for your advice so far. I’m getting to know this filter better. At > least I can now use it and save the result into a file. > > > > I updated my code which now uses a volume. > > The values of the input and output volume are equal except of the > z-dimension of the output which I changed to the number of projection > images. > > const unsigned int Dimension = 3; > > …Reading image… > > ImageType::Pointer inputImage = reader->GetOutput(); > > ImageType::Pointer outputImage = ImageType::New(); > > > > RegionType inputRegion = inputImage->GetLargestPossibleRegion(); > > > > const unsigned int NumberOfProjectionImages = 45; > > > > ImageType::IndexType start; > > start[0] = inputRegion.GetIndex()[0]; > > start[1] = inputRegion.GetIndex()[1]; > > start[2] = inputRegion.GetIndex()[2]; > > > > RegionType::SizeType size; > > size[0] = inputRegion.GetSize()[0]; > > size[1] = inputRegion.GetSize()[1]; > > size[2] = NumberOfProjectionImages; > > > > RegionType region; > > region.SetSize(size); > > region.SetIndex(start); > > > > ImageType::SpacingType spacing; > > spacing[0] = inputImage->GetSpacing()[0]; > > spacing[1] = inputImage->GetSpacing()[1]; > > spacing[2] = inputImage->GetSpacing()[2]; > > > > ImageType::PointType origin; > > origin[0] = inputImage->GetOrigin()[0]; > > origin[1] = inputImage->GetOrigin()[1]; > > origin[2] = inputImage->GetOrigin()[2]; > > > > outputImage->SetRegions(region); > > outputImage->SetSpacing(spacing); > > outputImage->SetOrigin(origin); > > outputImage->Allocate(); > > outputImage->FillBuffer(size[0] * size[1] * size[2]); > > > > // Geometry > > typedef rtk::ThreeDCircularProjectionGeometry GeometryType; > > GeometryType::Pointer geometry = GeometryType::New(); > > for (unsigned int i = 0; i < NumberOfProjectionImages; i++) > > geometry->AddProjection(510., 510., i*8.); > > > > typedef rtk::JosephForwardProjectionImageFilter<ImageType, > ImageType> ForwardType; > > ForwardType::Pointer forward = ForwardType::New(); > > forward->SetInput(1, inputImage); > > forward->SetInput(0, outputImage); > > forward->SetGeometry(geometry); > > forward->Update(); > > > > I’m not sure how I have to add projection for the geometry to get a result > which I’m familiar with. > > > > Can you describe how do I have to interpret the result of the > Joseph-Filter? I am familiar with sinograms and how they are computed, but > one sinogram is created by the projections of one image slice as far as I > understood it. What I want is to manipulate the image in the radon space > and then perform a filtered back projection. Typically the articles > describe techniques to reduce metal artifacts by using one 2D Slice. Is it > possible to achieve something like that with the JosephFilter? I’m guessing > since this filter won’t work well for 2D I have to use another filter. Can > you give me a hint if I can do something like that with RTK or another > ITK-based framework? Here <https://i.imgur.com/1yBz3o2.png> a diagram > which shows an example what I want to try. > > > > I tried also to use my code with one slice and added a zero padding like > that: > > typedef itk::ConstantPadImageFilter <ImageType, ImageType> > ConstantPadImageFilterType; > > ConstantPadImageFilterType::Pointer padFilter = > ConstantPadImageFilterType::New(); > > padFilter->SetInput(inputImage); > > > > ImageType::SizeType extendRegion; > > extendRegion[0] = 0; > > extendRegion[1] = 0; > > extendRegion[2] = 0; > > > > padFilter->SetPadLowerBound(extendRegion); > > padFilter->SetPadUpperBound(extendRegion); > > padFilter->SetConstant(1); > > padFilter->Update(); > > But the result of my png is a black square. > > Tobias > > > > > > *Von:* [email protected] [mailto:[email protected]] *Im Auftrag von *Simon > Rit > *Gesendet:* Freitag, 13. Mai 2016 18:15 > > *An:* Tobias Stein <[email protected]> > *Cc:* [email protected] > *Betreff:* Re: [Rtk-users] forward and back projection - MITK > > > > Dear Tobias, > > If you intialize the data of your itk::Image with FillBuffer with the > corresponding constant value, it is quite similar. The constant image value > allow us to stream the computation whereas you have to allocate the whole > image if you pass directly itk::Image. > > We only work in 3D and probably our projector don't work with 2D, that is > probably the problem here. You would have to correct the code for 2D but > that's a lot of work, we generally prefer pseudo 2D, see next point. > > I doubt that our forward projector will work well with one slice only if > you use 3D, we assume that the volume starts at the first pixel and ends at > the last pixel. You should probably add a zero border around with > itk::PadImage to see something. > > To avoid aliasing, you can cut frequencies in the ramp filter (look for > CutFrequency in the doxygen of rtk::FFTRampImageFilter). > > Simon > > > > On Fri, May 13, 2016 at 5:35 PM, Tobias Stein <[email protected]> > wrote: > > Thanks for the information. > > > > To keep it simple I want to project forward a single slice which I load as > a DICOM file. > > Now I want to compute the sinogram of that slice and write it in a new > file. > > I got a little confused about the ConstantImageSource which is used in the > test class. Do I need it to compute a sinogram or is just a replacement of > a itk::Image? > > > > Here is my code so far: > > > > const unsigned int Dimension = 2; > > typedef float PixelType; > > typedef itk::Image< PixelType, Dimension > ImageType; > > … > > (reading image by itk::ImageFileReader which works) > > … > > ImageType::Pointer inputImage = reader->GetOutput(); > > > > const unsigned int NumberOfProjectionImages = 1; > > > > ImageType::Pointer outputImage = ImageType::New(); > > > > typedef rtk::JosephForwardProjectionImageFilter<ImageType, > ImageType> ForwardType; > > ForwardType::Pointer forward = ForwardType::New(); > > forward->SetInput(1, inputImage); > > forward->SetInput(0, outputImage); > > forward->Update(); > > > > … > > (writing output image into a file) > > … > > > > The problem is that I get some errors when I instantiate the > forward-Filter. > > The errors are like that one: > > error C2784: "vnl_matrix<T> operator *(const vnl_diag_matrix<T> &,const > vnl_matrix<T> &)": template-Argument für "const vnl_diag_matrix<T> &" > konnte nicht von "vnl_matrix_fixed<T,4,4>" hergeleitet werden. > > > > Maybe there is something wrong with the PixelType that i used? > > > > If the forward projection works, which filter should i use to achieve a > filtered back projection without aliasing? > > > > Many thanks in advance. > > Tobias > > > > *Von:* [email protected] [mailto:[email protected]] *Im Auftrag von *Simon > Rit > *Gesendet:* Mittwoch, 4. Mai 2016 09:02 > *An:* Tobias Stein <[email protected]> > *Cc:* [email protected] > *Betreff:* Re: [Rtk-users] forward and back projection - MITK > > > > Dear Tobias, > > The forward projection has 2 inputs: > > - input 0: the stack of projections in which you wish to forward project, > > - input 1: the volume you wish to forward project. > > For the backprojection, it's exactly the same: > > - input 0: the volume in which you wish to backproject, > > - input 1: the stack of projections you wish to backproject. > > Be aware that JosephBackProjectionImageFilter is the transpose of the > forward projection and will have some aliasing (see, e.g., De Man and Basu > <http://iopscience.iop.org/article/10.1088/0031-9155/49/11/024/meta;jsessionid=87B598ABFDC2AA07D1DED2DEE607F0E7.c2.iopscience.cld.iop.org> > to see what I mean). > > I have personally never used MITK but don't hesitate to share your > experience on the mailing list. > > Simon > > > > On Wed, May 4, 2016 at 12:59 AM, Tobias Stein <[email protected]> > wrote: > > Hi all, > > > > i want to use the forward and backward projection to reduce metal > artefacts in ct images. I’ve seen so far that I may use the > JosephForwardProjectionImageFilter to perform the forward projection. I’ve > also seen the test for this class but I don’t get it, where should i put my > 2D slice as input to execute the transformation. About the back > transformation with the JosephBackProjectionImageFilter I also need more > information how I can use it to transform a sinogram back to a ct slice. > There is a test at the documentation, but the link is missing there. > > > > I’ve got another independent question. > > Are some of you familiar with MITK and know how to get the superbuild up > and running with RTK? I am writing a MITK-Plugin and want to reduce the > metal artifacts before a segmentation. So if some of you have experience > with the combination of RTK and MITK it would be nice if you will share it > ;) > > > > Best regards, > > Tobias > > > > > _______________________________________________ > Rtk-users mailing list > [email protected] > http://public.kitware.com/mailman/listinfo/rtk-users > > > > >
_______________________________________________ Rtk-users mailing list [email protected] http://public.kitware.com/mailman/listinfo/rtk-users
