Extract object/shell volume for combined models (sphere@hardsphere)

[toqduj]

Hi,

Maybe I'm missing something.. I've used sasmodels in McSAS3 via the following call:

F, Fsq, R_eff, V_shell, V_ratio = sasmodels.direct_model.call_Fq(
                self.kernel, dict(self.staticParameters, **parameters)
            )

Now I want to use this for a kernel constructed from a sphere@hardsphere model. I can do this in a jupyter notebook, e.g.:

comboFunc = sasmodels.core.load_model('sphere@hardsphere')
HSSkernel = comboFunc.make_kernel(md['Q'])
HSSFsq = sasmodels.direct_model.call_kernel(
                HSSkernel, dict(comboFunc.info.parameters.defaults, 
                                **{
                                    'radius':radius, 
                                    'radius_effective_mode':1, # 1: radius_effective follows radius
                                    'structure_factor_mode':1, # beta approximation
                                    'volfraction': volFrac,
                                    'background': 0
                                })
            )

However, I'll need the V_shell as well. When I use call_Fq instead of call_kernel, I get a "not implemented" error. Is there another universal way besides call_Fq to extract V_shell which will work both for "simple" as well as compound models?

[pkienzle]

Shell volume and volume ratio are available in HSSkernel.results(): https://github.com/SasView/sasmodels/blob/26813653223cffe3939c791f1ec92a08d556e055/sasmodels/product.py#L287-L306

You will need to test that HSSkernel is a product.ProductKernel first. Normal kernels do not have a results attribute and mixture kernels do not collect Fq intermediates.

[toqduj]

Ah thanks. So for normal operation, I'll just need Fsq or I, and a measure of volume. Based on your information, to make this as universal as possible, I should (please correct if I'm misunderstanding):

1. check if the kernel is instance of product.ProductKernel, if so, use call_kernel and extract the volume using HSSkernel.results()['volume']
2. check if the kernel is instance of mixture.MixedKernel, if so, I can use call_Fq to get Fsq and volume?
3. Otherwise use call_Fq to get Fsq and volume... (as I see that ProductKernel is also an instance of kernel.Kernel so I can't test for that)

So I guess 2 and 3 would be identical then... (I haven't tried mixed kernels yet). Thanks for the help!

[pkienzle]

You cannot get volume data from MixtureKernel with the current implementation. It is not even well defined since in many cases the underlying model (e.g., polymer micelle) does not have fixed volume. Mixture kernels also allow multiplication by a structure factor, which again does not have a volume component.

You could modify the Iq loop, collecting volume data when it is available and putting it into the intermediates object, similar to what is done in the product kernel. I suppose ignore structure factors and use NaN for kernels without volume (polymer, fractal, broad peak, etc.). I'm not sure how to handle subtraction.

https://github.com/SasView/sasmodels/blob/26813653223cffe3939c791f1ec92a08d556e055/sasmodels/mixture.py#L199-L236

[toqduj]

Thanks for the input so far. I'm slowly getting there with allowing ProductKernels in McSAS3. Before that, though, I must delve into the differences between call_Fq and call_kernel outputs, and how they differ from SasView's output:

Leaving mixturekernel aside for the moment (I'll probably need a relative-volume calculator for that if that becomes an issue), I've spent a few hours trying to match the output of call_Fq with call_kernel, since the documentation on calc_kernel does not specify what it outputs. I've found that the output of call_kernel is still volume-weighted so needs another multiplication with the object volume before it's behaving correctly with multiple radii. Fsq out of call_Fq, however, is already scaling with V^2.

One question that I'm struggling with now though, is that I cannot get the output of either of these to match the output when I use the sasview GUI, so there's still a factor I'm missing. It seems I have to multiply the output Fsq from call_Fq, and the output FsqV_shell of call_kernel with 1/(4/3pi) to match the answer SasView is giving me for identical settings (r=1, sld=1, sld_solvent=0, scale=1, background=0). Just not sure where this factor 1/(4/3 pi) originates from, seems to have dropped out of a volume calculation somewhere.

Also good to know that when using a structure factor, that the overall volume fraction of particles will be the structure factor volume fraction * scaling. I assume this is in the documentation somewhere, but I haven't come across it yet.

Test code:

comboFunc = sasmodels.core.load_model('sphere@hardsphere')
sFunc = sasmodels.core.load_model('sphere')
HSSkernel = comboFunc.make_kernel(Q)
Skernel = sFunc.make_kernel(Q)
# remember parameter values
# cumInt = np.zeros(Q[0].shape) # cumulated intensities
# sCumInt = np.zeros(Q[0].shape) # cumulated intensities
HSSVset = np.zeros(radii.shape[0])
SVset = np.zeros(radii.shape[0])

radii = np.array([1, 2, 3])

fh, ah = plt.subplots(figsize=[8, 6])
# call the model for each parameter set explicitly
# otherwise the model gets complex for multiple params incl. fitting
for i, radius in enumerate(radii): # for each contribution
    print(f'{i}: {radius}')
    # result squared or not is model type dependent
    P, Psq, R_eff, V_shell, V_ratio = sasmodels.direct_model.call_Fq(
                    Skernel, dict(sFunc.info.parameters.defaults, **{
                        'radius':radius, 
                        'background': 0,
                        'sld': 1,
                        'sld_solvent': 0,                    
                        'scale': 1,
                    })
                )
    # sCumInt += Psq # * V_shell
    SphereInt = Psq / (4/3*np.pi)
    SVset[i] = V_shell
    # now for the HSS one:
    HSSFsq = sasmodels.direct_model.call_kernel(
                    HSSkernel, dict(comboFunc.info.parameters.defaults, 
                                    **{
                                        'radius':radius, 
                                        # 'radius_effective': 10, 
                                        'radius_effective_mode':1, # 1: radius_effective follows radius
                                        'structure_factor_mode':1, # beta approximation
                                        'volfraction': 0.01,
                                        'sld': 1,
                                        'sld_solvent': 0,
                                        'background': 0,
                                        'scale': 1
                                    })
                )
    V_shell = HSSkernel.results()["volume"]
    HSSVset[i] = V_shell

    # cumInt += HSSFsq * V_shell # vset[i]
    HardSphereInt = HSSFsq * V_shell / (4/3*np.pi)
    ah.plot(Q[0], Psq , color='blue', label = f'sphere, r={radius}')
    ah.plot(Q[0], HSSFsq * V_shell, color='red', label = f'HSS, r={radius}')
    # check the I0 values:
    print(f'r={radius} sphere I(0): {(SphereInt)[0]}, HSS I(0): {(HardSphereInt)[0]}')

plt.yscale('log')
plt.xscale('log')
plt.legend(loc=0)

output:

[pkienzle]

SasView should be using the following: https://github.com/SasView/sasmodels/blob/26813653223cffe3939c791f1ec92a08d556e055/sasmodels/kernel.py#L107-L111 In your code above I don't see you using scale*Psq/V_shell + background. The value of V_shell should be 4πR³/3, so for scale=1 R=1 background=0 that would explain the result.

The β approximation calculation is much more complicated. Docs are here. Feel free to add missing detail.

Sasview applies a resolution function to the returned Iq.