Fractal model not useable with dQ data due to being undefined for negative q

[smk78]

Reported by Rafal Kryszton and verified by smk78.

Affects Mac & Win.

The Fractal model in 4.1 builds and the 4.0.1 release is exhibiting unusual behaviour, not present in the 3.1.2 release, when dQ smearing is enabled.

An error is thrown: Plotting Error: Data has no positive values, and therefore cannot be log-scaled.

Depending on the root cause, this issue may affect other models.

To demonstrate:

Load the attached data. Send to Fitting. Select the Fractal model. Unselect all parameters.

With Instrumental Smearing = None, set fractal_dim to a number other than 2 and hit Enter. The model should recompute.

Now set Instrumental Smearing = Use dQ data and repeat. If fractal_dim is any value other than 2 the error is thrown and the model is not recomputed or disappears from view.

Reselecting Instrumental Smearing = None redisplays the model line.

Side Comment:

The .txt file is the users original data file. This reads into 4.¼.0.1 BUT NOT into 3.1.2. The .csv version reads into all.

[smk78]

Rafal data TXT

[smk78]

Rafal Data CSV

[smk78]

Update.

Following a suggestion on todays developers call, I have just repeated my tests above with the same dataset but with the interweaved data points removed.

It makes no difference.

So the error is coming from the dQ data.

There are no negative values in any of the data columns.

[smk78]

Rafal Data CSV but with interwoven points removed (4 of them)

[krzywon]

The further I look into this, the less I'm convinced this is related to #826. What I've noticed (in case this helps anyone):

(1) Any non-integer value for the fractal dimension when dQ data smearing is selected throws the "Plotting Error: Data has no positive values, and therefore cannot be log-scaled." error.
(2) Fitting the data with an initial integer value for the fractal dimension, the fractal dimension as a fittable parameter, and dQ data set results in a non-converging fit "Fitting did not converge!!!".
(3) Fitting the data as in (2), but with the fractal dimension as a fixed parameter results in a converging fit.
(4) Fitting the data as in (2), but with no dQ data set results in a non-integer value for the fractal dimension.

As far as I can tell, there is nothing forcing the fractal dimension to be an integer in sasmodels, but if this parameter is supposed to be an integer, this could be related to #815.

[ajj]

• The model doesn't work at all when any resolution is added and the fractal dimension is non integer.
• It is not dependent on dQ values from the data : putting in a custom pinhole on the model alone breaks it

[pkienzle]

The resolution calculator is extending into negative q, which for the fractal model, leads to (q*r)^-D ⇒ NaN, which then spreads itself to every point via the resolution matrix, hence no good values.

Options are to take the absolute value of q when computing resolution values, or cut off at some minimum q.

Allowing negative q makes some sort of sense: you are saying that for points near the direct beam, beam divergence will pick up values from the other side of the detector, so some rings will be double counted with different weights.

If we choose to cutoff then the cutoff can't be zero since some models diverge. A value of qmin*0.1 would be practical. We will need this even if we allow negative q since q=0 may be undefined in that case as well.

[richardh]

After discussion with Andrew, we think that the models ought to deal with both Q<0 and Q=0, not the smearing routine. If needs be all the models should do Q = abs(Q) when they start. Thus the fractal model in ticket #835 should be modified.
To avoid issues with the few models that actually diverge or are undefined at exactly Q=0, then the smearing routine could avoid Q=0 when extrapolating. How does the smearing choose values to use when extrapolating? [I'm not clear either whether the smearing only interpolates the theory curve at the user supplied Q points, or whether it calculates extra theory points on a finer grid, as it ought to if it is going to get decent results on curves with sharp dips & peaks. Really we ought to have options for "fast" smearing where there are no really sharp features or "slow but accurate".]

[pkienzle]

Putting the correction in each model means that some models will forget to implement it, and sasview will fail inconsistently over a very long period of time (the problem only shows up when the resolution is unusually large at low q AND the user happens to try the failing model on that data set AND that the particular model parameters happen to give !NaN results in that case. Even worse, some models will silently seem to work even though they give bad values for negative q. The user never sees the calculated I(q) before smearing (even if we add the option to show it, they won't necessarily know that they should be looking at it for their particular data set).

Resolution is the best place to put the correct since it already incorporates intimate details of the geometry. Compare Pinhole1D vs. Slit1D for example. Both of them make assumptions about how the data is being measured when determining how to weight the q values. They can also assume that q values on the opposite side of the beam center have the same circular average as the positive q values, and thus return abs(q).

We should not put the correction in the 1D modelling wrapper—-somebody may be using sasview to fit data which is not SANS, and which has true y for negative x. In this hypothetical case we would need to use a different resolution function appropriate for the data, which for gaussian resolution happens to look like Pinhole1D without the absolute q correction.

The resolution function does a simple extrapolation from the q range to determine the q values to calculate. It uses the spacing between the first pair and last pair as a guide. For pinhole geometry the extrapolation is linear. For slit geometric the extrapolation is geometric. See sasmodels/resolution.py for details.

Sometimes the extrapolated value may hit q=0 or an otherwise very low q value. Proper handling for models which are divergent at low q is tricky, and may, for example, require proper handling of multiple scattering for the specific instrument geometry. Including absurdly low q values in the resolution calculation will do nothing to help this condition. One option is to set a |q| cutoff based on the minimum measured q so that models are well behaved. Another option is to simply reject points whose resolution width goes below 0.1*qmin, with the argument that correctly including them requires too much information about the instrument geometry and the scattering model. Better to not model it all than to pretend that we can model it when we know that the result will be incorrect. Since option 2 is a lot more work, I suggest we do option 1 for 4.1 and save option 2 for 4.2.

Beyond 4.1, the correct approach to resolution would be to check each measured q value to make sure that it has a sufficient basis of support in the calculated q values to get a reasonable estimate of the resolution, oversampling the data throughout the range as necessary (or possibly undersampling if the data density is higher than necessary for the resolution). Selecting the number of points required for support allows us to vary between "fast" and "slow but accurate" as suggested by Richard. There may already be a ticket for this; it will require numerous changes throughout the gui and fitting code.

USANS data for a model with Bragg peaks at high q are problematic: linear sampling over the entirety of the q range is impractical, but the geometric extrapolation may miss the peaks. Choosing the appropriate values would require intimate knowledge of the d-spacing in the model to know which q values should be sampled. Perhaps we need a way for the user to set the q calculated values so they can have this precise control?

[richardh]

So if for large dQ values at small Q the required Q has gone negative in the smearing, we are likely going to be giving the "wrong" answer as we should be say doing full 2d smearing rather than just 1d smearing.

As a simple compromise fix for this could we persuade the resolution routine itself to do Q=abs(Q) before it calls the model function?

[richardh]

Sorry, Paul K had previously said, in email -

Use abs(q) when computing negative q values for resolution calculations.

This still doesn't address the problem where a sample q point happens to lie at q=0 for a model which diverges at q=0. Perhaps adjust any points where |q| < Qmin to q=Qmin, with Qmin = 0.1*min(q) to protect against this?

I think this would be OK, might go smaller still, say 0.02*min(Q) ?
Richard

[pkienzle]

The algorithm now excludes q values in (-0.02 qmin, 0.02 qmin). There may still be points at negative q if the delta q at qmin is large. In the rafal dataset, which extrapolates from the q range (0.001547381, 0.4673276) to the q range (-0.0004811488, 0.58678529), the extrapolation uses the following points:

-0.00048115 -0.00022758 --- 0.00027955 0.00053312 0.00078668 0.00104025 0.00129381

with the point near zero excluded. From this set of q centers, it computes mid points between each and uses the error function to integrate the gaussian value between those midpoints. It then takes the absolute values of the q centers since q at -0.00048115 is the same as q at 0.00048115, but still associates them with the distribution weights for the negative q values.

[pkienzle]

Question: can USANS resolution include negative q values at the lowest q? Currently the algorithm uses logarithmic spacing for extrapolating, but it would be easy enough to use linear extrapolation at low q with:

diff --git a/sasmodels/resolution.py b/sasmodels/resolution.py
index 73cc553..c4592df 100644
--- a/sasmodels/resolution.py
+++ b/sasmodels/resolution.py
@@ -354,9 +354,12 @@ def slit_extend_q(q, width, height):
     that each point I(q) has sufficient support from the underlying
     function.
     """
-    q_min, q_max = np.min(q-height), np.max(np.sqrt((q+height)**2 + width**2))
+    data_min = np.min(q)
+    q_min, q_max = data_min-height, np.max(np.sqrt((q+height)**2 + width**2))
 
-    return geometric_extrapolation(q, q_min, q_max)
+    # geometric extrapolation above, linear extrapolation below
+    q = geometric_extrapolation(q, data_min, q_max)
+    q = linear_extrapolation(q, q_min, q_max)

[GitHub <noreply@…>]

In 1aa7990916fb1aa2b527637492552f3d694bba5b/sasmodels:

Merge pull request #28 from SasView?/ticket-835

ticket-835: allow negative q in resolution calculation. Closes #835