[17bbadd] | 1 | """ |
---|
| 2 | Product model |
---|
| 3 | ------------- |
---|
| 4 | |
---|
| 5 | The product model multiplies the structure factor by the form factor, |
---|
| 6 | modulated by the effective radius of the form. The resulting model |
---|
| 7 | has a attributes of both the model description (with parameters, etc.) |
---|
| 8 | and the module evaluator (with call, release, etc.). |
---|
| 9 | |
---|
| 10 | To use it, first load form factor P and structure factor S, then create |
---|
| 11 | *ProductModel(P, S)*. |
---|
| 12 | """ |
---|
| 13 | import numpy as np |
---|
| 14 | |
---|
| 15 | from .core import call_ER_VR |
---|
| 16 | from .generate import process_parameters |
---|
| 17 | |
---|
| 18 | SCALE=0 |
---|
| 19 | BACKGROUND=1 |
---|
[3bcd03d] | 20 | RADIUS_EFFECTIVE=2 |
---|
[17bbadd] | 21 | VOLFRACTION=3 |
---|
| 22 | |
---|
| 23 | def make_product_info(p_info, s_info): |
---|
| 24 | """ |
---|
| 25 | Create info block for product model. |
---|
| 26 | """ |
---|
| 27 | p_id, p_name, p_pars = p_info['id'], p_info['name'], p_info['parameters'] |
---|
| 28 | s_id, s_name, s_pars = s_info['id'], s_info['name'], s_info['parameters'] |
---|
| 29 | # We require models to start with scale and background |
---|
[fcd7bbd] | 30 | assert s_pars[SCALE].name == 'scale' |
---|
| 31 | assert s_pars[BACKGROUND].name == 'background' |
---|
[17bbadd] | 32 | # We require structure factors to start with effect radius and volfraction |
---|
[3bcd03d] | 33 | assert s_pars[RADIUS_EFFECTIVE].name == 'radius_effective' |
---|
[fcd7bbd] | 34 | assert s_pars[VOLFRACTION].name == 'volfraction' |
---|
[17bbadd] | 35 | # Combine the parameter sets. We are skipping the first three |
---|
| 36 | # parameters of S since scale, background are defined in P and |
---|
| 37 | # effect_radius is set from P.ER(). |
---|
| 38 | pars = p_pars + s_pars[3:] |
---|
| 39 | # check for duplicates; can't use assertion since they may never be checked |
---|
[fcd7bbd] | 40 | if len(set(p.name for p in pars)) != len(pars): |
---|
[17bbadd] | 41 | raise ValueError("Duplicate parameters in %s and %s"%(p_id)) |
---|
| 42 | # For comparison with sasview, determine the old parameters. |
---|
| 43 | oldname = [p_info['oldname'], s_info['oldname']] |
---|
| 44 | oldpars = {'scale':'scale_factor'} |
---|
| 45 | oldpars.update(p_info['oldpars']) |
---|
| 46 | oldpars.update(s_info['oldpars']) |
---|
| 47 | |
---|
| 48 | model_info = {} |
---|
| 49 | model_info['id'] = '*'.join((p_id, s_id)) |
---|
| 50 | model_info['name'] = ' X '.join((p_name, s_name)) |
---|
| 51 | model_info['filename'] = None |
---|
| 52 | model_info['title'] = 'Product of %s and structure factor %s'%(p_name, s_name) |
---|
| 53 | model_info['description'] = model_info['title'] |
---|
[fcd7bbd] | 54 | model_info['docs'] = model_info['title'] |
---|
[17bbadd] | 55 | model_info['category'] = "custom" |
---|
| 56 | model_info['parameters'] = pars |
---|
[fcd7bbd] | 57 | #model_info['single'] = p_info['single'] and s_info['single'] |
---|
| 58 | model_info['structure_factor'] = False |
---|
| 59 | model_info['variant_info'] = None |
---|
| 60 | #model_info['tests'] = [] |
---|
| 61 | #model_info['source'] = [] |
---|
| 62 | # Iq, Iqxy, form_volume, ER, VR and sesans |
---|
[17bbadd] | 63 | model_info['oldname'] = oldname |
---|
| 64 | model_info['oldpars'] = oldpars |
---|
[fcd7bbd] | 65 | model_info['composition'] = ('product', [p_info, s_info]) |
---|
[17bbadd] | 66 | process_parameters(model_info) |
---|
| 67 | return model_info |
---|
| 68 | |
---|
| 69 | class ProductModel(object): |
---|
[72a081d] | 70 | def __init__(self, model_info, P, S): |
---|
| 71 | self.info = model_info |
---|
[17bbadd] | 72 | self.P = P |
---|
| 73 | self.S = S |
---|
| 74 | |
---|
| 75 | def __call__(self, q_vectors): |
---|
| 76 | # Note: may be sending the q_vectors to the GPU twice even though they |
---|
| 77 | # are only needed once. It would mess up modularity quite a bit to |
---|
| 78 | # handle this optimally, especially since there are many cases where |
---|
| 79 | # separate q vectors are needed (e.g., form in python and structure |
---|
| 80 | # in opencl; or both in opencl, but one in single precision and the |
---|
| 81 | # other in double precision). |
---|
| 82 | p_kernel = self.P(q_vectors) |
---|
| 83 | s_kernel = self.S(q_vectors) |
---|
| 84 | return ProductKernel(self.info, p_kernel, s_kernel) |
---|
| 85 | |
---|
| 86 | def release(self): |
---|
| 87 | """ |
---|
| 88 | Free resources associated with the model. |
---|
| 89 | """ |
---|
| 90 | self.P.release() |
---|
| 91 | self.S.release() |
---|
| 92 | |
---|
| 93 | |
---|
| 94 | class ProductKernel(object): |
---|
| 95 | def __init__(self, model_info, p_kernel, s_kernel): |
---|
| 96 | dim = '2d' if p_kernel.q_input.is_2d else '1d' |
---|
| 97 | |
---|
| 98 | # Need to know if we want 2D and magnetic parameters when constructing |
---|
| 99 | # a parameter map. |
---|
| 100 | par_map = {} |
---|
| 101 | p_info = p_kernel.info['partype'] |
---|
| 102 | s_info = s_kernel.info['partype'] |
---|
| 103 | vol_pars = set(p_info['volume']) |
---|
| 104 | if dim == '2d': |
---|
| 105 | num_p_fixed = len(p_info['fixed-2d']) |
---|
| 106 | num_p_pd = len(p_info['pd-2d']) |
---|
| 107 | num_s_fixed = len(s_info['fixed-2d']) |
---|
[35b4c47] | 108 | num_s_pd = len(s_info['pd-2d']) - 1 # exclude effect_radius |
---|
[17bbadd] | 109 | # volume parameters are amongst the pd pars for P, not S |
---|
| 110 | vol_par_idx = [k for k,v in enumerate(p_info['pd-2d']) |
---|
| 111 | if v in vol_pars] |
---|
| 112 | else: |
---|
| 113 | num_p_fixed = len(p_info['fixed-1d']) |
---|
| 114 | num_p_pd = len(p_info['pd-1d']) |
---|
| 115 | num_s_fixed = len(s_info['fixed-1d']) |
---|
[35b4c47] | 116 | num_s_pd = len(s_info['pd-1d']) - 1 # exclude effect_radius |
---|
[17bbadd] | 117 | # volume parameters are amongst the pd pars for P, not S |
---|
| 118 | vol_par_idx = [k for k,v in enumerate(p_info['pd-1d']) |
---|
| 119 | if v in vol_pars] |
---|
| 120 | |
---|
| 121 | start = 0 |
---|
| 122 | par_map['p_fixed'] = np.arange(start, start+num_p_fixed) |
---|
| 123 | # User doesn't set scale, background or effect_radius for S in P*S, |
---|
| 124 | # so borrow values from end of p_fixed. This makes volfraction the |
---|
| 125 | # first S parameter. |
---|
[d5ba841] | 126 | start += num_p_fixed |
---|
| 127 | par_map['s_fixed'] = np.hstack(([start,start], |
---|
| 128 | np.arange(start, start+num_s_fixed-2))) |
---|
[17bbadd] | 129 | par_map['volfraction'] = num_p_fixed |
---|
[d5ba841] | 130 | start += num_s_fixed-2 |
---|
[17bbadd] | 131 | # vol pars offset from the start of pd pars |
---|
| 132 | par_map['vol_pars'] = [start+k for k in vol_par_idx] |
---|
| 133 | par_map['p_pd'] = np.arange(start, start+num_p_pd) |
---|
[d5ba841] | 134 | start += num_p_pd-1 |
---|
| 135 | par_map['s_pd'] = np.hstack((start, |
---|
| 136 | np.arange(start, start+num_s_pd-1))) |
---|
[17bbadd] | 137 | |
---|
| 138 | self.fixed_pars = model_info['partype']['fixed-' + dim] |
---|
| 139 | self.pd_pars = model_info['partype']['pd-' + dim] |
---|
| 140 | self.info = model_info |
---|
| 141 | self.p_kernel = p_kernel |
---|
| 142 | self.s_kernel = s_kernel |
---|
| 143 | self.par_map = par_map |
---|
| 144 | |
---|
| 145 | def __call__(self, fixed_pars, pd_pars, cutoff=1e-5): |
---|
| 146 | pars = fixed_pars + pd_pars |
---|
| 147 | scale = pars[SCALE] |
---|
| 148 | background = pars[BACKGROUND] |
---|
| 149 | s_volfraction = pars[self.par_map['volfraction']] |
---|
| 150 | p_fixed = [pars[k] for k in self.par_map['p_fixed']] |
---|
| 151 | s_fixed = [pars[k] for k in self.par_map['s_fixed']] |
---|
| 152 | p_pd = [pars[k] for k in self.par_map['p_pd']] |
---|
| 153 | s_pd = [pars[k] for k in self.par_map['s_pd']] |
---|
| 154 | vol_pars = [pars[k] for k in self.par_map['vol_pars']] |
---|
| 155 | |
---|
| 156 | effect_radius, vol_ratio = call_ER_VR(self.p_kernel.info, vol_pars) |
---|
| 157 | |
---|
| 158 | p_fixed[SCALE] = s_volfraction |
---|
| 159 | p_fixed[BACKGROUND] = 0.0 |
---|
| 160 | s_fixed[SCALE] = scale |
---|
| 161 | s_fixed[BACKGROUND] = 0.0 |
---|
[35b4c47] | 162 | s_fixed[2] = s_volfraction/vol_ratio |
---|
| 163 | s_pd[0] = [effect_radius], [1.0] |
---|
[17bbadd] | 164 | |
---|
[8e0d974] | 165 | p_res = self.p_kernel(p_fixed, p_pd, cutoff) |
---|
| 166 | s_res = self.s_kernel(s_fixed, s_pd, cutoff) |
---|
[35b4c47] | 167 | #print s_fixed, s_pd, p_fixed, p_pd |
---|
[17bbadd] | 168 | |
---|
| 169 | return p_res*s_res + background |
---|
| 170 | |
---|
| 171 | def release(self): |
---|
| 172 | self.p_kernel.release() |
---|
| 173 | self.q_kernel.release() |
---|
| 174 | |
---|