[BUG]: Issue with “missing channel” in Incoherent algorithm?

[willend]

Describe the bug
Bing Li at ESS came by with an interesting bug that occurs in the case of using Incoherent in a hollow, low-aspect cylinder:

COMPONENT iso_samp1 = Incoherent(
 geometry = 0, radius = 0.012,
 xwidth = 0, yheight = 0.004,
 zdepth = 0, thickness = 0.00025,
 concentric = 0, order=1)
AT (0, 0, 0) RELATIVE sample
ROTATED (0.0, 0.0, 0.0) RELATIVE sample
EXTEND %{
  if (!SCATTERED) ABSORB;
%}

Illuminating a banana-shaped detector around the sample, a band around the scattering plane gets systematically lower intensity:

Interestingly, using e.g. PowderN in the V (Vanadium) configuration but suppressing the powerlines (no reflection list) does not yield this stripe...

COMPONENT iso_samp = PowderN(
geometry = 0, radius = 0.012,sigma_abs=5.08, sigma_inc=5.08,Vc=13.827,
 xwidth = 0, yheight = 0.004,
 zdepth = 0, thickness = 0.00025,
concentric = 0, order=1)
AT (0, 0, 0) RELATIVE sample
ROTATED (0.0, 0.0, 0.0) RELATIVE sample
EXTEND %{
  if (!SCATTERED) ABSORB;
  %}

Document steps to reproduce the bug
Run the attached instrument (PowderN is commented out)

Add your code / dataset
Instr file attached LET_backend.instr.TXT

Speculation:
Problem with neutrons that (either?) enter / exit through the cylinder ‘caps’?

[g5t]

Does this extend to other shapes for the Incoherent component?
I've used it to approximate crystal electric field levels on BIFROST -- which only has detectors in the plane -- but probably used a solid cylinder for the shape.

[willend]

@g5t full extent of the issue remains to be studied, will report back on this issue later with a status and comparison across all geometries hollow/solid with multiple components in each case.

[willend]

Intermediate status report is that among the possible geometries, only the hollow cylinder really stands out. I.e. the solid cylinder seems OK.

That said, I will try to get a full(er) picture via some more rigorous testing across Incoherent/PowderN/Isotropic_sqw and write up new test instrs / a document…

[mads-bertelsen]

I think I see the problem. When forcing the multiple scattering loop to end, the attenuation for the last path out of the geometry is applied with the code below. That does not support the hollow geometries, so if a ray is scattered in a direction that would transverse the hollow part, it is still attenuated as if that was made of the same material.

  // Add attenuation of exit flight path for non-multiple scattering
  if (order && SCATTERED) {
    if (VarsInc.shape == 0)
      intersect = cylinder_intersect (&t0, &t3, x, y, z, vx, vy, vz, radius, yheight);
    else if (VarsInc.shape == 1)
      intersect = box_intersect (&t0, &t3, x, y, z, vx, vy, vz, xwidth, yheight, zdepth);
    else if (VarsInc.shape == 2)
      intersect = sphere_intersect (&t0, &t3, x, y, z, vx, vy, vz, radius);
    #ifdef USE_OFF
    else if (VarsInc.shape == 3)
      intersect = off_intersect (&t0, &t3, NULL, NULL, x, y, z, vx, vy, vz, 0, 0, 0, thread_offdata);
    #endif

    d_path = v * t3; /* Length of full path through sample */

    p_trans = exp (-my_t * d_path);
    p *= p_trans;
    PROP_DT (t3);
  }

[mads-bertelsen]

I spotted a difference in the Incoherent and PowderN interpretation of a hollow cylinder which makes direction comparison between the two components less useful.

For the Incoherent component, there is a top and bottom lid with the thickness specified, while the PowderN version effectively have the hollow part go all the way to the top/bottom.

[mads-bertelsen]

Updating the Incoherent components exit clause that performs attenuation in case of limited order seems to have done the trick, will make a pull request soon.

For a quick comparison, here are the cases:

The broken Incoherent version

The PowderN version referenced in this issue, does not have the line, which is because the PowderN hollow geometry does not have the "lids".

The fixed Incoherent version, notice there is still a line, there should be

A Union reference that only show first order scattering, the line is there as it should be with this geometry (done with twice the ncount as rays are lost to higher order scattering that can not be disabled). The intensities match the fixed version of Incoherent.

[willend]

@mads-bertelsen could you please detail why that stripe is meant to occur?

[mads-bertelsen]

@willend Sure, the stripe in the horizontal plane is from absorption in the "lid" part, so its the stripe at latitude = 0 and y = 0. The bugged version of the component also has inhomogeneity away from that line that was exaggerated by applying attenuation for the part of the path in the hollow part of the geometry.

[willend]

Then I guess that for a more ‘elongated’ cylinder and with sufficient detector resolution, there should in fact be two stripes?

Also, for full circle, could you potentially illustrate with a union-assembly that cylindrical V-surface (i.e. without ‘lids’) behaves like PowderN? (I guess this can be done with a hollow V-cylinder and ‘vacuum’ lids?)

[mads-bertelsen]

@willend Here is a comparison between a Union tube, so no lids, compared to the PowderN. No line in the center anymore. I am surprised to see the intensities do not match, perhaps it is because I only take a single order of the Union output, but thought that would correspond to PowderN supporting just order=1..

[willend]

@mads-bertelsen since @Lomholy did a quite thorough validation of PowderN less than one year ago (to my knowledge in connection with the infamous ’sample paper’, perhaps minor details did not make it across somehow?

Should we perhaps aim for revisiting all ‘comparison-instruments’ in the hope of ensuring agreement wherever we can? (And filling in gaps wherever something can be compared more rigorously that what we are currently doing.)

(One thing that suddenly comes to mind is that Isotropic_Sqw also ‘includes the lids’ - would be good to get a checkmark also there.)

[mads-bertelsen]

@willend Happy to get more cross comparisons. I don't know if @Lomholy checked hollow geometries? It looks like there are very few test cases that use the thickness parameter, perhaps we could add testcases to the incoherent test bunch with all the hollow geometries?

Should we enforce that they have the same choice of interpretation of the thickness parameter so that they should give the same intensity in that test?

[willend]

@mads-bertelsen you are probably right that @Lomholy worked on solid geometries only.

And yes I believe that thickness should mean the same across the samples, at least it was intended to.

[mads-bertelsen]

To help us in making a sensible decision on whether we should include top/bottom for cylinders with a given thickness, I checked the behavior of all the components with the parameters radius and thickness.

Only wall, no top / bottom
Tunneling_sample
Res_sample
PowderN
TOFRes_sample
V_sample

Top / bottom included
Incoherent
Isotropic_Sqw

Supports negative thickness that expand size
PowderN

Shape_simple just calls McDisplay with
cylinder (0, 0, 0, radius, yheight, thickness, 0, 1, 0)

Shape just calls McDisplay with
cylinder (0, 0, 0, radius, yheight, thickness, nx, ny, nz);

simple mcdisplay does not show anything different when adding thickness parameter
mcdisplay-webgl shows geometry without top / bottom as below