Reinvert spect axis

[danieldeidda]

Changes in this pull request

Testing performed

Related issues

Checklist before requesting a review

• [] I have performed a self-review of my code
• [] I have added docstrings/doxygen in line with the guidance in the developer guide
• [] I have implemented unit tests that cover any new or modified functionality (if applicable)
• [] The code builds and runs on my machine
• [] documentation/release_XXX.md has been updated with any functionality change (if applicable)

Contribution Notes

Please tick the following:

• [ v] The content of this Pull Request (the Contribution) is intentionally submitted for inclusion in STIR (the Work) under the terms and conditions of the Apache-2.0 License.
• [v ] I (or my institution) have signed the STIR Contribution License Agreement (not required for small changes).

[KrisThielemans]

Did you run the recon_test_pack SPECT tests?

It's possible that the SRT2DSPECT code needs adaptation. @Dimitra-Kyriakopoulou in our current SPECT projector, the "x" (i.e. horizontal) axis is flipped between emission and attenuation. This PR fixes this. Not sure if you had modified your code to be compatible with the STIR SPECTUB projector.

[danieldeidda]

no sorry i didn't have STIR installed. I'll check on another PC

[Dimitra-Kyriakopoulou]

Dear Professor @KrisThielemans,
simulate_data_for_tests.sh made up for the fliplr of SPECTUB.

I) As long as line 113 of the current simulate_data_for_test.sh in master
atten_input_file="my_atten_image_SPECT.hv"
is replaced by
atten_input_file="my_atten_image.hv"

the SRT2DSPECT error should not exist any more.

II) I cannot test it immediately now, as the inversion is not in master; but if the error is not corrected, i will add in the GRD2D-DDSR2D repository the inversion changes to test it.

III) What are further redundant and need to be removed -but do not affect the tests- are:

1. lines 93-94 of the current simulate_data_for_test.sh in master
2. the generate_atten_cylinder_SPECT.par in recon_test_pack folder.

I am really sorry for the late reply,
Dimitra

[danieldeidda]

Dear Professor @KrisThielemans, simulate_data_for_tests.sh made up for the fliplr of SPECTUB.

I) As long as line 113 of the current simulate_data_for_test.sh in master atten_input_file="my_atten_image_SPECT.hv" is replaced by atten_input_file="my_atten_image.hv"

the SRT2DSPECT error should not exist any more.

II) I cannot test it immediately now, as the inversion is not in master; but if the error is not corrected, i will add in the GRD2D-DDSR2D repository the inversion changes to test it.

III) What are further redundant and need to be removed -but do not affect the tests- are:

1. lines 93-94 of the current simulate_data_for_test.sh in master
2. the generate_atten_cylinder_SPECT.par in recon_test_pack folder.

I am really sorry for the late reply, Dimitra

OK thanks Dimitra, I will have a look

[danieldeidda]

@KrisThielemans is there a way to see this log file? "There were errors. Check my_SRT2DSPECT_test_sim.par_SPECT.log" I don't seem to be bale to download the artifacts

[Dimitra-Kyriakopoulou]

Dear @danieldeidda,
Would you mind perhaps trying the change I suggested?
Specifically, on line 113 of the current simulate_data_for_test.sh, could you please replace atten_input_file="my_atten_image_SPECT.hv"
by
atten_input_file="my_atten_image.hv" ?

Kind regards,
Dimitra

[KrisThielemans]

@KrisThielemans is there a way to see this log file? "There were errors. Check my_SRT2DSPECT_test_sim.par_SPECT.log" I don't seem to be bale to download the artifacts

go to the "summary" page of the run of the action (in this case, https://github.com/UCL/STIR/actions/runs/18588852662?pr=1631). scrolldown a bit. There's zip files with all logs.

[danieldeidda]

Hi,
so I have installed this and tried the recon_test_pack.
I changed the simulation according to Dimitra's suggestion but the SPECT recons have the attenuation image inverted between org and out (folders).
I am not sure I understand when the org files are created but they are created with the axis inverted and don't seem to update.
I don't find anywhere in the recon_test any occurrence of invert_axis.

[KrisThielemans]

recon_test_pack/SPECT/SPECTUB is a bit of a mess regarding flipping/rotation, but we should avoid checking in new binary files.

As this is a test that says "do I get the same as before",we should invert_axis the attenuation image, as that's used in

STIR/recon_test_pack/SPECT/SPECTUB/OSEM_2DPSF.par

Line 26 in 4e9666a

attenuation map := attMapRec.hv

so I think this needs

• modify run_SPECTUB_tests.sh to x-flip attMapRec.hv to a new file
• modify *.par to use that new file

This test doesn't check SRT2DSPECT, so no changes for Dimitra.

On the other hand, we have `run_test_simulate_and_recon.sh, which first simulates some data

STIR/recon_test_pack/run_test_simulate_and_recon.sh

Line 91 in 4e9666a

./simulate_data_for_tests.sh --SPECT --suffix "$SPECT_suffix"

and then reconstructs it. For OSEM, whatever flipping won't matter, as both simulation and recon use the same SPECTUB projector. However, SRT2DSPECT does not use SPECTUB, so if we flip SPECTUB, simulate with SPECTUB, the recon will be correct. Therefore, we should flip the attenuation image used for the simulation (the OSEM test won't care, but SRT2DSPECT will). I guess that's what @Dimitra-Kyriakopoulou is saying.

[Dimitra-Kyriakopoulou]

Dear Professor @KrisThielemans and @danieldeidda,
I) Apologies -I thought the inversion had already been finalized. The SRT2DSPECT change will come after that.

II) I agree with Professor Thielemans' proposal.

In addition, I would like to express a few thoughts that came to my mind regarding the 90 degree rotation and z-flip issues, which can appear as the current left–right (x) problem:

A) src/IO/InterfilePDFSHeaderSPECT.cxx

num_segments = 1;
num_views = -1;
add_key("number of projections", &num_views);
start_angle = 0;
add_key("start angle", &start_angle);
direction_of_rotation = "cw";
add_key("direction of rotation", &direction_of_rotation);
extent_of_rotation = double_value_not_set;
add_key("extent of rotation", &extent_of_rotation);

Different datasets may use different 0° definitions and CW/CCW sense. That creates a global rotation (often ~90°) which can be misread as “left/right” mismatch.

Hence, right after the add_key(...) calls, we must normalize angles to one convention (0°=+x, CCW, positive increment). If already canonical, it’s a no-op; otherwise it corrects once and removes rotation ambiguity.

B) src/recon_buildblock/ProjMatrixByBinSPECTUB.cxx

vol.Xcmd2 = vol.Ncold2 * vol.szcm; // half-size x (cm)
vol.Ycmd2 = vol.Nrowd2 * vol.szcm; // half-size y (cm)
vol.Zcmd2 = vol.Nslid2 * vol.thcm; // half-size z (cm)

vol.x0 = ((float)0.5 - vol.Ncold2) * vol.szcm; // x of first voxel
vol.y0 = ((float)0.5 - vol.Nrowd2) * vol.szcm; // y of first voxel
vol.z0 = ((float)0.5 - vol.Nslid2) * vol.thcm; // z of first voxel

Some images arrive with opposite z slice order . That appears as a z-mirror during comparison and can confuse left/right assessments.

At image load (or just before the projector binds the VoxelsOnCartesianGrid), we can add a one-time z-order normalization that reorders slices to canonical +z.

C) src/recon_buildblock/SPECTUB_Tools.cxx

const float deg = wmh.prj.ang0 + (float)i * wmh.prj.incr; // angle (deg)
ang[i].cos = cos(deg * dg2rd);
ang[i].sin = sin(deg * dg2rd);

// reduced angle to 0..45°, track quadrant
float angR = fabs(deg);
int   quad = (int)floor(angR / 90.f);
angR       = fabs(angR - 90.f * quad);
if (angR > 45.f) angR = fabs(90.f - angR);

// detector-line direction and first-bin position
ang[i].incx  = wmh.prj.szcm * ang[i].cos;  // Δx per bin
ang[i].incy  = wmh.prj.szcm * ang[i].sin;  // Δy per bin

ang[i].xbin0 = -ang[i].Rrad * ang[i].sin - wmh.prj.lngcmd2 * ang[i].cos;
ang[i].ybin0 =  ang[i].Rrad * ang[i].cos - wmh.prj.lngcmd2 * ang[i].sin;

If (A) isn’t normalized, deg carries incompatible angle conventions, making the geometry rotated (which may be perceived as left/right mismatch).

We must ensure deg = ang0 + i*incr is exactly the canonical CCW, 0°=+x delivered by (A):

• If A is in place → no code change usually needed in C.
• If any extra CW/CCW or handedness reinterpretation exists here → remove it so C simply uses the canonical angles.

D) src/recon_buildblock/SPECTUB_Weight3d.cxx

float ux = bin.x - vox.x;
float uy = bin.y - vox.y;
float uz = bin.z - vox.z;

int signx = SIGN(ux);
int signy = SIGN(uy);
int signz = SIGN(uz);

// ... next-plane distances and stepping ...
dx = (next_x - vox.x) / (ux + EPSILON);
dy = (next_y - vox.y) / (uy + EPSILON);
dz = (next_z - vox.z) / (uz + EPSILON);

If z-order is normalized once at (B), D is automatically correct for all cases.

[KrisThielemans]

Thanks @Dimitra-Kyriakopoulou. These are good points, but I think they belong in a new issue, not in a comment in this (intentionally very small) PR which attempts to fix the SPECTUB projector and hence all iterative algorithms. I've created #1635 for that.

[KrisThielemans]

run_SPECTUB_tests.sh, need to essentially revert ac9cc0f, so still need to remove this line

STIR/recon_test_pack/SPECT/SPECTUB/README.txt

Line 8 in c6ee75a

In addition, the reconstructed images are inverted along x. So the script inverts the reconstructed images along x.

[KrisThielemans]

Clean-up:

• remove line

  STIR/recon_test_pack/simulate_data_for_tests.sh

  Line 94 in c6ee75a

  generate_image generate_atten_cylinder_SPECT.par

  and associated generate_atten_cylinder_SPECT.par
• as per Dimitra's suggestion: modify

  STIR/recon_test_pack/simulate_data_for_tests.sh

  Line 113 in c6ee75a

  atten_input_file="my_atten_image_SPECT.hv"

  
  to use atten_input_file="my_atten_image.hv"
• rename my_atten_image_SPECT_modified.hv to my_atten_image_SPECT.hv (at least 2 places)

[KrisThielemans]

This is a bit confusing (and has a few typos). Are we now inverting the x-axis in the reconstructed emission image? I hope not. If not, does the following describe it better?

Suggested change

	Fixed a bug in the SPECTUB matrix (introduced in 4.0) that inverted the SPECT reconstruction along the axis.
	This had resulted with an inconsitent attenuation correction orientation compared to the reonctructed image.
	The axis inversion is now reverted, so if your scripts were using an inverted attenuation map they need to be
	updated with this release.
	See <a href="https://github.com/UCL/STIR/pull/1631"">PR #1631</a>" . </li>
	Fixed a bug in the SPECTUB matrix (introduced in 4.0) that effectively inverted the "x" (horizontal) axis for the
	attenuation image. This resulted in an inconsistent attenuation image orientation compared to the reconstructed emission image.
	This axis inversion is now reverted, so if your scripts are using an attenuation map, you need to ***remove**
	calls to `invert_axis x` (or code to that effect) to get correct results.<br>
	See <a href="https://github.com/UCL/STIR/pull/1631"">PR #1631</a>".
	</li>

[danieldeidda]

yes sorry

[KrisThielemans]

hmmm. I think we are indeed x-flipping the emission image. Worse than I thought.

[KrisThielemans]

We need to be 100% clear in our release notes, and I have doubts that my suggested edit is ok.

Also, you will need to run pre-commit run -all, commit, and then merge master.

[KrisThielemans]

sorry. forward projection now fails in run_test_simulate_and_recon.sh

[Dimitra-Kyriakopoulou]

Dear Professor @KrisThielemans,
my_atten_image_SPECT_modified.hv files are created by commenting out the line with the modality; this was done specifically so forward_project would not fail, as it expects PET.

Therefore, could you please change only line 113 of simulate_data_for_tests.sh in master from
atten_input_file="my_atten_image_SPECT.hv"
to
atten_input_file="my_atten_image.hv"?

[KrisThielemans]

@Dimitra-Kyriakopoulou we renamed the file in the script (dropping "_modified"), but we forgot to update the .par file. Should be fine now.

@danieldeidda I also did some more clean-up, and adjusted the release notes again according to my latest understanding. Please check!

This is a rather drastic change however. It's going to make lots of people unhappy...

[Dimitra-Kyriakopoulou]

Dear Professor @KrisThielemans,
Motivated by your idea of x-flipping attMapRec.hv, I noticed that the tests pass without changing simulate_data_for_tests.sh -provided we update run_SPECTUB_tests.sh beyond just removing the flip.

In hindsight, after the attenuation-related change did not resolve the issue, it became clear this wasn’t a numerical issue but a baseline mismatch: the org/ images still contained the previous x-flip, while the new outputs (with the flip removed) are in the true geometry.

After refreshing the baselines to the current, unflipped geometry, the tests passed. The attached updated script can run in two ways:

A) Refresh mode (run with --refresh-org)

It intentionally updates the references to the new convention, ie. after removing the previous x-flip:

• Reconstructs images into out/.
• Renames the current org/ to a timestamped backup org_legacy_YYYYMMDD_HHMMSS/ (so you don’t lose the previousbaselines).
• Copies out/*.hv and out/*.v into org/ (these become the new references).
• Stops (no compare on this run).
  Exit code: 0 if refresh succeeded; non-zero if a reconstruction failed.

-- Next normal run will compare against the new org/.

B) Normal mode (run with no flag)
Comparison test only — does not change references.

Run TWICE: . If the org/ baselines are the previous, x-flipped references, a plain compare will fail until you update them.
Hence, i) refresh org/ once (to replace the previous flipped baselines with the new, unflipped references), and ii)compare normally from then on.

THANK YOU SO MUCH FOR READING THIS MESSAGE!
run_SPECTUB_tests3.sh

[KrisThielemans]

Thanks @Dimitra-Kyriakopoulou I've resolved this in another way. (I won't to avoid modifying the org data as long as possible in any case). You can check the history to see what I changed.

We just have a problem due to #1639 which will be resolved when merging with master again.

[Dimitra-Kyriakopoulou]

Dear Professor @KrisThielemans,
I had seen the issue had been resolved, but wished to comment on the org/ change being equivalent to the simulate_data_for_tests.sh change — i.e., one or the other needs to change.

For the clean up, it may also be safe to remove the generate_atten_cylinder_SPECT.par in recon_test_pack folder.

THANK YOU WHOLEHEARTEDLY FOR YOUR COMMENT!!!

[danieldeidda]

@KrisThielemans i think the notes are OK. @Dimitra-Kyriakopoulou thanks I removed the parfile too.

[KrisThielemans]

@danieldeidda the effect of this change is rather dramatic. It is going to break @samdporter and @varzakis work for the STIR User's meeting. I think we should postpone it for after the release (i.e. until after MIC).

@Dimitra-Kyriakopoulou I then plan to add a note that SRT2DSPECT gives the correct orientation.

[KrisThielemans]

AppVeyor builds fails due to #1642. Please ignore.

[danieldeidda]

OK for me. Is there anything to do before then on this?

[KrisThielemans]

we'll have to merge master again, and move the text in the release notes, but let's do that after the release.

sorry