verify and document up_frac_i and up_frac_f calculations for magnetic models

[pkienzle]

The weights for each cross section computed from up_frac_i and up_frac_f need to be checked. The code uses:

    wuu = sqrt(up_i * up_f)
    wud = sqrt(up_i * (1-up_f))
    wdu = sqrt((1-up_i) * up_f)
    wdd = sqrt((1-up_i) * (1-up_f))
    I = wuu*Iuu + wud*Iud + wdu*Idu + wdd*Idd 

If either up_i or up_f is not 0 or 1, then the sum of the weights does not equal one which does not sound correct. That is, if up_i/up_f represents leakage between channels, the total number of neutrons should be preserved.

Consider up_i = 0.9, which lets through 90% up and 10% down from the guides (so up/(up+down)=0.9, and similarly up_f=0.9, which lets through 90% up and 10% down after the sample. Tracing the individual channels:

initial =>  sample   =>  final selection

 up=0.9 =>  0.9*Iuu  =>  0.9*0.9*Iuu = 0.81*Iuu   [true pos]
                         0.1*0.9*Iuu              [false neg]
            0.9*Iud  =>  0.1*0.9*Iud = 0.09*Iud   [false pos]
                         0.9*0.9*Iud              [true neg] 
 dn=0.1 =>  0.1*Idu  =>  0.9*0.1*Idu = 0.09*Idu   [true pos]
                         0.1*0.1*Idu              [false neg]
            0.1*Idd  =>  0.1*0.1*Idd = 0.01*Idd   [false pos]
                         0.9*0.1*Idd              [true neg]

So weights should not use sqrt. Note that the false positive and false negative rates may be independent, and further, with He3 polarizers, may be time dependent, but we keep things simple by having only one number per polarizer rather than two.

Doing the same analysis with a partially polarized or unpolarized beam requires an additional factor of two intensity correction. In this case, with 90% efficiency on the front end and 50% on the back end:

initial  => sample  => final selection
 up=0.9 =>  0.9*Iuu  =>  0.5*0.9*Iuu = 0.45*Iuu
                         0.5*0.9*Iuu
            0.9*Iud  =>  0.5*0.9*Iud = 0.45*Iud
                         0.5*0.9*Iud
 dn=0.1 =>  0.1*Idu  =>  0.5*0.1*Idu = 0.05*Idu
                         0.5*0.1*Idu
            0.1*Idd  =>  0.5*0.1*Idd = 0.05*Idd
                         0.5*0.1*Idd

For non-magnetic systems Iuu=Idd and Iud=Idu=0. In the first example, a correction of 1/(up_i*up*f + (1-up_i)*(1-up_f)) is needed so that I=Iuu=Iud. In the second example the correction is 2.

Update the user guide to contain complete details of the final implementation and how it is to be interpreted.

[pkienzle]

This code exists in sascalc/calculator/c_extensions/libfunc.c:cal_msld and sasmodels/kernel_iq.c:set_spin_weights

[pkienzle]

This analysis works regardless of whether selectivity is based on spin leakage in the polarizer/analyser or whether it is due to depolarization in the guide field.

On the front end, leakage of down through the front polarizer can give the same up/(up+down) ratio entering the sample as a perfect up selection in the polarizer followed by percentage depolarization. On the back end with non-spin-flip (NSF) samples, false positive down neutrons from Idd are indistinguishable from a portion of the down neutrons from Idd flipping to up and being accepted by a perfect analyser. For spin-flip (SF) samples with Iud==Idu (the usual case) depolarization is again indistinguishable from analyser selectivity. Only case that may cause problems is SF with Iud!=Idu since then it matters whether the sample sees up vs down neutrons, but this never occurs with the simple model of magnetism we are using (the imaginary sld is +mz for Iud and -mz for Idu so the squared contrast is always the same). Such samples will need specialized 2D models which allow magnetism that varies within the sample, and these can define their own notion of spin up fraction if they need it.

[awashington]

I spoke with Rob Dalgliesh about this earlier. We agree that we should be calculating intensities and not amplitudes, so we should not be taking the square roots. The square roots might make sense if we were calculating an amplitude and there was going to be interference between the spin states but

1. The initial unpolarised beam is a mixed state that is operated on by the polariser and analyser. These mixed states will not interfere with each other. Sasview 9 might consider allowing for misaligned flippers that do not perform a full pi rotation, but that is well outside of our current scope.
2. The weight is being multiplied onto the scattered intensity, not the scattered amplitude, so the treating the weights as amplitudes wouldn't make any sense.

[dirk]

The weights for the polarisation efficiencies of the instrument act on intensities, i.e. remove the square rooted weights from sascalc/calculator/c_extensions/libfunc.c:cal_msld and add in sasview\src\sas\sascalc\calculator\c_extensions\sld2i.c:genicomXY the proper weight for Iout derived as the sum of the spin-resolved scattering cross-sections.
Concerning the weights, make sure that these are polarisation efficiencies, i.e. for the incoming side eff_u=(1+up_i)/2 and eff_d=(1-up_i)/2 and not polarisations up_i or up_f. I assume that will also introduce proper weights such that one can reconstruct from POLARIS cross sections the SANSPOL and unpolarised cross sections with correct absolute intensities.

[awashington]

Fixed for sasmodels in pr 69.

I haven't touched this yet for sasview, because genicomXY is different from the primary scattering code and looks to be calculating an actual amplitude. I'm basing this mostly on the fact that it's calculating scattering from individual nuclei, so I'm guessing that we're on a length scale where interference will be important, so the amplitude will be necessary. Plus the fact that everything si complex and the amplitude is calculated at the end.

In this case, I believe that we do need to keep the square root terms, since the amplitude is squared at the end of the function. I think that we would need a sqrt term, but not perhaps the fourth root that we currently have?

[dirk]

Still there is a unconventional definition of the weights not beeing expressed as polarisation, but in terms of fraction of UP spins before the sample (in_spin or up_i in the code) and after sample (out_spin).
The polarisation P of a neutron beam is the average over the individual neutron polarisation with 0< P<1. One can then introduce the fraction up_i=½(1+P). For an unpolarised beam up_i=50% and P=0, and P=+1 or P=-1 for a fully polarised beam in one or the other of the two spin directions. With an additional analyser behind the sample with polarisation power A, one can distinguish 4 cross-sections.
The two non-spin-flip scattering cross section (i.e. uu and dd, depending on the sign of the incoming polarisation P) are given by
NSF= N² + P.(N* Q + N Q*) + (P.Q)²
and the spin-flip scattering ( i.e. ud and du)
SF=Q.Q-(P.Q)² -i P.(Q x Q*).
The scalar product W.X and W x X the vector product between the two vector quantities W and X, * denotes a complex-conjugated quantity. The vector quantities are correctly calculated in libfunc.c.

For half-polarised SANS (SANSPOL) without spin-state discrimination after the sample (A=0 or out_spin=0.5, I⁺=uu+ud and I^-=dd+du), you end up with the well known result:
I=N² + Q.Q + P.(N* Q + N Q*) -i P.(Q x Q*) and here the polarisation P not up_i appears as weight.

I have the feeling that the final spin-resolved scattering cross-section for non-perfect polariser and analyser is given by
I=1/(1+|A*P|)((1+A*P)*NSF+(1-A*P)*SF).

The weights are hence not up_i and up_f but P=2 up_i -1 and A=2 up_f -1 and work on intensities not amplitudes.

[awashington]

The use of spin fraction as the base value is certainly muddying the waters, but I think that I may still be misunderstanding something. I would have argued that the weights cannot be the polarisation because that would allow for negative weights when there is negative polarisation, which I would have thought to be non-physical.

[butler]

We probably should have a separate conf call to sort out this issue? At the moment it doesn't seem like it can be resolved for 4.2 so going to move it to 4.3 at this point. Should perhaps be noted in release notes for 4.2?

[dirk]

At the end of the day, the weights serve to construct various magnet scattering cross sections, which are linear combinations of the spin-resolved cross sections. The efficiency weights and work on intensities not amplitudes.
The parameters in_spin and out_spin can be easily to polarisation efficiencies (of the instrument) ranging from 0.5 (unpolarised) beam to 1 (perfect optics).
For in_spin or out_spin <0.5 one assumes a cross section, that has the spin reversed/flipped with respect to the initial supermirror polariser. The actual polarisation efficiency in this case is however e_in/out = 1-in/out_spin.
Using in/out_spin (which are not quite right efficiencies), a norm is needed to make sure that the scattering cross sections are correctly weighted with incoming/outgoing flux, such that the sum of spin-resolved measurements adds up to the unpolarised or half-polarised scattering cross section.
For out_spin, the norm is

if (out_spin < 0.5){norm=1-out_spin;}
  else{norm=out_spin;}

and the intensity weights look like

 wuu = (1.0-in_spin) * (1.0-out_spin) / norm;
    wud = (1.0-in_spin) * out_spin / norm
    wdu = in_spin * (1.0-out_spin) / norm
    wdd = in_spin * out_spin / norm
    I = wuu*Iuu + wud*Iud + wdu*Idu + wdd*Idd 

No intensity weighting is needed on the incoming polariser side taking for granted that a user has normalised to the incoming flux with polariser in for SANSPOL and unpolarised beam, respectively.

Tested and works for me.

[dirk]

In 5e1875c4a542fb1755c98100e92d33fb7f9b43dc/sasmodels:

Correct polarisation efficiency weights for magnetic models. Addresses #1086

[dirk]

I done some testing and I think that the code in sascalc/calculator/c_extensions/libfunc.c:cal_msld is obsolete together with sascalc/calculator/c_extensions/sld2i.c:genicomXY. Will be commented out everywhere, let's see if it is not breaking anything.

[dirk]

In b36e7c726c1824fc66e2a5a76abb7351c49d66e6/sasview:

clean up obsolete? code of magnetic scattering calculation. Addresses #1086.