Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

values/pages from SH12A bdo file are not pr particle, bug? #400

Open
kkorsgaard opened this issue Jun 7, 2021 · 10 comments
Open

values/pages from SH12A bdo file are not pr particle, bug? #400

kkorsgaard opened this issue Jun 7, 2021 · 10 comments
Assignees

Comments

@kkorsgaard
Copy link
Contributor

Hi
After I should compare some different versions of Shieldhit12A I did notice that different identical setup had different values for dose depending on how many simulated particles I did use.

The problem

I did find out that the reason for the different values are that pymchelper does output value from Shieldhit12A multiplied with how many particles I simulate pr simulation.

This means that even though I combine 200 simulations each with 100k primaries it only gives the total for ONE simulation. Very unintuitive if it is a feature and not a bug. If it should give the total I would like it to be the TOTAL for all simulations (primaries) I ran and independent of how I do arrive to this number.

Goal

But anyway I would rather have it to output value (e.g. dose) pr simulated particle as I thought it did. This is what Shieldhit12A do if I output as a txt file.

Reproduce the problem

I can easy reproduce the problem. Just by running a simulation and having both a bdo file and txt file saved. This is also the reason that I do say pymchelper multiply all values with number of simulated primaries as the txt and bdo files that way get identical (I multiply the txt values with number of primaries). Furthermore I am using the method fromfile(). Convertmc image gives the same plots/ values. I have attached a txt file and bdo file which are from the same run and setup.
files.zip

Tested values

I have tested dose, Nkerma, Energy and materiel. It should at least not multiply material and zone at it does not depend on the number of simulated particles, so it indicates that it is a bug.

Version

I am using pymchelper version 1.5.1 and Shieldhit version: v0.9.1-6-gc0599d8.

@kkorsgaard
Copy link
Contributor Author

As mentioned in post #395 when doing parallelization not all jobs will finish for a given time limit. So I guess which could lead to some potential errors - if pymchelper in 1/10 of the cases gives the total dose for 100k particles and in 9/10 gives the total doses for 200k which would lead to some large uncentainties.
However I don't think it is something I have seen, but the question remains how does pymchelper handle this?

I cannot think of a very intuitive way (I would say it did something wrong if it e.g. multiplied all values with 200k even if some only have run with 100k primaries). So dose pr simulated particle would give more sense in these cases.

@nbassler
Copy link
Member

nbassler commented Jun 8, 2021

@grzanka it might be helpful for pymchelper to add a normalization tag on SH side.
SH_NORM_NONE
SH_NORM_PER_PARTICLE
... ?

@nbassler
Copy link
Member

nbassler commented Jun 9, 2021

Not directrly related, but SHBDO_PAG_NORMALIZE tag was not set previously on SH12A.

I propose to use it, then only SH12A needs to keep track of which detectors need to be normalized and which one do not.
This is the MR which activates the tag: https://gitlab.au.dk/shieldhit/shieldhit/-/merge_requests/259

@nbassler
Copy link
Member

nbassler commented Jun 10, 2021

@grzanka following our discussion, I realize we have 3 cases how to deal with results in .bdo.

  • Assuming independent simulation j
  • The number of primaris in this simulation is I_j
  • The raw result (not per particle) is in the BDO file as x_j
  • The user requested results is X

so, we need these averaging cases:

1) SUM:          X = sum_j x_j                    # this is e.g for COUNT, COUNTPOINT and other artificial stuff.
2) NORMALIZED:   X = sum_j x_j / sum_j I_j        # this changes unit into "per primary". Dose, Fluence etc.
3) AVERAGED:     X = (sum_j x_j * sum_j I_j) / sum_j I_j      #  unit is unchanged. dLET, tLET, ...z2/beta2

I do not see any difference in fluence or dose averaged here, since the feld is identical for the voxel of the constiuent compoents, so
fluence = kf * I_j
and
dose = kd * I_j

kf and kd are constants and independent of j (assuming sufficent statistics). Therefore they cancels out in case 3) if you do dose or fluence averaged, respectively.

SHBDO_PAG_NORMALIZE is a char, which can carry the modus.

@nbassler
Copy link
Member

We can even add a 4th case, where just the first result is used.

0) X = x_0    # for GEOMAP scorers, they are always the same for all j.

@nbassler
Copy link
Member

nbassler commented Jun 14, 2021

https://gitlab.au.dk/shieldhit/shieldhit/-/merge_requests/259

new tag behaviour (was not set at all before)

    SHBDO_PAG_NORMALIZE,             /* Flags for page->postproc on how to postprocess the data in SHBDO_PAG_DATA

                                        Given:
                                        - the data in page->data as x_j
                                        - for j instances of this simulation
                                        - which was done with I_j number of paritcles.

                                        The resulting data will be termed X and has the units given by SHBDO_PAG_DATA_UNIT

                                        0: X = x_1                                for GEOMAP type scorers
                                        1: X = sum_j x_j                          COUNT, ...
                                        2: X = (sum_j x_j) / (sum_j I_j)          NORMCOUNT, ...
                                        3: X = (sum_j x_j * I_j) / (sum_j I_j)    LET, ...

                                        (Defines are in score/sh_scoredef.h)
                                      */

@grzanka grzanka self-assigned this Jan 7, 2022
@grzanka grzanka pinned this issue Jan 7, 2022
@grzanka grzanka linked a pull request Jan 7, 2022 that will close this issue
@grzanka grzanka removed a link to a pull request Jan 14, 2022
@grzanka
Copy link
Contributor

grzanka commented Jan 14, 2022

Part of the scorers (code 2) will be fixed in PR #521

@nbassler
Copy link
Member

related: #628 which adds another averaging tag.

@reviewpad
Copy link
Contributor

reviewpad bot commented Apr 19, 2023

AI-Generated Summary: This issue describes a problem with the output values of pymchelper when processing Shieldhit12A simulations. The user reports that different identical setups produce varied dose values depending on the number of simulated particles. It seems that pymchelper multiplies the output values by the number of simulated particles, which leads to incorrect results. The goal is to output the value (e.g., dose) per simulated particle as Shieldhit12A does when outputting as a text file. The issue can be reproduced by running a simulation and having both a BDO and TXT file saved. The user has tested this for dose, Nkerma, Energy, and material parameters. The issue is observed in pymchelper version 1.5.1 and Shieldhit version v0.9.1-6-gc0599d8.

@grzanka
Copy link
Contributor

grzanka commented Mar 26, 2024

Should be fixed in PR #723

Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
None yet
Projects
None yet
Development

No branches or pull requests

3 participants