<div dir="ltr">Hi all,<div><br></div><div>I've been trying for a while to use RTK in my own software to perform reconstructions. I've looked at the applications and examples to understand how to set it up but I can't get it working.</div><div><br></div><div>I have my own data object which is stored in a single array. I used ITK's Import Image Filter to import the data in RTK. Right now I'm trying to get everything working with my own Shepp-Logan projections. So my pipeline looks like this for now:</div><div>- Create the geometry</div><div>- Create a Constant Image Source</div><div>- Import my projection data</div><div>- Set up the backprojection</div><div>- Take the output of the back-projection and export it to my own data object</div><div><br></div><div>The problem is that all that I get as the output is zeros everywhere, on the whole volume. I've tried skipping the back-projection to see if it was my data that was the problem, but then it works fine and I just get back my projections as the output. I've tried a few different back-projection methods, but it's always giving me zeros as the output.</div><div><br></div><div>Has anybody experienced this before? What did I miss in my setup? I'll put a bit of my code below.
<span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">My guess is that I didn't properly set the input image with the projections.</span>
</div><div><br></div><div>Thanks a lot for your time!</div><div><br></div><div>Adam</div><div><br></div><div>Code:</div><div><br></div><div>// Geometry</div><div><div>rtk::<wbr>ThreeDCircularProjectionGeomet<wbr>ry::Pointer geometry;</div><div>geometry = rtk::<wbr>ThreeDCircularProjectionGeomet<wbr>ry::New();</div><div>for (uint32_t p = 0; p < numberOfProjections; p++)</div><div>{</div><div><span style="white-space:pre-wrap"> </span>double angle = beginAngle + (double)p * (endAngle - beginAngle) / (double)numberOfProjections;</div><div><span style="white-space:pre-wrap"> </span>geometry->AddProjection(sid, sdd, angle, isox, isoy, outAngle, inAngle, sourceOffsetX, sourceOffsetY);</div><div>}</div></div><div>TRY_AND_EXIT_ON_ITK_EXCEPTION(<wbr>geometry->Update());<br></div><div><br></div><div>// Constant Image Source</div><div><div>ConstantImageSourceType::<wbr>SizeType reconstructionSize;</div><div>reconstructionSize[0] = outputChannelSizeVector.getX()<wbr>;</div><div>reconstructionSize[1] = outputChannelSizeVector.getY()<wbr>;</div><div>reconstructionSize[2] = outputChannelSizeVector.getZ()<wbr>;</div><div><br></div><div>ConstantImageSourceType::<wbr>PointType origin;</div><div>origin[0] = outputOriginVector.getX();</div><div>origin[1] = outputOriginVector.getY();</div><div>origin[2] = outputOriginVector.getZ();</div><div><br></div><div>ConstantImageSourceType::<wbr>SpacingType reconstructionSpacing;</div><div>reconstructionSpacing[0] = outputSpacingVector.getX();</div><div>reconstructionSpacing[1] = outputSpacingVector.getY();</div><div>reconstructionSpacing[2] = outputSpacingVector.getZ();</div><div><br></div><div>ConstantImageSourceType::<wbr>Pointer constantImageSource = ConstantImageSourceType::New()<wbr>;</div><div>constantImageSource-><wbr>SetOrigin(origin);</div><div>constantImageSource-><wbr>SetSpacing(<wbr>reconstructionSpacing);</div><div>constantImageSource->SetSize(<wbr>reconstructionSize);</div><div>constantImageSource-><wbr>SetConstant(0.0);</div></div><div><br></div><div>// ITK Import Image Filter</div><div><div>ImportFilterType::SizeType projectionsSize;</div><div>projectionsSize[0] = projectionsX<span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">Size</span>;</div><div>projectionsSize[1] =
<span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">projectionsY</span><span style="font-size:small;text-decoration-style:initial;text-decoration-color:initial;background-color:rgb(255,255,255);float:none;display:inline">Size</span>;</div><div>projectionsSize[2] =
<span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">projectionsZ</span><span style="font-size:small;text-decoration-style:initial;text-decoration-color:initial;background-color:rgb(255,255,255);float:none;display:inline">Size</span><span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">;</span></div><div><br></div><div>ImportFilterType::SpacingType projectionsSpacing;</div><div>projectionsSpacing[0] = projectionsXSpacing;</div><div>projectionsSpacing[1] =
<span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">projectionsYSpacing</span>;</div><div>projectionsSpacing[2] =
<span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">projectionsZSpacing</span>; // not that it matters much...</div><div><br></div><div>const bool importImageFilterWillOwnTheBuf<wbr>fer = false;</div><div><br></div><div>double projectionsOrigin[3] = { 0.0 };</div><div><br></div><div>ImportFilterType::Pointer importFilter = ImportFilterType::New();</div><div>ImportFilterType::IndexType importFilterStartIndex;</div><div>importFilterStartIndex.Fill(0)<wbr>;</div></div><div><br></div><div><div>ImportFilterType::RegionType region;</div><div>region.SetIndex(<wbr>importFilterStartIndex);</div><div>region.SetSize(<wbr>projectionsSize);</div><div><br></div><div>importFilter->SetRegion(<wbr>region);</div><div>importFilter->SetOrigin(<wbr>projectionsOrigin);</div><div>importFilter->SetSpacing(<wbr>projectionsSpacing);</div></div><div><br></div><div><div>PixelType *projectionsDataBuffer = (PixelType *)projectionsChannelFloat.<wbr>getTRawDataChunk(0);</div><div>importFilter-><wbr>SetImportPointer(<wbr>projectionsDataBuffer, projectionsSize[0]* projectionsSize[1] * projectionsSize[2], importImageFilterWillOwnTheBuf<wbr>fer);</div></div><div><br></div><div>TRY_AND_EXIT_ON_ITK_EXCEPTION(<wbr>importFilter->Update());<br></div><div><br></div><div>// Backrprojection setup</div><div><div>rtk::<wbr>BackProjectionImageFilter<<wbr>OutputImageType, OutputImageType>::Pointer bp;</div><div>bp = rtk::<wbr>FDKBackProjectionImageFilter<<wbr>OutputImageType, OutputImageType>::New();</div></div><div><br></div><div><div>bp->SetInput(0, constantImageSource-><wbr>GetOutput());</div><div>bp->SetInput(1, importFilter->GetOutput());</div><div>bp->SetGeometry(geometry);</div><div>TRY_AND_EXIT_ON_ITK_EXCEPTION(<wbr>bp->Update());</div></div><div><br></div><div>OutputImageType::<wbr>PixelContainer::Element * iterator;<br></div><div>iterator = bp->GetOutput()-><wbr>GetBufferPointer();<br></div><div><br></div><div><div>PixelType *outputBuffer = (PixelType *)outputChannel.<wbr>getTRawDataChunk(0);</div><div>uint64_t numberOfPixels = xSize * ySize * zSize;</div><div>for (uint64_t i = 0; i < numberOfPixels; i++)</div><div>{</div><div><span style="white-space:pre-wrap"> </span>
<span style="font-size:small;background-color:rgb(255,255,255);text-decoration-style:initial;text-decoration-color:initial;float:none;display:inline">outputBuffer<span> </span></span>[i] = *iterator;</div><div><span style="white-space:pre-wrap"> </span>iterator++;</div><div>}</div></div></div>