Changeset 975ec8e in sasview for sansmodels/src/sans
- Timestamp:
- Aug 14, 2009 5:58:58 PM (15 years ago)
- Branches:
- master, ESS_GUI, ESS_GUI_Docs, ESS_GUI_batch_fitting, ESS_GUI_bumps_abstraction, ESS_GUI_iss1116, ESS_GUI_iss879, ESS_GUI_iss959, ESS_GUI_opencl, ESS_GUI_ordering, ESS_GUI_sync_sascalc, costrafo411, magnetic_scatt, release-4.1.1, release-4.1.2, release-4.2.2, release_4.0.1, ticket-1009, ticket-1094-headless, ticket-1242-2d-resolution, ticket-1243, ticket-1249, ticket885, unittest-saveload
- Children:
- d5bd424
- Parents:
- cad821b
- Location:
- sansmodels/src/sans/models
- Files:
-
- 4 added
- 36 edited
Legend:
- Unmodified
- Added
- Removed
-
sansmodels/src/sans/models/BinaryHSModel.py
r9bd69098 r975ec8e 32 32 for details of the model. 33 33 List of default parameters: 34 l_radius = 1 60.0 [A]34 l_radius = 100.0 [A] 35 35 s_radius = 25.0 [A] 36 vol_frac_ls = 0. 236 vol_frac_ls = 0.1 37 37 vol_frac_ss = 0.2 38 38 ls_sld = 3.5e-006 [1/A²] … … 53 53 self.name = "BinaryHSModel" 54 54 ## Model description 55 self.description =""" 56 Model parameters: 57 58 l_radius : large radius of the binary hard sphere 59 s_radius : small radius of the binary hard sphere 55 self.description =""" Model parameters: l_radius : large radius of binary hard sphere 56 s_radius : small radius of binary hard sphere 60 57 vol_frac_ls : volume fraction of large spheres 61 58 vol_frac_ss : volume fraction of small spheres -
sansmodels/src/sans/models/EllipticalCylinderModel.py
r9bd69098 r975ec8e 69 69 70 70 ## fittable parameters 71 self.fixed=['cyl_phi.width', 'cyl_theta.width', 'cyl_psi.width', 'length.width', 'r_minor.width', 'r_ratio.width </text>']71 self.fixed=['cyl_phi.width', 'cyl_theta.width', 'cyl_psi.width', 'length.width', 'r_minor.width', 'r_ratio.width'] 72 72 73 73 ## parameters with orientation -
sansmodels/src/sans/models/LamellarModel.py
r9bd69098 r975ec8e 33 33 List of default parameters: 34 34 scale = 1.0 35 delta= 50.0 [A]36 s igma = 0.1537 contrast = 5.3e-006 [1/A²]35 bi_thick = 50.0 [A] 36 sld_bi = 1e-006 [1/A²] 37 sld_sol = 6.3e-006 [1/A²] 38 38 background = 0.0 [1/cm] 39 39 … … 52 52 self.description ="""[Dilute Lamellar Form Factor](from a lyotropic lamellar phase) 53 53 I(q)= 2*pi*P(q)/(delta *q^(2)), where 54 P(q)=2*(contrast/q)^(2)*(1-cos(q*delta) 55 *e^(1/2*(q*sigma)^(2)).56 delta = bilayer thickness57 s igma = variation in bilayer thickness58 = delta*polydispersity59 contrast = SLD_solvent - SLD_bilayer60 Note: the polydispersity in delta is included."""54 P(q)=2*(contrast/q)^(2)*(1-cos(q*delta))^(2)) 55 bi_thick = bilayer thickness 56 sld_bi = SLD of bilayer 57 sld_sol = SLD of solvent 58 background = Incoherent background 59 scale = scale factor 60 """ 61 61 62 62 ## Parameter details [units, min, max] 63 63 self.details = {} 64 64 self.details['scale'] = ['', None, None] 65 self.details[' delta'] = ['[A]', None, None]66 self.details['s igma'] = ['', None, None]67 self.details[' contrast'] = ['[1/A²]', None, None]65 self.details['bi_thick'] = ['[A]', None, None] 66 self.details['sld_bi'] = ['[1/A²]', None, None] 67 self.details['sld_sol'] = ['[1/A²]', None, None] 68 68 self.details['background'] = ['[1/cm]', None, None] 69 69 -
sansmodels/src/sans/models/OblateModel.py
r9bd69098 r975ec8e 90 90 91 91 ## parameters with orientation 92 self.orientation_params =[ ]92 self.orientation_params =['axis_phi', 'axis_theta', 'axis_phi.width', 'axis_theta.width'] 93 93 94 94 def clone(self): -
sansmodels/src/sans/models/ParallelepipedModel.py
r9bd69098 r975ec8e 38 38 contrast = 5.3e-006 [1/A²] 39 39 background = 0.0 [1/cm] 40 parallel_theta = 1.0 [rad] 41 parallel_phi = 1.0 [rad] 40 parallel_theta = 0.0 [rad] 41 parallel_phi = 0.0 [rad] 42 parallel_psi = 0.0 [rad] 42 43 43 44 """ … … 72 73 self.details['parallel_theta'] = ['[rad]', None, None] 73 74 self.details['parallel_phi'] = ['[rad]', None, None] 75 self.details['parallel_psi'] = ['[rad]', None, None] 74 76 75 77 ## fittable parameters 76 self.fixed=['short_edgeA.width', 'longer_edgeB.width', 'longuest_edgeC.width', 'parallel_phi.width', 'parallel_ theta.width']78 self.fixed=['short_edgeA.width', 'longer_edgeB.width', 'longuest_edgeC.width', 'parallel_phi.width', 'parallel_psi.width', 'parallel_theta.width'] 77 79 78 80 ## parameters with orientation 79 self.orientation_params =['parallel_phi', 'parallel_ theta', 'parallel_phi.width', 'parallel_theta.width']81 self.orientation_params =['parallel_phi', 'parallel_psi', 'parallel_theta', 'parallel_phi.width', 'parallel_psi.width', 'parallel_theta.width'] 80 82 81 83 def clone(self): -
sansmodels/src/sans/models/StackedDisksModel.py
r9bd69098 r975ec8e 34 34 scale = 0.01 35 35 radius = 3000.0 [A] 36 length= 10.0 [A]37 thickness= 15.0 [A]36 core_thick = 10.0 [A] 37 layer_thick = 15.0 [A] 38 38 core_sld = 4e-006 [1/A²] 39 39 layer_sld = -4e-007 [1/A²] 40 40 solvent_sld = 5e-006 [1/A²] 41 n layers= 1.042 s pacing= 0.041 n_stacking = 1.0 42 sigma_d = 0.0 43 43 background = 0.001 [1/cm] 44 axis_theta = 1.0 [rad]45 axis_phi = 1.0 [rad]44 axis_theta = 0.0 [rad] 45 axis_phi = 0.0 [rad] 46 46 47 47 """ … … 63 63 self.details['scale'] = ['', None, None] 64 64 self.details['radius'] = ['[A]', None, None] 65 self.details[' length'] = ['[A]', None, None]66 self.details[' thickness'] = ['[A]', None, None]65 self.details['core_thick'] = ['[A]', None, None] 66 self.details['layer_thick'] = ['[A]', None, None] 67 67 self.details['core_sld'] = ['[1/A²]', None, None] 68 68 self.details['layer_sld'] = ['[1/A²]', None, None] 69 69 self.details['solvent_sld'] = ['[1/A²]', None, None] 70 self.details['n layers'] = ['', None, None]71 self.details['s pacing'] = ['', None, None]70 self.details['n_stacking'] = ['', None, None] 71 self.details['sigma_d'] = ['', None, None] 72 72 self.details['background'] = ['[1/cm]', None, None] 73 73 self.details['axis_theta'] = ['[rad]', None, None] … … 75 75 76 76 ## fittable parameters 77 self.fixed=[' length.width', 'radius.width', 'axis_theta.width', 'axis_phi.width']77 self.fixed=['core_thick.width', 'layer_thick.width', 'radius.width', 'axis_theta.width', 'axis_phi.width'] 78 78 79 79 ## parameters with orientation -
sansmodels/src/sans/models/TriaxialEllipsoidModel.py
r9bd69098 r975ec8e 38 38 contrast = 5.3e-006 [1/A²] 39 39 background = 0.0 [1/cm] 40 axis_theta = 1.0 [rad] 41 axis_phi = 1.0 [rad] 40 axis_theta = 0.0 [rad] 41 axis_phi = 0.0 [rad] 42 axis_psi = 0.0 [rad] 42 43 43 44 """ … … 67 68 self.details['axis_theta'] = ['[rad]', None, None] 68 69 self.details['axis_phi'] = ['[rad]', None, None] 70 self.details['axis_psi'] = ['[rad]', None, None] 69 71 70 72 ## fittable parameters 71 self.fixed=['axis_p hi.width', 'axis_theta.width', 'semi_axisA.width', 'semi_axisB.width', 'semi_axisC.width']73 self.fixed=['axis_psi.width', 'axis_phi.width', 'axis_theta.width', 'semi_axisA.width', 'semi_axisB.width', 'semi_axisC.width'] 72 74 73 75 ## parameters with orientation 74 self.orientation_params =['axis_p hi', 'axis_theta', 'axis_phi.width', 'axis_theta.width']76 self.orientation_params =['axis_psi', 'axis_phi', 'axis_theta', 'axis_psi.width', 'axis_phi.width', 'axis_theta.width'] 75 77 76 78 def clone(self): -
sansmodels/src/sans/models/c_extensions/binaryHs.h
r3d25331f r975ec8e 7 7 [PYTHONCLASS] = BinaryHSModel 8 8 [DISP_PARAMS] = l_radius,s_radius 9 [DESCRIPTION] =<text> 10 Model parameters: 11 12 l_radius : large radius of the binary hard sphere 13 s_radius : small radius of the binary hard sphere 14 vol_frac_ls : volume fraction of large spheres 15 vol_frac_ss : volume fraction of small spheres 16 ls_sld: large sphere scattering length density 17 ss_sld: small sphere scattering length density 18 solvent_sld: solvent scattering length density 19 background: incoherent background 9 [DESCRIPTION] =<text> Model parameters: l_radius : large radius of binary hard sphere 10 s_radius : small radius of binary hard sphere 11 vol_frac_ls : volume fraction of large spheres 12 vol_frac_ss : volume fraction of small spheres 13 ls_sld: large sphere scattering length density 14 ss_sld: small sphere scattering length density 15 solvent_sld: solvent scattering length density 16 background: incoherent background 20 17 </text> 21 18 [FIXED]= l_radius.width;s_radius.width … … 23 20 */ 24 21 typedef struct { 25 22 26 23 /// large radius of the binary hard sphere [A] 27 // [DEFAULT]=l_radius= 1 60.0 [A]24 // [DEFAULT]=l_radius= 100.0 [A] 28 25 double l_radius; 29 26 … … 32 29 double s_radius; 33 30 34 /// volume fraction of large spheres 35 // [DEFAULT]=vol_frac_ls= 0. 231 /// volume fraction of large spheres 32 // [DEFAULT]=vol_frac_ls= 0.1 36 33 double vol_frac_ls; 37 34 38 /// volume fraction of small spheres 39 // [DEFAULT]=vol_frac_ss= 0.2 35 /// volume fraction of small spheres 36 // [DEFAULT]=vol_frac_ss= 0.2 40 37 double vol_frac_ss; 41 38 -
sansmodels/src/sans/models/c_extensions/elliptical_cylinder.c
r3fe701a r975ec8e 31 31 } 32 32 33 double elliptical_cylinder_kernel(EllipticalCylinderParameters *pars, double q, double alpha, double psi, doublenu) {33 double elliptical_cylinder_kernel(EllipticalCylinderParameters *pars, double q, double alpha, double nu) { 34 34 double qr; 35 35 double qL; 36 double Be,Si; 36 37 double r_major; 37 38 double kernel; … … 42 43 qL = q*pars->length*cos(alpha)/2.0; 43 44 44 kernel = 2.0*NR_BessJ1(qr)/qr * sin(qL)/qL; 45 if (qr==0){ 46 Be = 0.5; 47 }else{ 48 Be = NR_BessJ1(qr)/qr; 49 } 50 if (qL==0){ 51 Si = 1.0; 52 }else{ 53 Si = sin(qL)/qL; 54 } 55 56 57 kernel = 2.0*Be * Si; 45 58 return kernel*kernel; 46 59 } … … 129 142 } 130 143 131 answer = elliptical_cylinder_kernel(pars, q, alpha, pars->cyl_psi,nu);144 answer = elliptical_cylinder_kernel(pars, q, alpha,nu); 132 145 133 146 // Multiply by contrast^2 -
sansmodels/src/sans/models/c_extensions/lamellar.c
r42f193a r975ec8e 22 22 // Fill paramater array 23 23 dp[0] = pars->scale; 24 dp[1] = pars-> delta;25 dp[2] = pars->s igma;26 dp[3] = pars-> contrast;24 dp[1] = pars->bi_thick; 25 dp[2] = pars->sld_bi; 26 dp[3] = pars->sld_sol; 27 27 dp[4] = pars->background; 28 28 29 29 30 30 // Call library function to evaluate model 31 return LamellarFF(dp, q);31 return lamellar_kernel(dp, q); 32 32 } 33 34 33 35 /** 34 36 * Function to evaluate 2D scattering function -
sansmodels/src/sans/models/c_extensions/lamellar.h
r42f193a r975ec8e 3 3 /** Structure definition for lamellar parameters 4 4 * [PYTHONCLASS] = LamellarModel 5 * [DISP_PARAMS] = bi_thick 5 6 [DESCRIPTION] = <text>[Dilute Lamellar Form Factor](from a lyotropic lamellar phase) 6 7 P(q)=2*(contrast/q)^(2)*(1-cos(q*delta) 8 *e^(1/2*(q*sigma)^(2)).9 delta = bilayer thickness10 s igma = variation in bilayer thickness11 = delta*polydispersity12 contrast = SLD_solvent - SLD_bilayer13 Note: the polydispersity in delta is included. 7 I(q)= 2*pi*P(q)/(delta *q^(2)), where 8 P(q)=2*(contrast/q)^(2)*(1-cos(q*delta))^(2)) 9 bi_thick = bilayer thickness 10 sld_bi = SLD of bilayer 11 sld_sol = SLD of solvent 12 background = Incoherent background 13 scale = scale factor 14 14 15 </text> 16 [FIXED]= <text>bi_thick.width</text> 15 17 **/ 16 18 typedef struct { … … 19 21 double scale; 20 22 /// delta bilayer thickness [A] 21 // [DEFAULT]= delta=50.0 [A]22 double delta;23 /// variation in bilayer thickness24 // [DEFAULT]=s igma=0.1525 double s igma;26 /// Contrast [1/A²]27 // [DEFAULT]= contrast=5.3e-6 [1/A²]28 double contrast;23 // [DEFAULT]=bi_thick=50.0 [A] 24 double bi_thick; 25 /// SLD of bilayer [1/A²] 26 // [DEFAULT]=sld_bi=1.0e-6 [1/A²] 27 double sld_bi; 28 /// SLD of solvent [1/A²] 29 // [DEFAULT]=sld_sol=6.3e-6 [1/A²] 30 double sld_sol; 29 31 /// Incoherent Background [1/cm] 0.00 30 32 // [DEFAULT]=background=0.0 [1/cm] -
sansmodels/src/sans/models/c_extensions/oblate.c
r96b59384 r975ec8e 66 66 * @return: function value 67 67 */ 68 /* 68 69 69 double oblate_analytical_2D_scaled(OblateParameters *pars, double q, double q_x, double q_y) { 70 70 71 71 return 1.0; 72 72 } 73 */ 73 -
sansmodels/src/sans/models/c_extensions/oblate.h
r96b59384 r975ec8e 24 24 25 25 [FIXED] = <text>major_core.width;minor_core.width; major_shell.width; minor_shell.width</text> 26 [ORIENTATION_PARAMS] =26 [ORIENTATION_PARAMS]= <text>axis_phi; axis_theta; axis_phi.width; axis_theta.width</text> 27 27 28 28 **/ … … 52 52 // [DEFAULT]=background=0.001 [1/cm] 53 53 double background; 54 / *//Disable for now54 //Disable for now 55 55 /// Orientation of the oblate axis w/respect incoming beam [rad] 56 56 // [DEFAULT]=axis_theta=1.0 [rad] … … 59 59 // [DEFAULT]=axis_phi=1.0 [rad] 60 60 double axis_phi; 61 */ 61 62 62 } OblateParameters; 63 63 -
sansmodels/src/sans/models/c_extensions/parallelepiped.c
r5068697 r975ec8e 19 19 double parallelepiped_analytical_1D(ParallelepipedParameters *pars, double q) { 20 20 double dp[6]; 21 21 22 22 // Fill paramater array 23 23 dp[0] = pars->scale; … … 29 29 30 30 // Call library function to evaluate model 31 return Parallelepiped(dp, q); 31 return Parallelepiped(dp, q); 32 32 } 33 34 35 double pkernel(double a, double b,double c, double ala, double alb, double alc){ 36 // mu passed in is really mu*sqrt(1-sig^2) 37 double argA,argB,argC,tmp1,tmp2,tmp3; //local variables 38 39 //handle arg=0 separately, as sin(t)/t -> 1 as t->0 40 argA = a*ala; 41 argB = b*alb; 42 argC = c*alc; 43 if(argA==0.0) { 44 tmp1 = 1.0; 45 } else { 46 tmp1 = sin(argA)*sin(argA)/argA/argA; 47 } 48 if (argB==0.0) { 49 tmp2 = 1.0; 50 } else { 51 tmp2 = sin(argB)*sin(argB)/argB/argB; 52 } 53 54 if (argC==0.0) { 55 tmp3 = 1.0; 56 } else { 57 tmp3 = sin(argC)*sin(argC)/argC/argC; 58 } 59 60 return (tmp1*tmp2*tmp3); 61 62 }//Function pkernel() 63 64 65 66 33 67 /** 34 68 * Function to evaluate 2D scattering function … … 41 75 q = sqrt(qx*qx+qy*qy); 42 76 return parallelepiped_analytical_2D_scaled(pars, q, qx/q, qy/q); 43 } 77 } 44 78 45 79 … … 53 87 double parallelepiped_analytical_2D(ParallelepipedParameters *pars, double q, double phi) { 54 88 return parallelepiped_analytical_2D_scaled(pars, q, cos(phi), sin(phi)); 55 } 56 89 } 90 57 91 /** 58 92 * Function to evaluate 2D scattering function … … 64 98 */ 65 99 double parallelepiped_analytical_2D_scaled(ParallelepipedParameters *pars, double q, double q_x, double q_y) { 66 double parallel_x, parallel_y, parallel_z;100 double cparallel_x, cparallel_y, cparallel_z, bparallel_x, bparallel_y, parallel_x, parallel_y, parallel_z; 67 101 double q_z; 68 double alpha, vol, cos_val ;69 double aa, mu, uu; 102 double alpha, vol, cos_val_c, cos_val_b, cos_val_a, edgeA, edgeB, edgeC; 103 70 104 double answer; 71 72 105 double pi = 4.0*atan(1.0); 73 106 74 // parallelepiped orientation 75 parallel_x = sin(pars->parallel_theta) * cos(pars->parallel_phi); 76 parallel_y = sin(pars->parallel_theta) * sin(pars->parallel_phi); 77 parallel_z = cos(pars->parallel_theta); 78 107 edgeA = pars->short_edgeA; 108 edgeB = pars->longer_edgeB; 109 edgeC = pars->longuest_edgeC; 110 111 112 // parallelepiped c axis orientation 113 cparallel_x = sin(pars->parallel_theta) * cos(pars->parallel_phi); 114 cparallel_y = sin(pars->parallel_theta) * sin(pars->parallel_phi); 115 cparallel_z = cos(pars->parallel_theta); 116 79 117 // q vector 80 118 q_z = 0; 81 119 82 120 // Compute the angle btw vector q and the 83 121 // axis of the parallelepiped 84 cos_val = parallel_x*q_x + parallel_y*q_y + parallel_z*q_z; 85 122 cos_val_c = cparallel_x*q_x + cparallel_y*q_y + cparallel_z*q_z; 123 alpha = acos(cos_val_c); 124 125 // parallelepiped a axis orientation 126 parallel_x = (1-sin(pars->parallel_theta)*sin(pars->parallel_phi))*sin(pars->parallel_psi);//cos(pars->parallel_theta) * sin(pars->parallel_phi)*sin(pars->parallel_psi); 127 parallel_y = (1-sin(pars->parallel_theta)*sin(pars->parallel_phi))*cos(pars->parallel_psi);//cos(pars->parallel_theta) * cos(pars->parallel_phi)*cos(pars->parallel_psi); 128 cos_val_a = (parallel_x*q_x + parallel_y*q_y); 129 130 131 132 // parallelepiped b axis orientation 133 bparallel_x = (1-sin(pars->parallel_theta)*cos(pars->parallel_phi))*cos(pars->parallel_psi);//cos(pars->parallel_theta) * cos(pars->parallel_phi)* cos(pars->parallel_psi); 134 bparallel_y = (1-sin(pars->parallel_theta)*cos(pars->parallel_phi))*sin(pars->parallel_psi);//cos(pars->parallel_theta) * sin(pars->parallel_phi)* sin(pars->parallel_psi); 135 // axis of the parallelepiped 136 cos_val_b = (bparallel_x*q_x + bparallel_y*q_y) ; 137 138 139 86 140 // The following test should always pass 87 if (fabs(cos_val )>1.0) {141 if (fabs(cos_val_c)>1.0) { 88 142 printf("parallel_ana_2D: Unexpected error: cos(alpha)>1\n"); 89 143 return 0; 90 144 } 91 92 // Note: cos(alpha) = 0 and 1 will get an93 // undefined value from PPKernel94 alpha = acos( cos_val );95 145 96 aa = pars->short_edgeA/pars->longer_edgeB;97 mu = 1.0;98 uu = 1.0;99 100 146 // Call the IGOR library function to get the kernel 101 answer = PPKernel( aa, mu, uu);102 147 answer = pkernel( q*edgeA, q*edgeB, q*edgeC, cos_val_a,cos_val_b,cos_val_c); 148 103 149 // Multiply by contrast^2 104 150 answer *= pars->contrast*pars->contrast; 105 151 106 152 //normalize by cylinder volume 107 153 //NOTE that for this (Fournet) definition of the integral, one must MULTIPLY by Vparallel 108 vol = pars->short_edgeA* pars->longer_edgeB * pars->longuest_edgeC;154 vol = 8*edgeA* edgeB * edgeC; 109 155 answer *= vol; 110 156 111 157 //convert to [cm-1] 112 158 answer *= 1.0e8; 113 159 114 160 //Scale 115 161 answer *= pars->scale; 116 162 117 163 // add in the background 118 164 answer += pars->background; 119 165 120 166 return answer; 121 167 } 122 123 168 -
sansmodels/src/sans/models/c_extensions/parallelepiped.h
r5068697 r975ec8e 7 7 /** Structure definition for Parallelepiped parameters 8 8 * [PYTHONCLASS] = ParallelepipedModel 9 * [DISP_PARAMS] = short_edgeA, longer_edgeB, longuest_edgeC,parallel_phi, parallel_theta9 * [DISP_PARAMS] = short_edgeA, longer_edgeB, longuest_edgeC,parallel_phi,parallel_psi, parallel_theta 10 10 [DESCRIPTION] = <text> Calculates the form factor for a rectangular solid with uniform scattering length density. 11 11 12 12 scale:Scale factor 13 13 short_edgeA: Shortest edge of the parallelepiped [A] … … 15 15 longuest_edgeC: Longuest edge of the parallelepiped [A] 16 16 constrast: particle_sld - solvent_sld 17 background:Incoherent Background [1/cm] 17 background:Incoherent Background [1/cm] 18 18 </text> 19 [FIXED]= <text>short_edgeA.width; longer_edgeB.width; longuest_edgeC.width;parallel_phi.width; parallel_theta.width</text>20 [ORIENTATION_PARAMS]= <text>parallel_phi; parallel_theta; parallel_phi.width; parallel_theta.width</text>19 [FIXED]= <text>short_edgeA.width; longer_edgeB.width; longuest_edgeC.width;parallel_phi.width;parallel_psi.width; parallel_theta.width</text> 20 [ORIENTATION_PARAMS]= <text>parallel_phi;parallel_psi; parallel_theta; parallel_phi.width;parallel_psi.width; parallel_theta.width</text> 21 21 22 22 … … 38 38 // [DEFAULT]=contrast=5.3e-6 [1/A²] 39 39 double contrast; 40 /// Incoherent Background [1/cm] 40 /// Incoherent Background [1/cm] 41 41 // [DEFAULT]=background=0.0 [1/cm] 42 42 double background; 43 43 /// Orientation of the parallelepiped axis w/respect incoming beam [rad] 44 // [DEFAULT]=parallel_theta= 1.0 [rad]44 // [DEFAULT]=parallel_theta=0.0 [rad] 45 45 double parallel_theta; 46 /// Orientation of the parallelepiped in the plane of the detector [rad]47 // [DEFAULT]=parallel_phi= 1.0 [rad]46 /// Orientation of the longitudinal axis of the parallelepiped in the plane of the detector [rad] 47 // [DEFAULT]=parallel_phi=0.0 [rad] 48 48 double parallel_phi; 49 /// Orientation of the cross-sectional minor axis of the parallelepiped in the plane of the detector [rad] 50 // [DEFAULT]=parallel_psi=0.0 [rad] 51 double parallel_psi; 49 52 50 53 … … 60 63 double parallelepiped_analytical_2DXY(ParallelepipedParameters *pars, double qx, double qy); 61 64 double parallelepiped_analytical_2D_scaled(ParallelepipedParameters *pars, double q, double q_x, double q_y); 62 63 65 #endif -
sansmodels/src/sans/models/c_extensions/stacked_disks.h
r5068697 r975ec8e 7 7 /** Structure definition for stacked disks parameters 8 8 * [PYTHONCLASS] = StackedDisksModel 9 * [DISP_PARAMS] = length, radius, axis_theta, axis_phi9 * [DISP_PARAMS] = core_thick, layer_thick, radius, axis_theta, axis_phi 10 10 [DESCRIPTION] = <text> 11 11 </text> 12 [FIXED]= <text> length.width; radius.width; axis_theta.width; axis_phi.width</text>12 [FIXED]= <text>core_thick.width;layer_thick.width; radius.width; axis_theta.width; axis_phi.width</text> 13 13 [ORIENTATION_PARAMS]= <text>axis_phi; axis_theta; axis_phi.width; axis_theta.width</text> 14 14 … … 22 22 // [DEFAULT]=radius=3000 [A] 23 23 double radius; 24 /// Length of the stakeddisk [A]25 // [DEFAULT]= length=10.0 [A]26 double length;24 /// Thickness of the core disk [A] 25 // [DEFAULT]=core_thick=10.0 [A] 26 double core_thick; 27 27 /// Thickness of the staked disk [A] 28 // [DEFAULT]= thickness=15.0 [A]29 double thickness;28 // [DEFAULT]=layer_thick=15.0 [A] 29 double layer_thick; 30 30 /// Core scattering length density[1/A²] 31 31 // [DEFAULT]=core_sld=4e-6 [1/A²] … … 37 37 // [DEFAULT]=solvent_sld=5.0e-6 [1/A²] 38 38 double solvent_sld; 39 /// number of layers40 // [DEFAULT]=n layers=141 double n layers;42 /// GSD of disks s pacing43 // [DEFAULT]=s pacing=044 double s pacing;45 /// Incoherent Background [1/cm] 39 /// number of stacking 40 // [DEFAULT]=n_stacking=1 41 double n_stacking; 42 /// GSD of disks sigma_d 43 // [DEFAULT]=sigma_d=0 44 double sigma_d; 45 /// Incoherent Background [1/cm] 46 46 // [DEFAULT]=background=0.001 [1/cm] 47 47 double background; 48 48 /// Orientation of the staked disk axis w/respect incoming beam [rad] 49 // [DEFAULT]=axis_theta= 1.0 [rad]49 // [DEFAULT]=axis_theta=0.0 [rad] 50 50 double axis_theta; 51 51 /// Orientation of the staked disk in the plane of the detector [rad] 52 // [DEFAULT]=axis_phi= 1.0 [rad]52 // [DEFAULT]=axis_phi=0.0 [rad] 53 53 double axis_phi; 54 54 -
sansmodels/src/sans/models/c_extensions/triaxial_ellipsoid.c
r34c3020 r975ec8e 19 19 double triaxial_ellipsoid_analytical_1D(TriaxialEllipsoidParameters *pars, double q) { 20 20 double dp[6]; 21 21 22 22 // Fill paramater array 23 23 dp[0] = pars->scale; … … 27 27 dp[4] = pars->contrast; 28 28 dp[5] = pars->background; 29 29 30 30 // Call library function to evaluate model 31 return TriaxialEllipsoid(dp, q); 31 return TriaxialEllipsoid(dp, q); 32 32 } 33 34 double triaxial_ellipsoid_kernel(TriaxialEllipsoidParameters *pars, double q, double alpha, double nu) { 35 double t,a,b,c; 36 double kernel; 37 double pi = acos(-1.0); 38 39 a = pars->semi_axisA ; 40 b = pars->semi_axisB ; 41 c = pars->semi_axisC ; 42 43 t = q * sqrt(a*a*cos(nu)*cos(nu)+b*b*sin(nu)*sin(nu)*sin(alpha)*sin(alpha)+c*c*cos(alpha)*cos(alpha)); 44 if (t==0){ 45 kernel = 1.0; 46 }else{ 47 kernel = 3*(sin(t)-t*cos(t))/(t*t*t); 48 } 49 return kernel*kernel; 50 } 51 33 52 34 53 /** … … 42 61 q = sqrt(qx*qx+qy*qy); 43 62 return triaxial_ellipsoid_analytical_2D_scaled(pars, q, qx/q, qy/q); 44 } 63 } 45 64 46 65 … … 54 73 double triaxial_ellipsoid_analytical_2D(TriaxialEllipsoidParameters *pars, double q, double phi) { 55 74 return triaxial_ellipsoid_analytical_2D_scaled(pars, q, cos(phi), sin(phi)); 56 } 57 75 } 76 58 77 /** 59 78 * Function to evaluate 2D scattering function … … 65 84 */ 66 85 double triaxial_ellipsoid_analytical_2D_scaled(TriaxialEllipsoidParameters *pars, double q, double q_x, double q_y) { 67 double cyl_x, cyl_y, cyl_z ;86 double cyl_x, cyl_y, cyl_z, ell_x, ell_y; 68 87 double q_z; 69 double dx, dy;88 double cos_nu,nu; 70 89 double alpha, vol, cos_val; 71 90 double answer; … … 75 94 cyl_y = sin(pars->axis_theta) * sin(pars->axis_phi); 76 95 cyl_z = cos(pars->axis_theta); 77 96 78 97 // q vector 79 98 q_z = 0; 80 81 dx = 1.0;82 dy = 1.0;99 100 //dx = 1.0; 101 //dy = 1.0; 83 102 // Compute the angle btw vector q and the 84 103 // axis of the cylinder 85 104 cos_val = cyl_x*q_x + cyl_y*q_y + cyl_z*q_z; 86 105 87 106 // The following test should always pass 88 107 if (fabs(cos_val)>1.0) { … … 90 109 return 0; 91 110 } 92 111 93 112 // Note: cos(alpha) = 0 and 1 will get an 94 113 // undefined value from CylKernel 95 114 alpha = acos( cos_val ); 96 115 116 //ellipse orientation: 117 // the elliptical corss section was transformed and projected 118 // into the detector plane already through sin(alpha)and furthermore psi remains as same 119 // on the detector plane. 120 // So, all we need is to calculate the angle (nu) of the minor axis of the ellipse wrt 121 // the wave vector q. 122 123 //x- y- component on the detector plane. 124 ell_x = cos(pars->axis_psi); 125 ell_y = sin(pars->axis_psi); 126 127 // calculate the axis of the ellipse wrt q-coord. 128 cos_nu = ell_x*q_x + ell_y*q_y; 129 nu = acos(cos_nu); 130 97 131 // Call the IGOR library function to get the kernel 98 answer = TriaxialKernel(q,pars->semi_axisA, pars->semi_axisB, pars->semi_axisC, dx, dy);99 132 answer = triaxial_ellipsoid_kernel(pars, q, alpha, nu); 133 100 134 // Multiply by contrast^2 101 135 answer *= pars->contrast*pars->contrast; 102 136 103 137 //normalize by cylinder volume 104 138 //NOTE that for this (Fournet) definition of the integral, one must MULTIPLY by Vcyl 105 139 vol = 4/3 * pi * pars->semi_axisA * pars->semi_axisB * pars->semi_axisC; 106 140 answer *= vol; 107 141 108 142 //convert to [cm-1] 109 143 answer *= 1.0e8; 110 144 111 145 //Scale 112 146 answer *= pars->scale; 113 147 114 148 // add in the background 115 149 answer += pars->background; 116 150 117 151 return answer; 118 152 } 119 -
sansmodels/src/sans/models/c_extensions/triaxial_ellipsoid.h
r5068697 r975ec8e 3 3 /** Structure definition for cylinder parameters 4 4 * [PYTHONCLASS] = TriaxialEllipsoidModel 5 * [DISP_PARAMS] = axis_theta, axis_phi 5 * [DISP_PARAMS] = axis_theta, axis_phi, axis_psi 6 6 [DESCRIPTION] = <text> Note: 7 7 Constraints must be applied during fitting to ensure that the inequality a<b<c is not 8 8 violated. The calculation will not report an error, but the results will not be correct. 9 9 </text> 10 [FIXED]= <text>axis_p hi.width; axis_theta.width; semi_axisA.width; semi_axisB.width; semi_axisC.width </text>11 [ORIENTATION_PARAMS]= <text>axis_p hi; axis_theta; axis_phi.width; axis_theta.width</text>10 [FIXED]= <text>axis_psi.width; axis_phi.width; axis_theta.width; semi_axisA.width; semi_axisB.width; semi_axisC.width </text> 11 [ORIENTATION_PARAMS]= <text>axis_psi; axis_phi; axis_theta; axis_psi.width; axis_phi.width; axis_theta.width</text> 12 12 13 13 … … 33 33 double background; 34 34 /// Orientation of the triaxial_ellipsoid axis w/respect incoming beam [rad] 35 // [DEFAULT]=axis_theta= 1.0 [rad]35 // [DEFAULT]=axis_theta=0.0 [rad] 36 36 double axis_theta; 37 37 /// Orientation of the triaxial_ellipsoid in the plane of the detector [rad] 38 // [DEFAULT]=axis_phi= 1.0 [rad]38 // [DEFAULT]=axis_phi=0.0 [rad] 39 39 double axis_phi; 40 /// Orientation of the cross section of the triaxial_ellipsoid in the plane of the detector [rad] 41 // [DEFAULT]=axis_psi=0.0 [rad] 42 double axis_psi; 40 43 41 44 } TriaxialEllipsoidParameters; -
sansmodels/src/sans/models/c_models/CBinaryHSModel.cpp
r9bd69098 r975ec8e 86 86 87 87 // Initialize parameter dictionary 88 PyDict_SetItemString(self->params,"vol_frac_ls",Py_BuildValue("d",0. 200000));88 PyDict_SetItemString(self->params,"vol_frac_ls",Py_BuildValue("d",0.100000)); 89 89 PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.001000)); 90 90 PyDict_SetItemString(self->params,"vol_frac_ss",Py_BuildValue("d",0.200000)); … … 93 93 PyDict_SetItemString(self->params,"ss_sld",Py_BuildValue("d",0.000000)); 94 94 PyDict_SetItemString(self->params,"s_radius",Py_BuildValue("d",25.000000)); 95 PyDict_SetItemString(self->params,"l_radius",Py_BuildValue("d",1 60.000000));95 PyDict_SetItemString(self->params,"l_radius",Py_BuildValue("d",100.000000)); 96 96 // Initialize dispersion / averaging parameter dict 97 97 DispersionVisitor* visitor = new DispersionVisitor(); … … 170 170 return PyArray_Return(result); 171 171 } 172 /** 173 * Function to call to evaluate model 174 * @param args: input numpy array [q[],phi[]] 175 * @return: numpy array object 176 */ 177 static PyObject * evaluateTwoDim( BinaryHSModel* model, 178 PyArrayObject *q, PyArrayObject *phi) 179 { 180 PyArrayObject *result; 181 //check validity of input vectors 182 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 183 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 184 || phi->dimensions[0] != q->dimensions[0]){ 185 186 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 187 PyErr_SetString(PyExc_ValueError ,"wrong input"); 188 return NULL; 189 } 190 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 191 192 if (result == NULL){ 193 const char * message= "Could not create result "; 194 PyErr_SetString(PyExc_RuntimeError , message); 195 return NULL; 196 } 197 198 for (int i = 0; i < q->dimensions[0]; i++) { 199 double q_value = *(double *)(q->data + i*q->strides[0]); 200 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 201 double *result_value = (double *)(result->data + i*result->strides[0]); 202 if (q_value == 0) 203 *result_value = 0.0; 204 else 205 *result_value = model->evaluate_rphi(q_value, phi_value); 206 } 207 return PyArray_Return(result); 208 } 172 209 173 /** 210 174 * Function to call to evaluate model -
sansmodels/src/sans/models/c_models/CCylinderModel.cpp
r9bd69098 r975ec8e 175 175 return PyArray_Return(result); 176 176 } 177 /** 178 * Function to call to evaluate model 179 * @param args: input numpy array [q[],phi[]] 180 * @return: numpy array object 181 */ 182 static PyObject * evaluateTwoDim( CylinderModel* model, 183 PyArrayObject *q, PyArrayObject *phi) 184 { 185 PyArrayObject *result; 186 //check validity of input vectors 187 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 188 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 189 || phi->dimensions[0] != q->dimensions[0]){ 190 191 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 192 PyErr_SetString(PyExc_ValueError ,"wrong input"); 193 return NULL; 194 } 195 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 196 197 if (result == NULL){ 198 const char * message= "Could not create result "; 199 PyErr_SetString(PyExc_RuntimeError , message); 200 return NULL; 201 } 202 203 for (int i = 0; i < q->dimensions[0]; i++) { 204 double q_value = *(double *)(q->data + i*q->strides[0]); 205 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 206 double *result_value = (double *)(result->data + i*result->strides[0]); 207 if (q_value == 0) 208 *result_value = 0.0; 209 else 210 *result_value = model->evaluate_rphi(q_value, phi_value); 211 } 212 return PyArray_Return(result); 213 } 177 214 178 /** 215 179 * Function to call to evaluate model -
sansmodels/src/sans/models/c_models/CEllipsoidModel.cpp
r9bd69098 r975ec8e 175 175 return PyArray_Return(result); 176 176 } 177 /** 178 * Function to call to evaluate model 179 * @param args: input numpy array [q[],phi[]] 180 * @return: numpy array object 181 */ 182 static PyObject * evaluateTwoDim( EllipsoidModel* model, 183 PyArrayObject *q, PyArrayObject *phi) 184 { 185 PyArrayObject *result; 186 //check validity of input vectors 187 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 188 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 189 || phi->dimensions[0] != q->dimensions[0]){ 190 191 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 192 PyErr_SetString(PyExc_ValueError ,"wrong input"); 193 return NULL; 194 } 195 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 196 197 if (result == NULL){ 198 const char * message= "Could not create result "; 199 PyErr_SetString(PyExc_RuntimeError , message); 200 return NULL; 201 } 202 203 for (int i = 0; i < q->dimensions[0]; i++) { 204 double q_value = *(double *)(q->data + i*q->strides[0]); 205 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 206 double *result_value = (double *)(result->data + i*result->strides[0]); 207 if (q_value == 0) 208 *result_value = 0.0; 209 else 210 *result_value = model->evaluate_rphi(q_value, phi_value); 211 } 212 return PyArray_Return(result); 213 } 177 214 178 /** 215 179 * Function to call to evaluate model -
sansmodels/src/sans/models/c_models/CEllipticalCylinderModel.cpp
r9bd69098 r975ec8e 183 183 return PyArray_Return(result); 184 184 } 185 /** 186 * Function to call to evaluate model 187 * @param args: input numpy array [q[],phi[]] 188 * @return: numpy array object 189 */ 190 static PyObject * evaluateTwoDim( EllipticalCylinderModel* model, 191 PyArrayObject *q, PyArrayObject *phi) 192 { 193 PyArrayObject *result; 194 //check validity of input vectors 195 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 196 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 197 || phi->dimensions[0] != q->dimensions[0]){ 198 199 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 200 PyErr_SetString(PyExc_ValueError ,"wrong input"); 201 return NULL; 202 } 203 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 204 205 if (result == NULL){ 206 const char * message= "Could not create result "; 207 PyErr_SetString(PyExc_RuntimeError , message); 208 return NULL; 209 } 210 211 for (int i = 0; i < q->dimensions[0]; i++) { 212 double q_value = *(double *)(q->data + i*q->strides[0]); 213 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 214 double *result_value = (double *)(result->data + i*result->strides[0]); 215 if (q_value == 0) 216 *result_value = 0.0; 217 else 218 *result_value = model->evaluate_rphi(q_value, phi_value); 219 } 220 return PyArray_Return(result); 221 } 185 222 186 /** 223 187 * Function to call to evaluate model -
sansmodels/src/sans/models/c_models/CLamellarModel.cpp
r9bd69098 r975ec8e 86 86 87 87 // Initialize parameter dictionary 88 PyDict_SetItemString(self->params,"sld_sol",Py_BuildValue("d",0.000006)); 88 89 PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000)); 89 PyDict_SetItemString(self->params," sigma",Py_BuildValue("d",0.150000));90 PyDict_SetItemString(self->params,"bi_thick",Py_BuildValue("d",50.000000)); 90 91 PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000)); 91 PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d",0.000005)); 92 PyDict_SetItemString(self->params,"delta",Py_BuildValue("d",50.000000)); 92 PyDict_SetItemString(self->params,"sld_bi",Py_BuildValue("d",0.000001)); 93 93 // Initialize dispersion / averaging parameter dict 94 94 DispersionVisitor* visitor = new DispersionVisitor(); 95 95 PyObject * disp_dict; 96 disp_dict = PyDict_New(); 97 self->model->bi_thick.dispersion->accept_as_source(visitor, self->model->bi_thick.dispersion, disp_dict); 98 PyDict_SetItemString(self->dispersion, "bi_thick", disp_dict); 96 99 97 100 … … 161 164 return PyArray_Return(result); 162 165 } 163 /** 164 * Function to call to evaluate model 165 * @param args: input numpy array [q[],phi[]] 166 * @return: numpy array object 167 */ 168 static PyObject * evaluateTwoDim( LamellarModel* model, 169 PyArrayObject *q, PyArrayObject *phi) 170 { 171 PyArrayObject *result; 172 //check validity of input vectors 173 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 174 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 175 || phi->dimensions[0] != q->dimensions[0]){ 176 177 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 178 PyErr_SetString(PyExc_ValueError ,"wrong input"); 179 return NULL; 180 } 181 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 182 183 if (result == NULL){ 184 const char * message= "Could not create result "; 185 PyErr_SetString(PyExc_RuntimeError , message); 186 return NULL; 187 } 188 189 for (int i = 0; i < q->dimensions[0]; i++) { 190 double q_value = *(double *)(q->data + i*q->strides[0]); 191 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 192 double *result_value = (double *)(result->data + i*result->strides[0]); 193 if (q_value == 0) 194 *result_value = 0.0; 195 else 196 *result_value = model->evaluate_rphi(q_value, phi_value); 197 } 198 return PyArray_Return(result); 199 } 166 200 167 /** 201 168 * Function to call to evaluate model … … 260 227 261 228 // Reader parameter dictionary 229 self->model->sld_sol = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sld_sol") ); 262 230 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 263 self->model-> sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );231 self->model->bi_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "bi_thick") ); 264 232 self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 265 self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") ); 266 self->model->delta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "delta") ); 233 self->model->sld_bi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sld_bi") ); 267 234 // Read in dispersion parameters 268 235 PyObject* disp_dict; 269 236 DispersionVisitor* visitor = new DispersionVisitor(); 237 disp_dict = PyDict_GetItemString(self->dispersion, "bi_thick"); 238 self->model->bi_thick.dispersion->accept_as_destination(visitor, self->model->bi_thick.dispersion, disp_dict); 270 239 271 240 … … 330 299 331 300 // Reader parameter dictionary 301 self->model->sld_sol = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sld_sol") ); 332 302 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 333 self->model-> sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );303 self->model->bi_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "bi_thick") ); 334 304 self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 335 self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") ); 336 self->model->delta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "delta") ); 305 self->model->sld_bi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sld_bi") ); 337 306 // Read in dispersion parameters 338 307 PyObject* disp_dict; 339 308 DispersionVisitor* visitor = new DispersionVisitor(); 309 disp_dict = PyDict_GetItemString(self->dispersion, "bi_thick"); 310 self->model->bi_thick.dispersion->accept_as_destination(visitor, self->model->bi_thick.dispersion, disp_dict); 340 311 341 312 … … 389 360 390 361 // Reader parameter dictionary 362 self->model->sld_sol = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sld_sol") ); 391 363 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 392 self->model-> sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") );364 self->model->bi_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "bi_thick") ); 393 365 self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 394 self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") ); 395 self->model->delta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "delta") ); 366 self->model->sld_bi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sld_bi") ); 396 367 // Read in dispersion parameters 397 368 PyObject* disp_dict; 398 369 DispersionVisitor* visitor = new DispersionVisitor(); 370 disp_dict = PyDict_GetItemString(self->dispersion, "bi_thick"); 371 self->model->bi_thick.dispersion->accept_as_destination(visitor, self->model->bi_thick.dispersion, disp_dict); 399 372 400 373 … … 452 425 // Ugliness necessary to go from python to C 453 426 // TODO: refactor this 454 { 427 if (!strcmp(par_name, "bi_thick")) { 428 self->model->bi_thick.dispersion = dispersion; 429 } else { 455 430 PyErr_SetString(CLamellarModelError, 456 431 "CLamellarModel.set_dispersion expects a valid parameter name."); -
sansmodels/src/sans/models/c_models/COblateModel.cpp
r9bd69098 r975ec8e 178 178 return PyArray_Return(result); 179 179 } 180 /** 181 * Function to call to evaluate model 182 * @param args: input numpy array [q[],phi[]] 183 * @return: numpy array object 184 */ 185 static PyObject * evaluateTwoDim( OblateModel* model, 186 PyArrayObject *q, PyArrayObject *phi) 187 { 188 PyArrayObject *result; 189 //check validity of input vectors 190 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 191 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 192 || phi->dimensions[0] != q->dimensions[0]){ 193 194 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 195 PyErr_SetString(PyExc_ValueError ,"wrong input"); 196 return NULL; 197 } 198 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 199 200 if (result == NULL){ 201 const char * message= "Could not create result "; 202 PyErr_SetString(PyExc_RuntimeError , message); 203 return NULL; 204 } 205 206 for (int i = 0; i < q->dimensions[0]; i++) { 207 double q_value = *(double *)(q->data + i*q->strides[0]); 208 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 209 double *result_value = (double *)(result->data + i*result->strides[0]); 210 if (q_value == 0) 211 *result_value = 0.0; 212 else 213 *result_value = model->evaluate_rphi(q_value, phi_value); 214 } 215 return PyArray_Return(result); 216 } 180 217 181 /** 218 182 * Function to call to evaluate model -
sansmodels/src/sans/models/c_models/CParallelepipedModel.cpp
r9bd69098 r975ec8e 88 88 PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000)); 89 89 PyDict_SetItemString(self->params,"longer_edgeB",Py_BuildValue("d",75.000000)); 90 PyDict_SetItemString(self->params,"parallel_psi",Py_BuildValue("d",0.000000)); 90 91 PyDict_SetItemString(self->params,"longuest_edgeC",Py_BuildValue("d",400.000000)); 91 PyDict_SetItemString(self->params,"parallel_phi",Py_BuildValue("d", 1.000000));92 PyDict_SetItemString(self->params,"parallel_theta",Py_BuildValue("d", 1.000000));92 PyDict_SetItemString(self->params,"parallel_phi",Py_BuildValue("d",0.000000)); 93 PyDict_SetItemString(self->params,"parallel_theta",Py_BuildValue("d",0.000000)); 93 94 PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000)); 94 95 PyDict_SetItemString(self->params,"short_edgeA",Py_BuildValue("d",35.000000)); … … 109 110 self->model->parallel_phi.dispersion->accept_as_source(visitor, self->model->parallel_phi.dispersion, disp_dict); 110 111 PyDict_SetItemString(self->dispersion, "parallel_phi", disp_dict); 112 disp_dict = PyDict_New(); 113 self->model->parallel_psi.dispersion->accept_as_source(visitor, self->model->parallel_psi.dispersion, disp_dict); 114 PyDict_SetItemString(self->dispersion, "parallel_psi", disp_dict); 111 115 disp_dict = PyDict_New(); 112 116 self->model->parallel_theta.dispersion->accept_as_source(visitor, self->model->parallel_theta.dispersion, disp_dict); … … 179 183 return PyArray_Return(result); 180 184 } 181 /** 182 * Function to call to evaluate model 183 * @param args: input numpy array [q[],phi[]] 184 * @return: numpy array object 185 */ 186 static PyObject * evaluateTwoDim( ParallelepipedModel* model, 187 PyArrayObject *q, PyArrayObject *phi) 188 { 189 PyArrayObject *result; 190 //check validity of input vectors 191 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 192 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 193 || phi->dimensions[0] != q->dimensions[0]){ 194 195 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 196 PyErr_SetString(PyExc_ValueError ,"wrong input"); 197 return NULL; 198 } 199 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 200 201 if (result == NULL){ 202 const char * message= "Could not create result "; 203 PyErr_SetString(PyExc_RuntimeError , message); 204 return NULL; 205 } 206 207 for (int i = 0; i < q->dimensions[0]; i++) { 208 double q_value = *(double *)(q->data + i*q->strides[0]); 209 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 210 double *result_value = (double *)(result->data + i*result->strides[0]); 211 if (q_value == 0) 212 *result_value = 0.0; 213 else 214 *result_value = model->evaluate_rphi(q_value, phi_value); 215 } 216 return PyArray_Return(result); 217 } 185 218 186 /** 219 187 * Function to call to evaluate model … … 280 248 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 281 249 self->model->longer_edgeB = PyFloat_AsDouble( PyDict_GetItemString(self->params, "longer_edgeB") ); 250 self->model->parallel_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "parallel_psi") ); 282 251 self->model->longuest_edgeC = PyFloat_AsDouble( PyDict_GetItemString(self->params, "longuest_edgeC") ); 283 252 self->model->parallel_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "parallel_phi") ); … … 297 266 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_phi"); 298 267 self->model->parallel_phi.dispersion->accept_as_destination(visitor, self->model->parallel_phi.dispersion, disp_dict); 268 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_psi"); 269 self->model->parallel_psi.dispersion->accept_as_destination(visitor, self->model->parallel_psi.dispersion, disp_dict); 299 270 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_theta"); 300 271 self->model->parallel_theta.dispersion->accept_as_destination(visitor, self->model->parallel_theta.dispersion, disp_dict); … … 363 334 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 364 335 self->model->longer_edgeB = PyFloat_AsDouble( PyDict_GetItemString(self->params, "longer_edgeB") ); 336 self->model->parallel_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "parallel_psi") ); 365 337 self->model->longuest_edgeC = PyFloat_AsDouble( PyDict_GetItemString(self->params, "longuest_edgeC") ); 366 338 self->model->parallel_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "parallel_phi") ); … … 380 352 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_phi"); 381 353 self->model->parallel_phi.dispersion->accept_as_destination(visitor, self->model->parallel_phi.dispersion, disp_dict); 354 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_psi"); 355 self->model->parallel_psi.dispersion->accept_as_destination(visitor, self->model->parallel_psi.dispersion, disp_dict); 382 356 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_theta"); 383 357 self->model->parallel_theta.dispersion->accept_as_destination(visitor, self->model->parallel_theta.dispersion, disp_dict); … … 435 409 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 436 410 self->model->longer_edgeB = PyFloat_AsDouble( PyDict_GetItemString(self->params, "longer_edgeB") ); 411 self->model->parallel_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "parallel_psi") ); 437 412 self->model->longuest_edgeC = PyFloat_AsDouble( PyDict_GetItemString(self->params, "longuest_edgeC") ); 438 413 self->model->parallel_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "parallel_phi") ); … … 452 427 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_phi"); 453 428 self->model->parallel_phi.dispersion->accept_as_destination(visitor, self->model->parallel_phi.dispersion, disp_dict); 429 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_psi"); 430 self->model->parallel_psi.dispersion->accept_as_destination(visitor, self->model->parallel_psi.dispersion, disp_dict); 454 431 disp_dict = PyDict_GetItemString(self->dispersion, "parallel_theta"); 455 432 self->model->parallel_theta.dispersion->accept_as_destination(visitor, self->model->parallel_theta.dispersion, disp_dict); … … 517 494 } else if (!strcmp(par_name, "parallel_phi")) { 518 495 self->model->parallel_phi.dispersion = dispersion; 496 } else if (!strcmp(par_name, "parallel_psi")) { 497 self->model->parallel_psi.dispersion = dispersion; 519 498 } else if (!strcmp(par_name, "parallel_theta")) { 520 499 self->model->parallel_theta.dispersion = dispersion; -
sansmodels/src/sans/models/c_models/CSphereModel.cpp
r9bd69098 r975ec8e 163 163 return PyArray_Return(result); 164 164 } 165 /** 166 * Function to call to evaluate model 167 * @param args: input numpy array [q[],phi[]] 168 * @return: numpy array object 169 */ 170 static PyObject * evaluateTwoDim( SphereModel* model, 171 PyArrayObject *q, PyArrayObject *phi) 172 { 173 PyArrayObject *result; 174 //check validity of input vectors 175 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 176 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 177 || phi->dimensions[0] != q->dimensions[0]){ 178 179 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 180 PyErr_SetString(PyExc_ValueError ,"wrong input"); 181 return NULL; 182 } 183 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 184 185 if (result == NULL){ 186 const char * message= "Could not create result "; 187 PyErr_SetString(PyExc_RuntimeError , message); 188 return NULL; 189 } 190 191 for (int i = 0; i < q->dimensions[0]; i++) { 192 double q_value = *(double *)(q->data + i*q->strides[0]); 193 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 194 double *result_value = (double *)(result->data + i*result->strides[0]); 195 if (q_value == 0) 196 *result_value = 0.0; 197 else 198 *result_value = model->evaluate_rphi(q_value, phi_value); 199 } 200 return PyArray_Return(result); 201 } 165 202 166 /** 203 167 * Function to call to evaluate model -
sansmodels/src/sans/models/c_models/CStackedDisksModel.cpp
r9bd69098 r975ec8e 86 86 87 87 // Initialize parameter dictionary 88 PyDict_SetItemString(self->params,"core_sld",Py_BuildValue("d",0.000004)); 89 PyDict_SetItemString(self->params,"core_thick",Py_BuildValue("d",10.000000)); 90 PyDict_SetItemString(self->params,"layer_thick",Py_BuildValue("d",15.000000)); 91 PyDict_SetItemString(self->params,"axis_theta",Py_BuildValue("d",0.000000)); 92 PyDict_SetItemString(self->params,"layer_sld",Py_BuildValue("d",-0.000000)); 93 PyDict_SetItemString(self->params,"axis_phi",Py_BuildValue("d",0.000000)); 94 PyDict_SetItemString(self->params,"solvent_sld",Py_BuildValue("d",0.000005)); 88 95 PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",0.010000)); 89 PyDict_SetItemString(self->params,"axis_theta",Py_BuildValue("d",1.000000));90 PyDict_SetItemString(self->params,"spacing",Py_BuildValue("d",0.000000));91 PyDict_SetItemString(self->params,"nlayers",Py_BuildValue("d",1.000000));92 PyDict_SetItemString(self->params,"layer_sld",Py_BuildValue("d",-0.000000));93 PyDict_SetItemString(self->params,"solvent_sld",Py_BuildValue("d",0.000005));94 PyDict_SetItemString(self->params,"thickness",Py_BuildValue("d",15.000000));95 PyDict_SetItemString(self->params,"axis_phi",Py_BuildValue("d",1.000000));96 PyDict_SetItemString(self->params,"length",Py_BuildValue("d",10.000000));97 PyDict_SetItemString(self->params,"core_sld",Py_BuildValue("d",0.000004));98 96 PyDict_SetItemString(self->params,"radius",Py_BuildValue("d",3000.000000)); 99 97 PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.001000)); 98 PyDict_SetItemString(self->params,"sigma_d",Py_BuildValue("d",0.000000)); 99 PyDict_SetItemString(self->params,"n_stacking",Py_BuildValue("d",1.000000)); 100 100 // Initialize dispersion / averaging parameter dict 101 101 DispersionVisitor* visitor = new DispersionVisitor(); 102 102 PyObject * disp_dict; 103 103 disp_dict = PyDict_New(); 104 self->model->length.dispersion->accept_as_source(visitor, self->model->length.dispersion, disp_dict); 105 PyDict_SetItemString(self->dispersion, "length", disp_dict); 104 self->model->core_thick.dispersion->accept_as_source(visitor, self->model->core_thick.dispersion, disp_dict); 105 PyDict_SetItemString(self->dispersion, "core_thick", disp_dict); 106 disp_dict = PyDict_New(); 107 self->model->layer_thick.dispersion->accept_as_source(visitor, self->model->layer_thick.dispersion, disp_dict); 108 PyDict_SetItemString(self->dispersion, "layer_thick", disp_dict); 106 109 disp_dict = PyDict_New(); 107 110 self->model->radius.dispersion->accept_as_source(visitor, self->model->radius.dispersion, disp_dict); … … 180 183 return PyArray_Return(result); 181 184 } 182 /** 183 * Function to call to evaluate model 184 * @param args: input numpy array [q[],phi[]] 185 * @return: numpy array object 186 */ 187 static PyObject * evaluateTwoDim( StackedDisksModel* model, 188 PyArrayObject *q, PyArrayObject *phi) 189 { 190 PyArrayObject *result; 191 //check validity of input vectors 192 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 193 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 194 || phi->dimensions[0] != q->dimensions[0]){ 195 196 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 197 PyErr_SetString(PyExc_ValueError ,"wrong input"); 198 return NULL; 199 } 200 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 201 202 if (result == NULL){ 203 const char * message= "Could not create result "; 204 PyErr_SetString(PyExc_RuntimeError , message); 205 return NULL; 206 } 207 208 for (int i = 0; i < q->dimensions[0]; i++) { 209 double q_value = *(double *)(q->data + i*q->strides[0]); 210 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 211 double *result_value = (double *)(result->data + i*result->strides[0]); 212 if (q_value == 0) 213 *result_value = 0.0; 214 else 215 *result_value = model->evaluate_rphi(q_value, phi_value); 216 } 217 return PyArray_Return(result); 218 } 185 219 186 /** 220 187 * Function to call to evaluate model … … 279 246 280 247 // Reader parameter dictionary 248 self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") ); 249 self->model->core_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_thick") ); 250 self->model->layer_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_thick") ); 251 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") ); 252 self->model->layer_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_sld") ); 253 self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") ); 254 self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") ); 281 255 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 282 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") );283 self->model->spacing = PyFloat_AsDouble( PyDict_GetItemString(self->params, "spacing") );284 self->model->nlayers = PyFloat_AsDouble( PyDict_GetItemString(self->params, "nlayers") );285 self->model->layer_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_sld") );286 self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") );287 self->model->thickness = PyFloat_AsDouble( PyDict_GetItemString(self->params, "thickness") );288 self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") );289 self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );290 self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") );291 256 self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") ); 292 257 self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 258 self->model->sigma_d = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma_d") ); 259 self->model->n_stacking = PyFloat_AsDouble( PyDict_GetItemString(self->params, "n_stacking") ); 293 260 // Read in dispersion parameters 294 261 PyObject* disp_dict; 295 262 DispersionVisitor* visitor = new DispersionVisitor(); 296 disp_dict = PyDict_GetItemString(self->dispersion, "length"); 297 self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict); 263 disp_dict = PyDict_GetItemString(self->dispersion, "core_thick"); 264 self->model->core_thick.dispersion->accept_as_destination(visitor, self->model->core_thick.dispersion, disp_dict); 265 disp_dict = PyDict_GetItemString(self->dispersion, "layer_thick"); 266 self->model->layer_thick.dispersion->accept_as_destination(visitor, self->model->layer_thick.dispersion, disp_dict); 298 267 disp_dict = PyDict_GetItemString(self->dispersion, "radius"); 299 268 self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict); … … 364 333 365 334 // Reader parameter dictionary 335 self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") ); 336 self->model->core_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_thick") ); 337 self->model->layer_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_thick") ); 338 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") ); 339 self->model->layer_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_sld") ); 340 self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") ); 341 self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") ); 366 342 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 367 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") );368 self->model->spacing = PyFloat_AsDouble( PyDict_GetItemString(self->params, "spacing") );369 self->model->nlayers = PyFloat_AsDouble( PyDict_GetItemString(self->params, "nlayers") );370 self->model->layer_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_sld") );371 self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") );372 self->model->thickness = PyFloat_AsDouble( PyDict_GetItemString(self->params, "thickness") );373 self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") );374 self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );375 self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") );376 343 self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") ); 377 344 self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 345 self->model->sigma_d = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma_d") ); 346 self->model->n_stacking = PyFloat_AsDouble( PyDict_GetItemString(self->params, "n_stacking") ); 378 347 // Read in dispersion parameters 379 348 PyObject* disp_dict; 380 349 DispersionVisitor* visitor = new DispersionVisitor(); 381 disp_dict = PyDict_GetItemString(self->dispersion, "length"); 382 self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict); 350 disp_dict = PyDict_GetItemString(self->dispersion, "core_thick"); 351 self->model->core_thick.dispersion->accept_as_destination(visitor, self->model->core_thick.dispersion, disp_dict); 352 disp_dict = PyDict_GetItemString(self->dispersion, "layer_thick"); 353 self->model->layer_thick.dispersion->accept_as_destination(visitor, self->model->layer_thick.dispersion, disp_dict); 383 354 disp_dict = PyDict_GetItemString(self->dispersion, "radius"); 384 355 self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict); … … 438 409 439 410 // Reader parameter dictionary 411 self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") ); 412 self->model->core_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_thick") ); 413 self->model->layer_thick = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_thick") ); 414 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") ); 415 self->model->layer_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_sld") ); 416 self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") ); 417 self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") ); 440 418 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 441 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") );442 self->model->spacing = PyFloat_AsDouble( PyDict_GetItemString(self->params, "spacing") );443 self->model->nlayers = PyFloat_AsDouble( PyDict_GetItemString(self->params, "nlayers") );444 self->model->layer_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "layer_sld") );445 self->model->solvent_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "solvent_sld") );446 self->model->thickness = PyFloat_AsDouble( PyDict_GetItemString(self->params, "thickness") );447 self->model->axis_phi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_phi") );448 self->model->length = PyFloat_AsDouble( PyDict_GetItemString(self->params, "length") );449 self->model->core_sld = PyFloat_AsDouble( PyDict_GetItemString(self->params, "core_sld") );450 419 self->model->radius = PyFloat_AsDouble( PyDict_GetItemString(self->params, "radius") ); 451 420 self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 421 self->model->sigma_d = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma_d") ); 422 self->model->n_stacking = PyFloat_AsDouble( PyDict_GetItemString(self->params, "n_stacking") ); 452 423 // Read in dispersion parameters 453 424 PyObject* disp_dict; 454 425 DispersionVisitor* visitor = new DispersionVisitor(); 455 disp_dict = PyDict_GetItemString(self->dispersion, "length"); 456 self->model->length.dispersion->accept_as_destination(visitor, self->model->length.dispersion, disp_dict); 426 disp_dict = PyDict_GetItemString(self->dispersion, "core_thick"); 427 self->model->core_thick.dispersion->accept_as_destination(visitor, self->model->core_thick.dispersion, disp_dict); 428 disp_dict = PyDict_GetItemString(self->dispersion, "layer_thick"); 429 self->model->layer_thick.dispersion->accept_as_destination(visitor, self->model->layer_thick.dispersion, disp_dict); 457 430 disp_dict = PyDict_GetItemString(self->dispersion, "radius"); 458 431 self->model->radius.dispersion->accept_as_destination(visitor, self->model->radius.dispersion, disp_dict); … … 516 489 // Ugliness necessary to go from python to C 517 490 // TODO: refactor this 518 if (!strcmp(par_name, "length")) { 519 self->model->length.dispersion = dispersion; 491 if (!strcmp(par_name, "core_thick")) { 492 self->model->core_thick.dispersion = dispersion; 493 } else if (!strcmp(par_name, "layer_thick")) { 494 self->model->layer_thick.dispersion = dispersion; 520 495 } else if (!strcmp(par_name, "radius")) { 521 496 self->model->radius.dispersion = dispersion; -
sansmodels/src/sans/models/c_models/CTriaxialEllipsoidModel.cpp
r9bd69098 r975ec8e 87 87 // Initialize parameter dictionary 88 88 PyDict_SetItemString(self->params,"scale",Py_BuildValue("d",1.000000)); 89 PyDict_SetItemString(self->params,"axis_theta",Py_BuildValue("d",1.000000)); 89 PyDict_SetItemString(self->params,"axis_psi",Py_BuildValue("d",0.000000)); 90 PyDict_SetItemString(self->params,"axis_theta",Py_BuildValue("d",0.000000)); 90 91 PyDict_SetItemString(self->params,"semi_axisA",Py_BuildValue("d",100.000000)); 91 92 PyDict_SetItemString(self->params,"semi_axisB",Py_BuildValue("d",35.000000)); 92 93 PyDict_SetItemString(self->params,"semi_axisC",Py_BuildValue("d",400.000000)); 93 PyDict_SetItemString(self->params,"axis_phi",Py_BuildValue("d", 1.000000));94 PyDict_SetItemString(self->params,"axis_phi",Py_BuildValue("d",0.000000)); 94 95 PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000)); 95 96 PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d",0.000005)); … … 103 104 self->model->axis_phi.dispersion->accept_as_source(visitor, self->model->axis_phi.dispersion, disp_dict); 104 105 PyDict_SetItemString(self->dispersion, "axis_phi", disp_dict); 106 disp_dict = PyDict_New(); 107 self->model->axis_psi.dispersion->accept_as_source(visitor, self->model->axis_psi.dispersion, disp_dict); 108 PyDict_SetItemString(self->dispersion, "axis_psi", disp_dict); 105 109 106 110 … … 170 174 return PyArray_Return(result); 171 175 } 172 /** 173 * Function to call to evaluate model 174 * @param args: input numpy array [q[],phi[]] 175 * @return: numpy array object 176 */ 177 static PyObject * evaluateTwoDim( TriaxialEllipsoidModel* model, 178 PyArrayObject *q, PyArrayObject *phi) 179 { 180 PyArrayObject *result; 181 //check validity of input vectors 182 if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 183 || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 184 || phi->dimensions[0] != q->dimensions[0]){ 185 186 //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 187 PyErr_SetString(PyExc_ValueError ,"wrong input"); 188 return NULL; 189 } 190 result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 191 192 if (result == NULL){ 193 const char * message= "Could not create result "; 194 PyErr_SetString(PyExc_RuntimeError , message); 195 return NULL; 196 } 197 198 for (int i = 0; i < q->dimensions[0]; i++) { 199 double q_value = *(double *)(q->data + i*q->strides[0]); 200 double phi_value = *(double *)(phi->data + i*phi->strides[0]); 201 double *result_value = (double *)(result->data + i*result->strides[0]); 202 if (q_value == 0) 203 *result_value = 0.0; 204 else 205 *result_value = model->evaluate_rphi(q_value, phi_value); 206 } 207 return PyArray_Return(result); 208 } 176 209 177 /** 210 178 * Function to call to evaluate model … … 270 238 // Reader parameter dictionary 271 239 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 240 self->model->axis_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_psi") ); 272 241 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") ); 273 242 self->model->semi_axisA = PyFloat_AsDouble( PyDict_GetItemString(self->params, "semi_axisA") ); … … 284 253 disp_dict = PyDict_GetItemString(self->dispersion, "axis_phi"); 285 254 self->model->axis_phi.dispersion->accept_as_destination(visitor, self->model->axis_phi.dispersion, disp_dict); 255 disp_dict = PyDict_GetItemString(self->dispersion, "axis_psi"); 256 self->model->axis_psi.dispersion->accept_as_destination(visitor, self->model->axis_psi.dispersion, disp_dict); 286 257 287 258 … … 347 318 // Reader parameter dictionary 348 319 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 320 self->model->axis_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_psi") ); 349 321 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") ); 350 322 self->model->semi_axisA = PyFloat_AsDouble( PyDict_GetItemString(self->params, "semi_axisA") ); … … 361 333 disp_dict = PyDict_GetItemString(self->dispersion, "axis_phi"); 362 334 self->model->axis_phi.dispersion->accept_as_destination(visitor, self->model->axis_phi.dispersion, disp_dict); 335 disp_dict = PyDict_GetItemString(self->dispersion, "axis_psi"); 336 self->model->axis_psi.dispersion->accept_as_destination(visitor, self->model->axis_psi.dispersion, disp_dict); 363 337 364 338 … … 413 387 // Reader parameter dictionary 414 388 self->model->scale = PyFloat_AsDouble( PyDict_GetItemString(self->params, "scale") ); 389 self->model->axis_psi = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_psi") ); 415 390 self->model->axis_theta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "axis_theta") ); 416 391 self->model->semi_axisA = PyFloat_AsDouble( PyDict_GetItemString(self->params, "semi_axisA") ); … … 427 402 disp_dict = PyDict_GetItemString(self->dispersion, "axis_phi"); 428 403 self->model->axis_phi.dispersion->accept_as_destination(visitor, self->model->axis_phi.dispersion, disp_dict); 404 disp_dict = PyDict_GetItemString(self->dispersion, "axis_psi"); 405 self->model->axis_psi.dispersion->accept_as_destination(visitor, self->model->axis_psi.dispersion, disp_dict); 429 406 430 407 … … 486 463 } else if (!strcmp(par_name, "axis_phi")) { 487 464 self->model->axis_phi.dispersion = dispersion; 465 } else if (!strcmp(par_name, "axis_psi")) { 466 self->model->axis_psi.dispersion = dispersion; 488 467 } else { 489 468 PyErr_SetString(CTriaxialEllipsoidModelError, -
sansmodels/src/sans/models/c_models/c_models.cpp
reba9885 r975ec8e 30 30 void addCLamellarPSModel(PyObject *module); 31 31 void addCLamellarPSHGModel(PyObject *module); 32 void addCCoreShellSpheroidModel(PyObject *module); 32 33 void addCOblateModel(PyObject *module); 33 34 void addCProlateModel(PyObject *module); … … 225 226 addCLamellarPSModel(m); 226 227 addCLamellarPSHGModel(m); 228 addCCoreShellSpheroidModel(m); 227 229 addCOblateModel(m); 228 230 addCProlateModel(m); -
sansmodels/src/sans/models/c_models/lamellar.cpp
r42f193a r975ec8e 34 34 LamellarModel :: LamellarModel() { 35 35 scale = Parameter(1.0); 36 delta = Parameter(50.0); 37 delta.set_min(0.0); 38 sigma = Parameter(0.15); 39 sigma.set_min(0.0); 40 contrast = Parameter(5.3e-6); 36 bi_thick = Parameter(50.0, true); 37 bi_thick.set_min(0.0); 38 sld_bi = Parameter(1.0e-6); 39 sld_sol = Parameter(6.3e-6); 41 40 background = Parameter(0.0); 42 41 … … 55 54 // Add the background after averaging 56 55 dp[0] = scale(); 57 dp[1] = delta();58 dp[2] = s igma();59 dp[3] = contrast();60 dp[4] = background();56 dp[1] = bi_thick(); 57 dp[2] = sld_bi(); 58 dp[3] = sld_sol(); 59 dp[4] = 0.0; 61 60 62 return LamellarFF(dp, q); 61 // Get the dispersion points for the bi_thick 62 vector<WeightPoint> weights_bi_thick; 63 bi_thick.get_weights(weights_bi_thick); 64 // Perform the computation, with all weight points 65 double sum = 0.0; 66 double norm = 0.0; 67 68 // Loop over short_edgeA weight points 69 for(int i=0; i< (int)weights_bi_thick.size(); i++) { 70 dp[1] = weights_bi_thick[i].value; 71 72 sum += weights_bi_thick[i].weight * lamellar_kernel(dp, q); 73 norm += weights_bi_thick[i].weight; 74 75 } 76 77 return sum/norm + background(); 63 78 } 64 79 -
sansmodels/src/sans/models/c_models/models.hh
r2c4b289 r975ec8e 46 46 class ParallelepipedModel{ 47 47 public: 48 // TODO: add 2D 48 // TODO: add 2D 49 49 // Model parameters 50 50 Parameter scale; … … 56 56 Parameter parallel_theta; 57 57 Parameter parallel_phi; 58 Parameter parallel_psi; 58 59 59 60 // Constructor … … 277 278 Parameter axis_theta; 278 279 Parameter axis_phi; 279 280 Parameter axis_psi; 281 280 282 // Constructor 281 283 TriaxialEllipsoidModel(); … … 298 300 Parameter axis_theta; 299 301 Parameter axis_phi; 300 302 301 303 // Constructor 302 304 FlexibleCylinderModel(); … … 312 314 // Model parameters 313 315 Parameter scale; 314 Parameter length;315 Parameter radius; 316 Parameter thickness;316 Parameter core_thick; 317 Parameter radius; 318 Parameter layer_thick; 317 319 Parameter core_sld; 318 320 Parameter layer_sld; 319 321 Parameter solvent_sld; 320 Parameter n layers;321 Parameter s pacing;322 Parameter background; 323 Parameter axis_theta; 324 Parameter axis_phi; 325 322 Parameter n_stacking; 323 Parameter sigma_d; 324 Parameter background; 325 Parameter axis_theta; 326 Parameter axis_phi; 327 326 328 // Constructor 327 329 StackedDisksModel(); … … 337 339 // Model parameters 338 340 Parameter scale; 339 Parameter delta;340 Parameter s igma;341 Parameter contrast;342 Parameter background; 343 341 Parameter bi_thick; 342 Parameter sld_bi; 343 Parameter sld_sol; 344 Parameter background; 345 344 346 // Constructor 345 347 LamellarModel(); … … 362 364 Parameter sld_solvent; 363 365 Parameter background; 364 366 365 367 // Constructor 366 368 LamellarFFHGModel(); … … 386 388 Parameter caille; 387 389 Parameter background; 388 390 389 391 // Constructor 390 392 LamellarPSModel(); … … 409 411 Parameter caille; 410 412 Parameter background; 411 413 412 414 // Constructor 413 415 LamellarPSHGModel(); 416 417 // Operators to get I(Q) 418 double operator()(double q); 419 double operator()(double qx, double qy); 420 double evaluate_rphi(double q, double phi); 421 }; 422 423 class CoreShellSpheroidModel{ 424 public: 425 // Model parameters 426 Parameter scale; 427 Parameter equat_core; 428 Parameter polar_core; 429 Parameter equat_shell; 430 Parameter polar_shell; 431 Parameter contrast; 432 Parameter sld_solvent; 433 Parameter background; 434 Parameter axis_theta; 435 Parameter axis_phi; 436 437 // Constructor 438 CoreShellSpheroidModel(); 414 439 415 440 // Operators to get I(Q) … … 432 457 Parameter axis_theta; 433 458 Parameter axis_phi; 434 459 435 460 // Constructor 436 461 OblateModel(); … … 454 479 Parameter axis_theta; 455 480 Parameter axis_phi; 456 481 457 482 // Constructor 458 483 ProlateModel(); … … 474 499 Parameter axis_theta; 475 500 Parameter axis_phi; 476 501 477 502 //Constructor 478 503 HollowCylinderModel(); 479 504 480 505 //Operators to get I(Q) 481 506 double operator()(double q); … … 495 520 Parameter n_pairs; 496 521 Parameter background; 497 522 498 523 //Constructor 499 524 MultiShellModel(); 500 525 501 526 //Operators to get I(Q) 502 527 double operator()(double q); … … 514 539 Parameter shell_sld; 515 540 Parameter background; 516 541 517 542 //Constructor 518 543 VesicleModel(); 519 544 520 545 //Operators to get I(Q) 521 546 double operator()(double q); … … 535 560 Parameter solvent_sld; 536 561 Parameter background; 537 562 538 563 //Constructor 539 564 BinaryHSModel(); 540 565 541 566 //Operators to get I(Q) 542 567 double operator()(double q); … … 556 581 Parameter solvent_sld; 557 582 Parameter background; 558 583 559 584 //Constructor 560 585 BinaryHSPSF11Model(); 561 586 562 587 //Operators to get I(Q) 563 588 double operator()(double q); -
sansmodels/src/sans/models/c_models/oblate.cpp
r96b59384 r975ec8e 121 121 * @return: function value 122 122 */ 123 123 /* 124 124 double OblateModel :: operator()(double qx, double qy) { 125 125 double q = sqrt(qx*qx + qy*qy); … … 127 127 return (*this).operator()(q); 128 128 } 129 129 */ 130 130 131 131 /** … … 140 140 } 141 141 142 /* disable below for now 143 double Oblate ShellModel :: operator()(double qx, double qy) {142 143 double OblateModel :: operator()(double qx, double qy) { 144 144 OblateParameters dp; 145 145 // Fill parameter array … … 233 233 return sum/norm + background(); 234 234 } 235 */// 235 -
sansmodels/src/sans/models/c_models/parallelepiped.cpp
r8a48713 r975ec8e 45 45 parallel_theta = Parameter(0.0, true); 46 46 parallel_phi = Parameter(0.0, true); 47 parallel_psi = Parameter(0.0, true); 47 48 } 48 49 … … 54 55 */ 55 56 double ParallelepipedModel :: operator()(double q) { 56 double dp[ 5];57 double dp[6]; 57 58 58 59 // Fill parameter array for IGOR library … … 63 64 dp[3] = longuest_edgeC(); 64 65 dp[4] = contrast(); 65 //dp[5] = background();66 66 dp[5] = 0.0; 67 67 … … 69 69 vector<WeightPoint> weights_short_edgeA; 70 70 short_edgeA.get_weights(weights_short_edgeA); 71 71 72 72 // Get the dispersion points for the longer_edgeB 73 73 vector<WeightPoint> weights_longer_edgeB; … … 83 83 double sum = 0.0; 84 84 double norm = 0.0; 85 85 86 86 // Loop over short_edgeA weight points 87 87 for(int i=0; i< (int)weights_short_edgeA.size(); i++) { … … 96 96 dp[3] = weights_longuest_edgeC[j].value; 97 97 98 sum += weights_short_edgeA[i].weight * weights_longer_edgeB[j].weight 98 sum += weights_short_edgeA[i].weight * weights_longer_edgeB[j].weight 99 99 * weights_longuest_edgeC[k].weight * Parallelepiped(dp, q); 100 100 … … 124 124 dp.parallel_theta = parallel_theta(); 125 125 dp.parallel_phi = parallel_phi(); 126 126 dp.parallel_psi = parallel_psi(); 127 127 128 128 129 // Get the dispersion points for the short_edgeA … … 146 147 parallel_phi.get_weights(weights_parallel_phi); 147 148 149 // Get angular averaging for psi 150 vector<WeightPoint> weights_parallel_psi; 151 parallel_psi.get_weights(weights_parallel_psi); 148 152 149 153 // Perform the computation, with all weight points … … 171 175 dp.parallel_phi = weights_parallel_phi[m].value; 172 176 173 double _ptvalue = weights_short_edgeA[i].weight 174 * weights_longer_edgeB[j].weight 175 * weights_longuest_edgeC[k].weight 176 * weights_parallel_theta[l].weight 177 * weights_parallel_phi[m].weight 178 * parallelepiped_analytical_2DXY(&dp, qx, qy); 179 if (weights_parallel_theta.size()>1) { 180 _ptvalue *= sin(weights_parallel_theta[l].value); 177 // Average over phi distribution 178 for(int n=0; n< (int)weights_parallel_psi.size(); n++) { 179 dp.parallel_psi = weights_parallel_psi[n].value; 180 181 double _ptvalue = weights_short_edgeA[i].weight 182 * weights_longer_edgeB[j].weight 183 * weights_longuest_edgeC[k].weight 184 * weights_parallel_theta[l].weight 185 * weights_parallel_phi[m].weight 186 * weights_parallel_psi[n].weight 187 * parallelepiped_analytical_2DXY(&dp, qx, qy); 188 if (weights_parallel_theta.size()>1) { 189 _ptvalue *= sin(weights_parallel_theta[l].value); 190 } 191 sum += _ptvalue; 192 193 norm += weights_short_edgeA[i].weight 194 * weights_longer_edgeB[j].weight 195 * weights_longuest_edgeC[k].weight 196 * weights_parallel_theta[l].weight 197 * weights_parallel_phi[m].weight 198 * weights_parallel_psi[n].weight; 181 199 } 182 sum += _ptvalue;183 184 norm += weights_short_edgeA[i].weight185 * weights_longer_edgeB[j].weight186 * weights_longuest_edgeC[k].weight187 * weights_parallel_theta[l].weight188 * weights_parallel_phi[m].weight;189 200 } 190 201 -
sansmodels/src/sans/models/c_models/stackeddisks.cpp
r9188cc1 r975ec8e 37 37 radius = Parameter(3000.0, true); 38 38 radius.set_min(0.0); 39 length= Parameter(10.0, true);40 length.set_min(0.0);41 thickness= Parameter(15.0);42 thickness.set_min(0.0);39 core_thick = Parameter(10.0, true); 40 core_thick.set_min(0.0); 41 layer_thick = Parameter(15.0); 42 layer_thick.set_min(0.0); 43 43 core_sld = Parameter(4.0e-6); 44 44 layer_sld = Parameter(-4.0e-7); 45 45 solvent_sld = Parameter(5.0e-6); 46 n layers= Parameter(1);47 s pacing= Parameter(0);46 n_stacking = Parameter(1); 47 sigma_d = Parameter(0); 48 48 background = Parameter(0.001); 49 49 axis_theta = Parameter(0.0, true); … … 64 64 dp[0] = scale(); 65 65 dp[1] = radius(); 66 dp[2] = length();67 dp[3] = thickness();66 dp[2] = core_thick(); 67 dp[3] = layer_thick(); 68 68 dp[4] = core_sld(); 69 69 dp[5] = layer_sld(); 70 70 dp[6] = solvent_sld(); 71 dp[7] = n layers();72 dp[8] = s pacing();71 dp[7] = n_stacking(); 72 dp[8] = sigma_d(); 73 73 dp[9] = 0.0; 74 75 // Get the dispersion points for the length76 vector<WeightPoint> weights_length;77 length.get_weights(weights_length);78 74 79 75 // Get the dispersion points for the radius … … 81 77 radius.get_weights(weights_radius); 82 78 83 // Get the dispersion points for the thickness 84 vector<WeightPoint> weights_thickness; 85 thickness.get_weights(weights_thickness); 79 // Get the dispersion points for the core_thick 80 vector<WeightPoint> weights_core_thick; 81 core_thick.get_weights(weights_core_thick); 82 83 // Get the dispersion points for the layer_thick 84 vector<WeightPoint> weights_layer_thick; 85 layer_thick.get_weights(weights_layer_thick); 86 86 87 87 // Perform the computation, with all weight points … … 90 90 91 91 // Loop over length weight points 92 for(int i=0; i< (int)weights_ length.size(); i++) {93 dp[1] = weights_ length[i].value;92 for(int i=0; i< (int)weights_radius.size(); i++) { 93 dp[1] = weights_radius[i].value; 94 94 95 95 // Loop over radius weight points 96 for(int j=0; j< (int)weights_ radius.size(); j++) {97 dp[2] = weights_ radius[j].value;96 for(int j=0; j< (int)weights_core_thick.size(); j++) { 97 dp[2] = weights_core_thick[j].value; 98 98 99 99 // Loop over thickness weight points 100 for(int k=0; k< (int)weights_ radius.size(); k++) {101 dp[3] = weights_ radius[k].value;102 103 sum += weights_ length[i].weight104 * weights_ radius[j].weight * weights_thickness[k].weight* StackedDiscs(dp, q);105 norm += weights_ length[i].weight106 * weights_ radius[j].weight* weights_thickness[k].weight;100 for(int k=0; k< (int)weights_layer_thick.size(); k++) { 101 dp[3] = weights_layer_thick[k].value; 102 103 sum += weights_radius[i].weight 104 * weights_core_thick[j].weight * weights_layer_thick[k].weight* StackedDiscs(dp, q); 105 norm += weights_radius[i].weight 106 * weights_core_thick[j].weight* weights_layer_thick[k].weight; 107 107 } 108 108 } … … 121 121 // Fill parameter array 122 122 dp.scale = scale(); 123 dp. length = length();123 dp.core_thick = core_thick(); 124 124 dp.radius = radius(); 125 dp. thickness = thickness();125 dp.core_thick = core_thick(); 126 126 dp.core_sld = core_sld(); 127 127 dp.layer_sld = layer_sld(); 128 128 dp.solvent_sld= solvent_sld(); 129 dp.n layers = nlayers();130 dp.s pacing = spacing();129 dp.n_stacking = n_stacking(); 130 dp.sigma_d = sigma_d(); 131 131 dp.background = 0.0; 132 132 dp.axis_theta = axis_theta(); … … 134 134 135 135 // Get the dispersion points for the length 136 vector<WeightPoint> weights_ length;137 length.get_weights(weights_length);136 vector<WeightPoint> weights_core_thick; 137 core_thick.get_weights(weights_core_thick); 138 138 139 139 // Get the dispersion points for the radius … … 142 142 143 143 // Get the dispersion points for the thickness 144 vector<WeightPoint> weights_ thickness;145 thickness.get_weights(weights_thickness);144 vector<WeightPoint> weights_layer_thick; 145 layer_thick.get_weights(weights_layer_thick); 146 146 147 147 // Get angular averaging for theta … … 158 158 159 159 // Loop over length weight points 160 for(int i=0; i< (int)weights_ length.size(); i++) {161 dp. length = weights_length[i].value;160 for(int i=0; i< (int)weights_core_thick.size(); i++) { 161 dp.core_thick = weights_core_thick[i].value; 162 162 163 163 // Loop over radius weight points … … 166 166 167 167 // Loop over thickness weight points 168 for(int k=0; k< (int)weights_ thickness.size(); k++) {169 dp. thickness = weights_thickness[k].value;168 for(int k=0; k< (int)weights_layer_thick.size(); k++) { 169 dp.layer_thick = weights_layer_thick[k].value; 170 170 171 171 for(int l=0; l< (int)weights_theta.size(); l++) { … … 176 176 dp.axis_phi = weights_phi[m].value; 177 177 178 double _ptvalue = weights_ length[i].weight178 double _ptvalue = weights_core_thick[i].weight 179 179 * weights_radius[j].weight 180 * weights_ thickness[k].weight180 * weights_layer_thick[k].weight 181 181 * weights_theta[l].weight 182 182 * weights_phi[m].weight … … 187 187 sum += _ptvalue; 188 188 189 norm += weights_ length[i].weight189 norm += weights_core_thick[i].weight 190 190 * weights_radius[j].weight 191 * weights_ thickness[k].weight191 * weights_layer_thick[k].weight 192 192 * weights_theta[l].weight 193 193 * weights_phi[m].weight; -
sansmodels/src/sans/models/c_models/triaxialellipsoid.cpp
r9188cc1 r975ec8e 44 44 axis_theta = Parameter(0.0, true); 45 45 axis_phi = Parameter(0.0, true); 46 axis_psi = Parameter(0.0, true); 46 47 } 47 48 … … 53 54 */ 54 55 double TriaxialEllipsoidModel :: operator()(double q) { 55 double dp[ 5];56 double dp[6]; 56 57 57 58 // Fill parameter array for IGOR library … … 119 120 dp.axis_theta = axis_theta(); 120 121 dp.axis_phi = axis_phi(); 122 dp.axis_psi = axis_psi(); 121 123 122 124 // Get the dispersion points for the semi_axis A … … 140 142 axis_phi.get_weights(weights_phi); 141 143 144 // Get angular averaging for psi 145 vector<WeightPoint> weights_psi; 146 axis_psi.get_weights(weights_psi); 147 142 148 // Perform the computation, with all weight points 143 149 double sum = 0.0; … … 163 169 for(int m=0; m <(int)weights_phi.size(); m++) { 164 170 dp.axis_phi = weights_phi[m].value; 165 166 double _ptvalue = weights_semi_axisA[i].weight 167 * weights_semi_axisB[j].weight 168 * weights_semi_axisC[k].weight 169 * weights_theta[l].weight 170 * weights_phi[m].weight 171 * triaxial_ellipsoid_analytical_2DXY(&dp, qx, qy); 172 if (weights_theta.size()>1) { 173 _ptvalue *= sin(weights_theta[k].value); 171 // Average over psi distribution 172 for(int n=0; n <(int)weights_psi.size(); n++) { 173 dp.axis_psi = weights_psi[n].value; 174 175 double _ptvalue = weights_semi_axisA[i].weight 176 * weights_semi_axisB[j].weight 177 * weights_semi_axisC[k].weight 178 * weights_theta[l].weight 179 * weights_phi[m].weight 180 * weights_psi[n].weight 181 * triaxial_ellipsoid_analytical_2DXY(&dp, qx, qy); 182 if (weights_theta.size()>1) { 183 _ptvalue *= sin(weights_theta[k].value); 184 } 185 sum += _ptvalue; 186 187 norm += weights_semi_axisA[i].weight 188 * weights_semi_axisB[j].weight 189 * weights_semi_axisC[k].weight 190 * weights_theta[l].weight 191 * weights_phi[m].weight 192 * weights_psi[m].weight; 174 193 } 175 sum += _ptvalue;176 177 norm += weights_semi_axisA[i].weight178 * weights_semi_axisB[j].weight179 * weights_semi_axisC[k].weight180 * weights_theta[l].weight181 * weights_phi[m].weight;182 194 } 183 195
Note: See TracChangeset
for help on using the changeset viewer.