Changeset c1c29b6 in sasview for sansmodels/src/sans/models


Ignore:
Timestamp:
Aug 21, 2009 12:34:33 PM (15 years ago)
Author:
Jae Cho <jhjcho@…>
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:
14c3887
Parents:
bda194e3
Message:

some corrections and removed polydispersity from inside of function and set dQ =0

Location:
sansmodels/src/sans/models
Files:
10 edited

Legend:

Unmodified
Added
Removed

  • sansmodels/src/sans/models/LamellarPSModel.py

    r9bd69098 rc1c29b6  
    3535         spacing         = 400.0 [A] 
    3636         delta           = 30.0 [A] 
    37          sigma           = 0.15  
    3837         contrast        = 5.3e-006 [1/A²] 
    3938         n_plates        = 20.0  
    … …  
    6665                *Parameters: spacing = repeat spacing, 
    6766                delta = bilayer thickness, 
    68                 sigma = variation in bilayer thickness 
    6967                contrast = SLD_solvent - SLD_bilayer 
    7068                n_plate = # of Lamellar plates 
    … …  
    7876        self.details['spacing'] = ['[A]', None, None] 
    7977        self.details['delta'] = ['[A]', None, None] 
    80         self.details['sigma'] = ['', None, None] 
    8178        self.details['contrast'] = ['[1/A²]', None, None] 
    8279        self.details['n_plates'] = ['', None, None] 
    … …  
    8582 
    8683                ## fittable parameters 
    87         self.fixed=['spacing.width'] 
     84        self.fixed=['delta.width', 'spacing.width'] 
    8885         
    8986        ## parameters with orientation 

  • sansmodels/src/sans/models/c_extensions/ellipsoid.h

    r1ed3834 rc1c29b6  
    2828 //             </text> 
    2929 //[FIXED]= <text> axis_phi.width; axis_theta.width;radius_a.width; 
    30  //radius_b.width; length.width; r_minor.width 
    31  //, r_ratio.width</text> 
     30 //radius_b.width; length.width; r_minor.width; 
     31 //r_ratio.width</text> 
    3232 //[ORIENTATION_PARAMS]=  axis_phi.width; axis_theta.width;axis_phi; axis_theta 
    3333 

  • sansmodels/src/sans/models/c_extensions/lamellar.c

    r0cfeff4 rc1c29b6  
    66#include "lamellar.h" 
    77#include <math.h> 
    8 #include "libCylinder.h" 
     8//#include "libCylinder.h" 
    99#include <stdio.h> 
    1010#include <stdlib.h> 

  • sansmodels/src/sans/models/c_extensions/lamellar.h

    r975ec8e rc1c29b6  
    3636} LamellarParameters; 
    3737 
    38  
    39  
     38/// kernel 
     39double lamellar_kernel(double dp[], double q); 
    4040/// 1D scattering function 
    4141double lamellar_analytical_1D(LamellarParameters *pars, double q); 

  • sansmodels/src/sans/models/c_extensions/lamellarPS.c

    rb4679de rc1c29b6  
    1010#include <stdlib.h> 
    1111 
     12/*LamellarPS_kernel() was moved from libigor to get rid of polydipersity in del(thickness) that we provide from control panel. 
     13/*      LamellarPSX  :  calculates the form factor of a lamellar structure - with S(q) effects included 
     14------- 
     15------- resolution effects ARE NOT included, but only a CONSTANT default value, not the real q-dependent resolution!! 
     16 
     17        */ 
     18double 
     19LamellarPS_kernel(double dp[], double q) 
     20{ 
     21        double scale,dd,del,contr,NN,Cp,bkg;            //local variables of coefficient wave 
     22        double inten, qval,Pq,Sq,alpha,temp,t1,t2,t3,dQ; 
     23        double Pi,Euler,dQDefault,fii; 
     24        int ii,NNint; 
     25        Euler = 0.5772156649;           // Euler's constant 
     26        dQDefault = 0.0;//0.0025;               //[=] 1/A, q-resolution, default value 
     27        dQ = dQDefault; 
     28 
     29        Pi = 4.0*atan(1.0); 
     30        qval = q; 
     31 
     32        scale = dp[0]; 
     33        dd = dp[1]; 
     34        del = dp[2]; 
     35        contr = dp[3]; 
     36        NN = trunc(dp[4]);              //be sure that NN is an integer 
     37        Cp = dp[5]; 
     38        bkg = dp[6]; 
     39 
     40        Pq = 2.0*contr*contr/qval/qval*(1.0-cos(qval*del)); 
     41 
     42        NNint = (int)NN;                //cast to an integer for the loop 
     43        ii=0; 
     44        Sq = 0.0; 
     45        for(ii=1;ii<(NNint-1);ii+=1) { 
     46 
     47                fii = (double)ii;               //do I really need to do this? 
     48 
     49                temp = 0.0; 
     50                alpha = Cp/4.0/Pi/Pi*(log(Pi*ii) + Euler); 
     51                t1 = 2.0*dQ*dQ*dd*dd*alpha; 
     52                t2 = 2.0*qval*qval*dd*dd*alpha; 
     53                t3 = dQ*dQ*dd*dd*ii*ii; 
     54 
     55                temp = 1.0-ii/NN; 
     56                temp *= cos(dd*qval*ii/(1.0+t1)); 
     57                temp *= exp(-1.0*(t2 + t3)/(2.0*(1.0+t1)) ); 
     58                temp /= sqrt(1.0+t1); 
     59 
     60                Sq += temp; 
     61        } 
     62 
     63        Sq *= 2.0; 
     64        Sq += 1.0; 
     65 
     66        inten = 2.0*Pi*scale*Pq*Sq/(dd*qval*qval); 
     67 
     68        inten *= 1.0e8;         // 1/A to 1/cm 
     69 
     70    return(inten+bkg); 
     71} 
    1272 
    1373/** 
    … …  
    1878 */ 
    1979double lamellarPS_analytical_1D(LamellarPSParameters *pars, double q) { 
    20         double dp[8]; 
     80        double dp[7]; 
    2181 
    2282        // Fill paramater array 
    … …  
    2484        dp[1] = pars->spacing; 
    2585        dp[2] = pars->delta; 
    26         dp[3] = pars->sigma; 
    27         dp[4] = pars->contrast; 
    28         dp[5] = pars->n_plates; 
    29         dp[6] = pars->caille; 
    30         dp[7] = pars->background; 
     86        dp[3] = pars->contrast; 
     87        dp[4] = pars->n_plates; 
     88        dp[5] = pars->caille; 
     89        dp[6] = pars->background; 
    3190 
    3291        // Call library function to evaluate model 
    33         return LamellarPS(dp, q); 
     92        return LamellarPS_kernel(dp, q); 
    3493} 
    3594/** 

  • sansmodels/src/sans/models/c_extensions/lamellarPS.h

    r96b59384 rc1c29b6  
    77/** Structure definition for concentrated lamellar form factor parameters 
    88 * [PYTHONCLASS] = LamellarPSModel 
    9  * [DISP_PARAMS] = spacing 
     9 * [DISP_PARAMS] = delta, spacing 
    1010   [DESCRIPTION] = <text>[Concentrated Lamellar Form Factor] Calculates the scattered 
    1111           intensity from a lyotropic lamellar phase. 
    … …  
    2121                *Parameters: spacing = repeat spacing, 
    2222                delta = bilayer thickness, 
    23                 sigma = variation in bilayer thickness 
    2423                contrast = SLD_solvent - SLD_bilayer 
    2524                n_plate = # of Lamellar plates 
    … …  
    2827                scale = scale factor 
    2928</text> 
    30    [FIXED]= spacing.width 
     29   [FIXED]= <text>delta.width; spacing.width</text> 
    3130   [ORIENTATION_PARAMS]= 
    3231 
    … …  
    4241    //  [DEFAULT]=delta=30 [A] 
    4342    double delta; 
    44         /// polydispersity of the bilayer thickness  [A] 
    45     //  [DEFAULT]=sigma=0.15 
    46     double sigma; 
    4743    /// Contrast [1/A²] 
    4844    //  [DEFAULT]=contrast=5.3e-6 [1/A²] 
    … …  
    6056} LamellarPSParameters; 
    6157 
    62  
    63  
     58/// kernel 
     59double LamellarPS_kernel(double dp[], double q); 
    6460/// 1D scattering function 
    6561double lamellarPS_analytical_1D(LamellarPSParameters *pars, double q); 

  • sansmodels/src/sans/models/c_models/CLamellarPSModel.cpp

    r9bd69098 rc1c29b6  
    9292        PyDict_SetItemString(self->params,"background",Py_BuildValue("d",0.000000)); 
    9393        PyDict_SetItemString(self->params,"delta",Py_BuildValue("d",30.000000)); 
    94         PyDict_SetItemString(self->params,"sigma",Py_BuildValue("d",0.150000)); 
    9594        PyDict_SetItemString(self->params,"contrast",Py_BuildValue("d",0.000005)); 
    9695        // Initialize dispersion / averaging parameter dict 
    9796        DispersionVisitor* visitor = new DispersionVisitor(); 
    9897        PyObject * disp_dict; 
     98        disp_dict = PyDict_New(); 
     99        self->model->delta.dispersion->accept_as_source(visitor, self->model->delta.dispersion, disp_dict); 
     100        PyDict_SetItemString(self->dispersion, "delta", disp_dict); 
    99101        disp_dict = PyDict_New(); 
    100102        self->model->spacing.dispersion->accept_as_source(visitor, self->model->spacing.dispersion, disp_dict); 
    … …  
    167169    return PyArray_Return(result);  
    168170 } 
    169 /** 
    170  * Function to call to evaluate model 
    171  * @param args: input numpy array  [q[],phi[]] 
    172  * @return: numpy array object  
    173  */ 
    174 static PyObject * evaluateTwoDim( LamellarPSModel* model,  
    175                               PyArrayObject *q, PyArrayObject *phi) 
    176  { 
    177     PyArrayObject *result; 
    178     //check validity of input vectors 
    179     if (q->nd != 1 || q->descr->type_num != PyArray_DOUBLE 
    180         || phi->nd != 1 || phi->descr->type_num != PyArray_DOUBLE 
    181         || phi->dimensions[0] != q->dimensions[0]){ 
    182       
    183         //const char * message= "Invalid array: q->nd=%d,type_num=%d\n",q->nd,q->descr->type_num; 
    184         PyErr_SetString(PyExc_ValueError ,"wrong input");  
    185         return NULL; 
    186     } 
    187         result= (PyArrayObject *)PyArray_FromDims(q->nd,(int*)(q->dimensions), PyArray_DOUBLE); 
    188  
    189         if (result == NULL){ 
    190             const char * message= "Could not create result "; 
    191         PyErr_SetString(PyExc_RuntimeError , message); 
    192             return NULL; 
    193         } 
    194          
    195     for (int i = 0; i < q->dimensions[0]; i++) { 
    196       double q_value = *(double *)(q->data + i*q->strides[0]); 
    197       double phi_value = *(double *)(phi->data + i*phi->strides[0]); 
    198       double *result_value = (double *)(result->data + i*result->strides[0]); 
    199       if (q_value == 0) 
    200           *result_value = 0.0; 
    201       else 
    202           *result_value = model->evaluate_rphi(q_value, phi_value); 
    203     } 
    204     return PyArray_Return(result);  
    205  } 
     171 
    206172 /** 
    207173 * Function to call to evaluate model 
    … …  
    272238    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 
    273239    self->model->delta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "delta") ); 
    274     self->model->sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") ); 
    275240    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") ); 
    276241    // Read in dispersion parameters 
    277242    PyObject* disp_dict; 
    278243    DispersionVisitor* visitor = new DispersionVisitor(); 
     244    disp_dict = PyDict_GetItemString(self->dispersion, "delta"); 
     245    self->model->delta.dispersion->accept_as_destination(visitor, self->model->delta.dispersion, disp_dict); 
    279246    disp_dict = PyDict_GetItemString(self->dispersion, "spacing"); 
    280247    self->model->spacing.dispersion->accept_as_destination(visitor, self->model->spacing.dispersion, disp_dict); 
    … …  
    347314    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 
    348315    self->model->delta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "delta") ); 
    349     self->model->sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") ); 
    350316    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") ); 
    351317    // Read in dispersion parameters 
    352318    PyObject* disp_dict; 
    353319    DispersionVisitor* visitor = new DispersionVisitor(); 
     320    disp_dict = PyDict_GetItemString(self->dispersion, "delta"); 
     321    self->model->delta.dispersion->accept_as_destination(visitor, self->model->delta.dispersion, disp_dict); 
    354322    disp_dict = PyDict_GetItemString(self->dispersion, "spacing"); 
    355323    self->model->spacing.dispersion->accept_as_destination(visitor, self->model->spacing.dispersion, disp_dict); 
    … …  
    411379    self->model->background = PyFloat_AsDouble( PyDict_GetItemString(self->params, "background") ); 
    412380    self->model->delta = PyFloat_AsDouble( PyDict_GetItemString(self->params, "delta") ); 
    413     self->model->sigma = PyFloat_AsDouble( PyDict_GetItemString(self->params, "sigma") ); 
    414381    self->model->contrast = PyFloat_AsDouble( PyDict_GetItemString(self->params, "contrast") ); 
    415382    // Read in dispersion parameters 
    416383    PyObject* disp_dict; 
    417384    DispersionVisitor* visitor = new DispersionVisitor(); 
     385    disp_dict = PyDict_GetItemString(self->dispersion, "delta"); 
     386    self->model->delta.dispersion->accept_as_destination(visitor, self->model->delta.dispersion, disp_dict); 
    418387    disp_dict = PyDict_GetItemString(self->dispersion, "spacing"); 
    419388    self->model->spacing.dispersion->accept_as_destination(visitor, self->model->spacing.dispersion, disp_dict); 
    … …  
    473442        // Ugliness necessary to go from python to C 
    474443            // TODO: refactor this 
    475     if (!strcmp(par_name, "spacing")) { 
     444    if (!strcmp(par_name, "delta")) { 
     445        self->model->delta.dispersion = dispersion; 
     446    } else    if (!strcmp(par_name, "spacing")) { 
    476447        self->model->spacing.dispersion = dispersion; 
    477448    } else { 

  • sansmodels/src/sans/models/c_models/lamellar.cpp

    r975ec8e rc1c29b6  
    2828 
    2929extern "C" { 
    30         #include "libCylinder.h" 
     30//      #include "libCylinder.h" 
    3131        #include "lamellar.h" 
    3232} 

  • sansmodels/src/sans/models/c_models/lamellarPS.cpp

    rb4679de rc1c29b6  
    3939        delta     = Parameter(30.0); 
    4040        delta.set_min(0.0); 
    41         sigma    = Parameter(0.15); 
    42         sigma.set_min(0.0); 
    4341        contrast   = Parameter(5.3e-6); 
    4442        n_plates     = Parameter(20.0); 
    … …  
    5553 */ 
    5654double LamellarPSModel :: operator()(double q) { 
    57         double dp[8]; 
     55        double dp[7]; 
    5856 
    5957        // Fill parameter array for IGOR library 
    … …  
    6260        dp[1] = spacing(); 
    6361        dp[2] = delta(); 
    64         dp[3] = sigma(); 
    65         dp[4] = contrast(); 
    66         dp[5] = n_plates(); 
    67         dp[6] = caille(); 
    68         dp[7] = 0.0; 
     62        dp[3] = contrast(); 
     63        dp[4] = n_plates(); 
     64        dp[5] = caille(); 
     65        dp[6] = 0.0; 
    6966 
    7067 
    71         // Get the dispersion points for (delta) thickness 
     68        // Get the dispersion points for spacing and delta (thickness) 
    7269        vector<WeightPoint> weights_spacing; 
    7370        spacing.get_weights(weights_spacing); 
     71        vector<WeightPoint> weights_delta; 
     72        delta.get_weights(weights_delta); 
    7473 
    7574        // Perform the computation, with all weight points 
    … …  
    8079        for(int i=0; i< (int)weights_spacing.size(); i++) { 
    8180                dp[1] = weights_spacing[i].value; 
     81                for(int j=0; j< (int)weights_spacing.size(); j++) { 
     82                        dp[2] = weights_delta[i].value; 
    8283 
    83                 sum += weights_spacing[i].weight * LamellarPS(dp, q); 
    84                 norm += weights_spacing[i].weight; 
    85  
     84                        sum += weights_spacing[i].weight * weights_delta[j].weight * LamellarPS_kernel(dp, q); 
     85                        norm += weights_spacing[i].weight * weights_delta[j].weight; 
     86                } 
    8687        } 
    8788        return sum/norm + background(); 
    … …  
    108109        return (*this).operator()(q); 
    109110} 
     111 

  • sansmodels/src/sans/models/c_models/models.hh

    r8e36cdd rc1c29b6  
    2020        #include "cylinder.h" 
    2121        #include "parallelepiped.h" 
     22        #include "lamellarPS.h" 
     23        #include "lamellar.h" 
    2224} 
    2325 
    … …  
    383385        Parameter spacing; 
    384386        Parameter delta; 
    385         Parameter sigma; 
    386387        Parameter contrast; 
    387388        Parameter n_plates; 
Note: See TracChangeset for help on using the changeset viewer.