[Rtk-users] Error Result by using FDK algorithm to Reconstruction 2D X- ray Scan projection images into 2D CT Slice image
Simon Rit
simon.rit at creatis.insa-lyon.fr
Fri Jan 13 10:05:45 CET 2023
Hi,
So this is the same thing for the projections, you probably want to center
them and set origin to (1152-1)*0.25*-0.5=-143.875. This can be set with
reader->SetOrigin. The third coordinate of the origin vector is not used in
the projections images, you can leave it at 0.
Next is the angular arc. It is less than 360°. You will need to use the
Parker short scan weighting to reduce these artefacts, see here
<https://github.com/SimonRit/RTK/blob/master/applications/rtkfdk/rtkfdk.cxx#L100-L115>
a code example in C++.
Best regards,
Simon
On Fri, Jan 13, 2023 at 9:30 AM 何明哲 <m10512067 at yuntech.org.tw> wrote:
> Hi Simon
>
> Sure, here is those information for reader->Getoutput()
>
> Size(XYZ)-> 1152,1152,270
> Spacing(XYZ)-> 0.25,025,0.25.
> Origin(XYZ)-> 143.875,143.875,143.875
>
> but your suggestion inspires me a lot ,
>
> I was use add TiffimageIO to set size spacing origin to projectionreader
> source
>
> finally ,I can get this (see attached file )
>
> Although the image still has a lot of artifact and is not clear enough,
>
> Did you have any suggestion for this?
>
> but finally it is a huge breakthrough for me
>
> thanks for help and reply again
>
>
>
>
> Simon Rit <simon.rit at creatis.insa-lyon.fr> 於 2023年1月12日 週四 下午10:57寫道:
>
>> Can you provide the size, origin and spacing of reader->GetOutput()? It's
>> a bit hard to help you without this information.
>> Simon
>>
>> On Thu, Jan 12, 2023 at 11:01 AM 何明哲 <m10512067 at yuntech.org.tw> wrote:
>>
>>>
>>> Dear Simon
>>>
>>> Thanks for your reply .Its great help for me!
>>>
>>> 1.)I was refer your suggestion and modified the code,
>>> but the result will be a large blank area and I don't understand why
>>>
>>> here is the code and result
>>>
>>> sizeOutput[0] = Imagewidth;
>>> sizeOutput[1] = 1;
>>> sizeOutput[2] = numberOfProjections;
>>> spacing[0] = 0.25;
>>> spacing[1] = 0.25;
>>> spacing[2] = 0.25;
>>> origin[0] = (Imagewidth-1)*0.5*spacing[0];
>>> origin[1] = (Imageheigh-1)*0.5*spacing[1];
>>> origin[2] = spacing[2]*SliceN*-1;
>>>
>>> recoVolume->SetSize(sizeOutput);
>>> recoVolume->SetSpacing(spacing);
>>> recoVolume->SetOrigin(origin);
>>> recoVolume->SetConstant(0.);
>>>
>>>
>>>
>>> 2.)if you're doing 2D reconstruction, be careful that your projections
>>> should not be 1D because the backprojection uses a 2D interpolation. It
>>> should have at least 2 lines.
>>> --> Yes , I hope I can eventually reconstruct a 2D CT Slice, if my
>>> read projection step or setting have any problem please let me know,
>>>
>>> thanks for your reply
>>>
>>> BR,
>>>
>>>
>>>
>>>
>>>
>>> Simon Rit <simon.rit at creatis.insa-lyon.fr> 於 2023年1月12日 週四 下午4:30寫道:
>>>
>>>> Hi,
>>>> I'm not sure I have enough information to answer... What I can say is:
>>>> - that the origin is not set to center your volume, for centering, you
>>>> should use (size-1)*0.5*spacing
>>>> - if you're doing 2D reconstruction, be careful that your projections
>>>> should not be 1D because the backprojection uses a 2D interpolation. It
>>>> should have at least 2 lines.
>>>> I hope it helps,
>>>> Simon
>>>>
>>>> On Thu, Jan 12, 2023 at 8:46 AM 何明哲 <m10512067 at yuntech.org.tw> wrote:
>>>>
>>>>>
>>>>> HI Everyone ,
>>>>>
>>>>> I have a question
>>>>> that When I was refer to the official example
>>>>> "Firstcudareconstruction.cpp",
>>>>> then I tried to modify the code and reconstruction 2D x ray
>>>>> projection image to 2D ct slice,
>>>>> the code is worked but I got the wrong result too,
>>>>> I don't understand what the problem I'm have ?
>>>>> it's possible that the setting of reconstruction output image size or
>>>>> origin is wrong?
>>>>>
>>>>> here is my code and result (I was tried to only reconstruction just a
>>>>> slice NO.100 of the whole volume)
>>>>>
>>>>> --------------------------------------------------------------------------------------------
>>>>> using GeometryType = rtk::ThreeDCircularProjectionGeometry;
>>>>> GeometryType::Pointer geometry = GeometryType::New();
>>>>> unsigned int numberOfProjections = 280;
>>>>> double firstAngle = 0;
>>>>> double angularArc = 280;
>>>>> unsigned int sid = 510;
>>>>> unsigned int sdd = 690;
>>>>>
>>>>> for (unsigned int noProj = 0; noProj < numberOfProjections; noProj++)
>>>>> {
>>>>> double angle = firstAngle + noProj * angularArc /
>>>>> numberOfProjections;
>>>>> geometry->AddProjection(sid, sdd, angle);
>>>>> }
>>>>>
>>>>> rtk::ThreeDCircularProjectionGeometryXMLFileWriter::Pointer
>>>>> xmlWriter;
>>>>> xmlWriter =
>>>>> rtk::ThreeDCircularProjectionGeometryXMLFileWriter::New();
>>>>> xmlWriter->SetFilename(".\\RTK_Geometry.xml");
>>>>> xmlWriter->SetObject(geometry);
>>>>> xmlWriter->WriteFile();
>>>>>
>>>>> using NameGeneratorType = itk::NumericSeriesFileNames;
>>>>> NameGeneratorType::Pointer nameGenerator = NameGeneratorType::New();
>>>>>
>>>>> nameGenerator->SetSeriesFormat(ProjectionInage +"\\%d.tif");
>>>>> nameGenerator->SetStartIndex(1);
>>>>> nameGenerator->SetEndIndex(280);
>>>>> nameGenerator->SetIncrementIndex(1);
>>>>>
>>>>> using IutputImageType = itk::CudaImage<float, 3>;
>>>>> using ReaderType = rtk::ProjectionsReader<IutputImageType>;
>>>>> ReaderType::Pointer reader = ReaderType::New();
>>>>> reader->SetFileNames(nameGenerator->GetFileNames());
>>>>> reader->Update();
>>>>>
>>>>> using ConstantImageSourceType =
>>>>> rtk::ConstantImageSource<OutputImageType>;
>>>>> ConstantImageSourceType::PointType origin;
>>>>> ConstantImageSourceType::SpacingType spacing;
>>>>> ConstantImageSourceType::SizeType sizeOutput;
>>>>> ConstantImageSourceType::Pointer recoVolume =
>>>>> ConstantImageSourceType::New();
>>>>>
>>>>> sizeOutput[0] = Imagewidth;
>>>>> sizeOutput[1] = Imageheigh;
>>>>> sizeOutput[2] = numofproj;
>>>>> spacing[0] = 0.25;
>>>>> spacing[1] = 0.25;
>>>>> spacing[2] = 0.25;
>>>>> origin[0] = Imagewidth*spacing[0]*0.5*-1;
>>>>> origin[1] = Imageheigh*spacing[0]*0.5*-1;
>>>>> origin[2] = spacing[2]*SliceN*-1;
>>>>>
>>>>> recoVolume->SetSize(sizeOutput);
>>>>> recoVolume->SetSpacing(spacing);
>>>>> recoVolume->SetOrigin(origin);
>>>>> recoVolume->SetConstant(0);
>>>>>
>>>>> using FDKGPUType = rtk::CudaFDKConeBeamReconstructionFilter;
>>>>> FDKGPUType::Pointer feldkamp = FDKGPUType::New();
>>>>>
>>>>> feldkamp->SetInput(0, recoVolume->GetOutput());
>>>>> feldkamp->SetInput(1, reader->GetOutput());
>>>>> feldkamp->SetGeometry(geometry);
>>>>> feldkamp->SetNumberOfThreads(7);
>>>>> feldkamp->Update();
>>>>>
>>>>> using WriterType = itk::ImageFileWriter<OutputImageType>;
>>>>> WriterType::Pointer writer = WriterType::New();
>>>>> writer->SetFileName("C:\\RTK.mhd");
>>>>> writer->SetInput(feldkamp->GetOutput());
>>>>> writer->Update();
>>>>>
>>>>>
>>>>>
>>>>>
>>>>> it should be look like this
>>>>>
>>>>>
>>>>> thanks for your reply
>>>>>
>>>>> BR,
>>>>>
>>>>>
>>>>>
>>>>> _______________________________________________
>>>>> Rtk-users mailing list
>>>>> rtk-users at openrtk.org
>>>>> https://www.creatis.insa-lyon.fr/mailman/listinfo/rtk-users
>>>>>
>>>>
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://www.creatis.insa-lyon.fr/pipermail/rtk-users/attachments/20230113/91de5e87/attachment.htm>
More information about the Rtk-users
mailing list