Changes in / [d9b7197:2fe721be] in sasview


Ignore:
Location:
src/sas/sascalc/calculator/c_extensions
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • src/sas/sascalc/calculator/c_extensions/libfunc.c

    rf54e82cf rb6c8abe  
    9595// spintheta: angle (anti-clock-wise) between neutron spin(up) and x axis 
    9696// Note: all angles are in degrees. 
    97 void cal_msld(polar_sld *p_sld, int isangle, double qx, double qy, double bn, 
     97polar_sld cal_msld(int isangle, double qx, double qy, double bn, 
    9898                                double m01, double mtheta1, double mphi1, 
    9999                                double spinfraci, double spinfracf, double spintheta) 
    … …  
    124124        double my = 0.0; 
    125125        double mz = 0.0; 
    126         double uu = sld; 
    127         double dd = sld; 
    128         double re_ud = 0.0; 
    129         double im_ud = 0.0; 
    130         double re_du = 0.0; 
    131         double im_du = 0.0; 
     126        polar_sld p_sld; 
     127        p_sld.uu = sld; 
     128        p_sld.dd = sld; 
     129        p_sld.re_ud = 0.0; 
     130        p_sld.im_ud = 0.0; 
     131        p_sld.re_du = 0.0; 
     132        p_sld.im_du = 0.0; 
    132133 
    133134        //No mag means no further calculation 
    134         if (isangle>0) { 
     135        if (isangle>0){ 
    135136                if (m_max < 1.0e-32){ 
    136                         uu = sqrt(sqrt(in_spin * out_spin)) * uu; 
    137                         dd = sqrt(sqrt((1.0 - in_spin) * (1.0 - out_spin))) * dd; 
    138                 } 
    139         } 
    140         else if (fabs(m_max)< 1.0e-32 && fabs(m_phi)< 1.0e-32 && fabs(m_theta)< 1.0e-32){ 
    141                         uu = sqrt(sqrt(in_spin * out_spin)) * uu; 
    142                         dd = sqrt(sqrt((1.0 - in_spin) * (1.0 - out_spin))) * dd; 
    143         } else { 
    144  
    145                 //These are needed because of the precision of inputs 
    146                 if (in_spin < 0.0) in_spin = 0.0; 
    147                 if (in_spin > 1.0) in_spin = 1.0; 
    148                 if (out_spin < 0.0) out_spin = 0.0; 
    149                 if (out_spin > 1.0) out_spin = 1.0; 
    150  
    151                 if (q_x == 0.0) q_angle = pi / 2.0; 
    152                 else q_angle = atan(q_y/q_x); 
    153                 if (q_y < 0.0 && q_x < 0.0) q_angle -= pi; 
    154                 else if (q_y > 0.0 && q_x < 0.0) q_angle += pi; 
    155  
    156                 q_angle = pi/2.0 - q_angle; 
    157                 if (q_angle > pi) q_angle -= 2.0 * pi; 
    158                 else if (q_angle < -pi) q_angle += 2.0 * pi; 
    159  
    160                 if (fabs(q_x) < 1.0e-16 && fabs(q_y) < 1.0e-16){ 
    161                         m_perp = 0.0; 
    162                         } 
    163                 else { 
    164                         m_perp = m_max; 
    165                         } 
    166                 if (is_angle > 0){ 
    167                         m_phi *= pi/180.0; 
    168                         m_theta *= pi/180.0; 
    169                         mx = m_perp * cos(m_theta) * cos(m_phi); 
    170                         my = m_perp * sin(m_theta); 
    171                         mz = -(m_perp * cos(m_theta) * sin(m_phi)); 
    172                 } 
    173                 else{ 
    174                         mx = m_perp; 
    175                         my = m_phi; 
    176                         mz = m_theta; 
    177                 } 
    178                 //ToDo: simplify these steps 
    179                 // m_perp1 -m_perp2 
    180                 m_perp_x = (mx) *  cos(q_angle); 
    181                 m_perp_x -= (my) * sin(q_angle); 
    182                 m_perp_y = m_perp_x; 
    183                 m_perp_x *= cos(-q_angle); 
    184                 m_perp_y *= sin(-q_angle); 
    185                 m_perp_z = mz; 
    186  
    187                 m_sigma_x = (m_perp_x * cos(-s_theta) - m_perp_y * sin(-s_theta)); 
    188                 m_sigma_y = (m_perp_x * sin(-s_theta) + m_perp_y * cos(-s_theta)); 
    189                 m_sigma_z = (m_perp_z); 
    190  
    191                 //Find b 
    192                 uu -= m_sigma_x; 
    193                 dd += m_sigma_x; 
    194                 re_ud = m_sigma_y; 
    195                 re_du = m_sigma_y; 
    196                 im_ud = m_sigma_z; 
    197                 im_du = -m_sigma_z; 
    198  
    199                 uu = sqrt(sqrt(in_spin * out_spin)) * uu; 
    200                 dd = sqrt(sqrt((1.0 - in_spin) * (1.0 - out_spin))) * dd; 
    201  
    202                 re_ud = sqrt(sqrt(in_spin * (1.0 - out_spin))) * re_ud; 
    203                 im_ud = sqrt(sqrt(in_spin * (1.0 - out_spin))) * im_ud; 
    204                 re_du = sqrt(sqrt((1.0 - in_spin) * out_spin)) * re_du; 
    205                 im_du = sqrt(sqrt((1.0 - in_spin) * out_spin)) * im_du; 
    206         } 
    207         p_sld->uu = uu; 
    208         p_sld->dd = dd; 
    209         p_sld->re_ud = re_ud; 
    210         p_sld->im_ud = im_ud; 
    211         p_sld->re_du = re_du; 
    212         p_sld->im_du = im_du; 
     137                        p_sld.uu = sqrt(sqrt(in_spin * out_spin)) * p_sld.uu; 
     138                        p_sld.dd = sqrt(sqrt((1.0 - in_spin) * (1.0 - out_spin))) * p_sld.dd; 
     139                        return p_sld; 
     140                } 
     141        } 
     142        else{ 
     143                if (fabs(m_max)< 1.0e-32 && fabs(m_phi)< 1.0e-32 && fabs(m_theta)< 1.0e-32){ 
     144                        p_sld.uu = sqrt(sqrt(in_spin * out_spin)) * p_sld.uu; 
     145                        p_sld.dd = sqrt(sqrt((1.0 - in_spin) * (1.0 - out_spin))) * p_sld.dd; 
     146                        return p_sld; 
     147                } 
     148        } 
     149 
     150        //These are needed because of the precision of inputs 
     151        if (in_spin < 0.0) in_spin = 0.0; 
     152        if (in_spin > 1.0) in_spin = 1.0; 
     153        if (out_spin < 0.0) out_spin = 0.0; 
     154        if (out_spin > 1.0) out_spin = 1.0; 
     155 
     156        if (q_x == 0.0) q_angle = pi / 2.0; 
     157        else q_angle = atan(q_y/q_x); 
     158        if (q_y < 0.0 && q_x < 0.0) q_angle -= pi; 
     159        else if (q_y > 0.0 && q_x < 0.0) q_angle += pi; 
     160 
     161        q_angle = pi/2.0 - q_angle; 
     162        if (q_angle > pi) q_angle -= 2.0 * pi; 
     163        else if (q_angle < -pi) q_angle += 2.0 * pi; 
     164 
     165        if (fabs(q_x) < 1.0e-16 && fabs(q_y) < 1.0e-16){ 
     166                m_perp = 0.0; 
     167                } 
     168        else { 
     169                m_perp = m_max; 
     170                } 
     171        if (is_angle > 0){ 
     172                m_phi *= pi/180.0; 
     173                m_theta *= pi/180.0; 
     174                mx = m_perp * cos(m_theta) * cos(m_phi); 
     175                my = m_perp * sin(m_theta); 
     176                mz = -(m_perp * cos(m_theta) * sin(m_phi)); 
     177        } 
     178        else{ 
     179                mx = m_perp; 
     180                my = m_phi; 
     181                mz = m_theta; 
     182        } 
     183        //ToDo: simplify these steps 
     184        // m_perp1 -m_perp2 
     185        m_perp_x = (mx) *  cos(q_angle); 
     186        m_perp_x -= (my) * sin(q_angle); 
     187        m_perp_y = m_perp_x; 
     188        m_perp_x *= cos(-q_angle); 
     189        m_perp_y *= sin(-q_angle); 
     190        m_perp_z = mz; 
     191 
     192        m_sigma_x = (m_perp_x * cos(-s_theta) - m_perp_y * sin(-s_theta)); 
     193        m_sigma_y = (m_perp_x * sin(-s_theta) + m_perp_y * cos(-s_theta)); 
     194        m_sigma_z = (m_perp_z); 
     195 
     196        //Find b 
     197        p_sld.uu -= m_sigma_x; 
     198        p_sld.dd += m_sigma_x; 
     199        p_sld.re_ud = m_sigma_y; 
     200        p_sld.re_du = m_sigma_y; 
     201        p_sld.im_ud = m_sigma_z; 
     202        p_sld.im_du = -m_sigma_z; 
     203 
     204        p_sld.uu = sqrt(sqrt(in_spin * out_spin)) * p_sld.uu; 
     205        p_sld.dd = sqrt(sqrt((1.0 - in_spin) * (1.0 - out_spin))) * p_sld.dd; 
     206 
     207        p_sld.re_ud = sqrt(sqrt(in_spin * (1.0 - out_spin))) * p_sld.re_ud; 
     208        p_sld.im_ud = sqrt(sqrt(in_spin * (1.0 - out_spin))) * p_sld.im_ud; 
     209        p_sld.re_du = sqrt(sqrt((1.0 - in_spin) * out_spin)) * p_sld.re_du; 
     210        p_sld.im_du = sqrt(sqrt((1.0 - in_spin) * out_spin)) * p_sld.im_du; 
     211 
     212        return p_sld; 
    213213} 
    214214 

  • src/sas/sascalc/calculator/c_extensions/libfunc.h

    rf54e82cf rb6c8abe  
    1818//double gamln(double x); 
    1919 
    20 void cal_msld(polar_sld*, int isangle, double qx, double qy, double bn, double m01, double mtheta1,  
     20polar_sld cal_msld(int isangle, double qx, double qy, double bn, double m01, double mtheta1,  
    2121                        double mphi1, double spinfraci, double spinfracf, double spintheta); 
    2222 

  • src/sas/sascalc/calculator/c_extensions/librefl.c

    rf54e82cf r4c29e4d  
    103103#endif // NEED_ERF 
    104104 
    105 void cassign(Cplx *x, double real, double imag) 
    106 { 
    107         x->re = real; 
    108         x->im = imag; 
    109 } 
    110  
    111  
    112 void cplx_add(Cplx *z, Cplx x, Cplx y) 
    113 { 
    114         z->re = x.re + y.re; 
    115         z->im = x.im + y.im; 
    116 } 
    117  
    118 void rcmult(Cplx *z, double x, Cplx y) 
    119 { 
    120         z->re = x*y.re; 
    121         z->im = x*y.im; 
    122 } 
    123  
    124 void cplx_sub(Cplx *z, Cplx x, Cplx y) 
    125 { 
    126         z->re = x.re - y.re; 
    127         z->im = x.im - y.im; 
    128 } 
    129  
    130  
    131 void cplx_mult(Cplx *z, Cplx x, Cplx y) 
    132 { 
    133         z->re = x.re*y.re - x.im*y.im; 
    134         z->im = x.re*y.im + x.im*y.re; 
    135 } 
    136  
    137 void cplx_div(Cplx *z, Cplx x, Cplx y) 
    138 { 
    139         z->re = (x.re*y.re + x.im*y.im)/(y.re*y.re + y.im*y.im); 
    140         z->im = (x.im*y.re - x.re*y.im)/(y.re*y.re + y.im*y.im); 
    141 } 
    142  
    143 void cplx_exp(Cplx *z, Cplx b) 
    144 { 
     105complex cassign(real, imag) 
     106        double real, imag; 
     107{ 
     108        complex x; 
     109        x.re = real; 
     110        x.im = imag; 
     111        return x; 
     112} 
     113 
     114 
     115complex cplx_add(x,y) 
     116        complex x,y; 
     117{ 
     118        complex z; 
     119        z.re = x.re + y.re; 
     120        z.im = x.im + y.im; 
     121        return z; 
     122} 
     123 
     124complex rcmult(x,y) 
     125        double x; 
     126    complex y; 
     127{ 
     128        complex z; 
     129        z.re = x*y.re; 
     130        z.im = x*y.im; 
     131        return z; 
     132} 
     133 
     134complex cplx_sub(x,y) 
     135        complex x,y; 
     136{ 
     137        complex z; 
     138        z.re = x.re - y.re; 
     139        z.im = x.im - y.im; 
     140        return z; 
     141} 
     142 
     143 
     144complex cplx_mult(x,y) 
     145        complex x,y; 
     146{ 
     147        complex z; 
     148        z.re = x.re*y.re - x.im*y.im; 
     149        z.im = x.re*y.im + x.im*y.re; 
     150        return z; 
     151} 
     152 
     153complex cplx_div(x,y) 
     154        complex x,y; 
     155{ 
     156        complex z; 
     157        z.re = (x.re*y.re + x.im*y.im)/(y.re*y.re + y.im*y.im); 
     158        z.im = (x.im*y.re - x.re*y.im)/(y.re*y.re + y.im*y.im); 
     159        return z; 
     160} 
     161 
     162complex cplx_exp(b) 
     163        complex b; 
     164{ 
     165        complex z; 
    145166        double br,bi; 
    146167        br=b.re; 
    147168        bi=b.im; 
    148         z->re = exp(br)*cos(bi); 
    149         z->im = exp(br)*sin(bi); 
    150 } 
    151  
    152  
    153 void cplx_sqrt(Cplx *c, Cplx z)    //see Schaum`s Math Handbook p. 22, 6.6 and 6.10 
    154 { 
     169        z.re = exp(br)*cos(bi); 
     170        z.im = exp(br)*sin(bi); 
     171        return z; 
     172} 
     173 
     174 
     175complex cplx_sqrt(z)    //see Schaum`s Math Handbook p. 22, 6.6 and 6.10 
     176        complex z; 
     177{ 
     178        complex c; 
    155179        double zr,zi,x,y,r,w; 
    156180 
    … …  
    160184        if (zr==0.0 && zi==0.0) 
    161185        { 
    162                 c->re=0.0; 
    163                 c->im=0.0; 
    164         } else { 
     186    c.re=0.0; 
     187        c.im=0.0; 
     188    return c; 
     189        } 
     190        else 
     191        { 
    165192                x=fabs(zr); 
    166193                y=fabs(zi); 
    … …  
    169196                        r=y/x; 
    170197                        w=sqrt(x)*sqrt(0.5*(1.0+sqrt(1.0+r*r))); 
    171                 } else { 
     198                } 
     199                else 
     200                { 
    172201                        r=x/y; 
    173202                        w=sqrt(y)*sqrt(0.5*(r+sqrt(1.0+r*r))); 
    … …  
    175204                if (zr >=0.0) 
    176205                { 
    177                         c->re=w; 
    178                         c->im=zi/(2.0*w); 
    179                 } else { 
    180                         c->im=(zi >= 0) ? w : -w; 
    181                         c->re=zi/(2.0*c->im); 
     206                        c.re=w; 
     207                        c.im=zi/(2.0*w); 
    182208                } 
    183         } 
    184 } 
    185  
    186 void cplx_cos(Cplx *z, Cplx b) 
    187 { 
    188         // cos(b) = (e^bi + e^-bi)/2 
    189         //        = (e^b.im e^-i bi.re) + e^-b.im e^i b.re)/2 
    190         //        = (e^b.im cos(-b.re) + e^b.im sin(-b.re) i)/2 + (e^-b.im cos(b.re) + e^-b.im sin(b.re) i)/2 
    191         //        = e^b.im cos(b.re)/2 - e^b.im sin(b.re)/2 i + 1/e^b.im cos(b.re)/2 + 1/e^b.im sin(b.re)/2 i 
    192         //        = (e^b.im + 1/e^b.im)/2 cos(b.re) + (-e^b.im + 1/e^b.im)/2 sin(b.re) i 
    193         //        = cosh(b.im) cos(b.re) - sinh(b.im) sin(b.re) i 
    194         double exp_b_im = exp(b.im); 
    195         z->re = 0.5*(+exp_b_im + 1.0/exp_b_im) * cos(b.re); 
    196         z->im = -0.5*(exp_b_im - 1.0/exp_b_im) * sin(b.re); 
     209                else 
     210                { 
     211                        c.im=(zi >= 0) ? w : -w; 
     212                        c.re=zi/(2.0*c.im); 
     213                } 
     214                return c; 
     215        } 
     216} 
     217 
     218complex cplx_cos(b) 
     219        complex b; 
     220{ 
     221        complex zero,two,z,i,bi,negbi; 
     222        zero = cassign(0.0,0.0); 
     223        two = cassign(2.0,0.0); 
     224        i = cassign(0.0,1.0); 
     225        bi = cplx_mult(b,i); 
     226        negbi = cplx_sub(zero,bi); 
     227        z = cplx_div(cplx_add(cplx_exp(bi),cplx_exp(negbi)),two); 
     228        return z; 
    197229} 
    198230 

  • src/sas/sascalc/calculator/c_extensions/librefl.h

    rf54e82cf r9e531f2  
    55        double re; 
    66        double im; 
    7 } Cplx; 
     7} complex; 
    88 
    99typedef struct { 
    10         Cplx a; 
    11         Cplx b; 
    12         Cplx c; 
    13         Cplx d; 
     10        complex a; 
     11        complex b; 
     12        complex c; 
     13        complex d; 
    1414} matrix; 
    1515 
    16 void cassign(Cplx*, double real, double imag); 
     16complex cassign(double real, double imag); 
    1717 
    18 void cplx_add(Cplx*, Cplx x,Cplx y); 
     18complex cplx_add(complex x,complex y); 
    1919 
    20 void rcmult(Cplx*, double x,Cplx y); 
     20complex rcmult(double x,complex y); 
    2121 
    22 void cplx_sub(Cplx*, Cplx x,Cplx y); 
     22complex cplx_sub(complex x,complex y); 
    2323 
    24 void cplx_mult(Cplx*, Cplx x,Cplx y); 
     24complex cplx_mult(complex x,complex y); 
    2525 
    26 void cplx_div(Cplx*, Cplx x,Cplx y); 
     26complex cplx_div(complex x,complex y); 
    2727 
    28 void cplx_exp(Cplx*, Cplx b); 
     28complex cplx_exp(complex b); 
    2929 
    30 void cplx_sqrt(Cplx*, Cplx z); 
     30complex cplx_sqrt(complex z); 
    3131 
    32 void cplx_cos(Cplx*, Cplx b); 
     32complex cplx_cos(complex b); 
    3333 
    3434double intersldfunc(int fun_type, double n_sub, double i, double nu, double sld_l, double sld_r); 

  • src/sas/sascalc/calculator/c_extensions/sld2i.c

    rf54e82cf re6f2009  
    5353        polar_sld b_sld; 
    5454        double qr = 0.0; 
    55         Cplx iqr; 
    56         Cplx ephase; 
    57         Cplx comp_sld; 
    58  
    59         Cplx sumj_uu; 
    60         Cplx sumj_ud; 
    61         Cplx sumj_du; 
    62         Cplx sumj_dd; 
    63         Cplx temp_fi; 
    64  
    65         int i, j; 
     55        complex iqr = cassign(0.0, 0.0); 
     56        complex ephase = cassign(0.0, 0.0); 
     57        complex comp_sld = cassign(0.0, 0.0); 
     58 
     59        complex sumj_uu; 
     60        complex sumj_ud; 
     61        complex sumj_du; 
     62        complex sumj_dd; 
     63        complex temp_fi; 
    6664 
    6765        double count = 0.0; 
    6866        //check if this computation is for averaging 
    69  
    70         cassign(&iqr, 0.0, 0.0); 
    71         cassign(&ephase, 0.0, 0.0); 
    72         cassign(&comp_sld, 0.0, 0.0); 
    7367 
    7468        //Assume that pixel volumes are given in vol_pix in A^3 unit 
    … …  
    7872 
    7973        // Loop over q-values and multiply apply matrix 
     74        int i; 
    8075        for(i=0; i<npoints; i++){ 
    8176                //I_out[i] = 0.0; 
    82                 cassign(&sumj_uu, 0.0, 0.0); 
    83                 cassign(&sumj_ud, 0.0, 0.0); 
    84                 cassign(&sumj_du, 0.0, 0.0); 
    85                 cassign(&sumj_dd, 0.0, 0.0); 
    86                 //printf("i: %d\n", i); 
     77                sumj_uu = cassign(0.0, 0.0); 
     78                sumj_ud = cassign(0.0, 0.0); 
     79                sumj_du = cassign(0.0, 0.0); 
     80                sumj_dd = cassign(0.0, 0.0); 
     81                //printf ("%d ", i); 
    8782                //q = sqrt(qx[i]*qx[i] + qy[i]*qy[i]); // + qz[i]*qz[i]); 
     83                int j; 
    8884                for(j=0; j<this->n_pix; j++){ 
    8985                        if (this->sldn_val[j]!=0.0 
    … …  
    9288                                ||this->mz_val[j]!=0.0) 
    9389                        { 
    94                             // printf("i,j: %d,%d\n", i,j); 
    9590                                //anisotropic 
    96                                 cassign(&temp_fi, 0.0, 0.0); 
    97                                 cal_msld(&b_sld, 0, qx[i], qy[i], this->sldn_val[j], 
     91                                temp_fi = cassign(0.0, 0.0); 
     92                                b_sld = cal_msld(0, qx[i], qy[i], this->sldn_val[j], 
    9893                                                         this->mx_val[j], this->my_val[j], this->mz_val[j], 
    9994                                                         this->inspin, this->outspin, this->stheta); 
    10095                                qr = (qx[i]*this->x_val[j] + qy[i]*this->y_val[j]); 
    101                                 cassign(&iqr, 0.0, qr); 
    102                                 cplx_exp(&ephase, iqr); 
     96                                iqr = cassign(0.0, qr); 
     97                                ephase = cplx_exp(iqr); 
    10398 
    10499                                //Let's multiply pixel(atomic) volume here 
    105                                 rcmult(&ephase, this->vol_pix[j], ephase); 
     100                                ephase = rcmult(this->vol_pix[j], ephase); 
    106101                                //up_up 
    107102                                if (this->inspin > 0.0 && this->outspin > 0.0){ 
    108                                         cassign(&comp_sld, b_sld.uu, 0.0); 
    109                                         cplx_mult(&temp_fi, comp_sld, ephase); 
    110                                         cplx_add(&sumj_uu, sumj_uu, temp_fi); 
     103                                        comp_sld = cassign(b_sld.uu, 0.0); 
     104                                        temp_fi = cplx_mult(comp_sld, ephase); 
     105                                        sumj_uu = cplx_add(sumj_uu, temp_fi); 
    111106                                } 
    112107                                //down_down 
    113108                                if (this->inspin < 1.0 && this->outspin < 1.0){ 
    114                                         cassign(&comp_sld, b_sld.dd, 0.0); 
    115                                         cplx_mult(&temp_fi, comp_sld, ephase); 
    116                                         cplx_add(&sumj_dd, sumj_dd, temp_fi); 
     109                                        comp_sld = cassign(b_sld.dd, 0.0); 
     110                                        temp_fi = cplx_mult(comp_sld, ephase); 
     111                                        sumj_dd = cplx_add(sumj_dd, temp_fi); 
    117112                                } 
    118113                                //up_down 
    119114                                if (this->inspin > 0.0 && this->outspin < 1.0){ 
    120                                         cassign(&comp_sld, b_sld.re_ud, b_sld.im_ud); 
    121                                         cplx_mult(&temp_fi, comp_sld, ephase); 
    122                                         cplx_add(&sumj_ud, sumj_ud, temp_fi); 
     115                                        comp_sld = cassign(b_sld.re_ud, b_sld.im_ud); 
     116                                        temp_fi = cplx_mult(comp_sld, ephase); 
     117                                        sumj_ud = cplx_add(sumj_ud, temp_fi); 
    123118                                } 
    124119                                //down_up 
    125120                                if (this->inspin < 1.0 && this->outspin > 0.0){ 
    126                                         cassign(&comp_sld, b_sld.re_du, b_sld.im_du); 
    127                                         cplx_mult(&temp_fi, comp_sld, ephase); 
    128                                         cplx_add(&sumj_du, sumj_du, temp_fi); 
    129                                 } 
     121                                        comp_sld = cassign(b_sld.re_du, b_sld.im_du); 
     122                                        temp_fi = cplx_mult(comp_sld, ephase); 
     123                                        sumj_du = cplx_add(sumj_du, temp_fi); 
     124                                } 
     125 
    130126 
    131127                                if (i == 0){ 
    … …  
    162158        //Assume that pixel volumes are given in vol_pix in A^3 unit 
    163159        // Loop over q-values and multiply apply matrix 
    164     int i, j, k; 
     160        int i;    
    165161        for(i=0; i<npoints; i++){ 
    166162                sumj =0.0; 
     163                int j; 
    167164                for(j=0; j<n_pix; j++){ 
    168165                        //Isotropic: Assumes all slds are real (no magnetic) 
    … …  
    182179                                //full calculation 
    183180                                //pragma omp parallel for 
     181                                int k; 
    184182                                for(k=0; k<n_pix; k++){ 
    185183                                        sld_j =  this->sldn_val[j] * this->sldn_val[k] * this->vol_pix[j] * this->vol_pix[k]; 

  • src/sas/sascalc/calculator/c_extensions/sld2i_module.c

    rf54e82cf r3010f68  
    44#include <Python.h> 
    55#include <stdio.h> 
    6 #include "sld2i.h" 
     6#include <sld2i.h> 
    77 
    88#if PY_MAJOR_VERSION < 3 
    … …  
    6868        double outspin; 
    6969        double stheta; 
    70         GenI *sld2i; 
    7170 
    7271        if (!PyArg_ParseTuple(args, "iOOOOOOOOddd", &n_pix, &x_val_obj, &y_val_obj, &z_val_obj, &sldn_val_obj, &mx_val_obj, &my_val_obj, &mz_val_obj, &vol_pix_obj, &inspin, &outspin, &stheta)) return NULL; 
    … …  
    7978        OUTVECTOR(mz_val_obj, mz_val, n_x); 
    8079        OUTVECTOR(vol_pix_obj, vol_pix, n_x); 
    81         sld2i =  PyMem_Malloc(sizeof(GenI)); 
     80        GenI* sld2i =  PyMem_Malloc(sizeof(GenI)); 
    8281        if (sld2i != NULL) { 
    8382                initGenI(sld2i, n_pix,x_val,y_val,z_val,sldn_val,mx_val,my_val,mz_val,vol_pix,inspin,outspin,stheta); 
    … …  
    9998        double *I_out; 
    10099        PyObject *gen_obj; 
    101         GenI *sld2i; 
    102100 
    103101        if (!PyArg_ParseTuple(args, "OiOOO",  &gen_obj, &npoints, &qx_obj, &qy_obj, &I_out_obj)) return NULL; 
    … …  
    110108 
    111109        // Set the array pointers 
    112         sld2i = (GenI *)PyCapsule_GetPointer(gen_obj, "GenI"); 
     110        GenI* sld2i = (GenI *)PyCapsule_GetPointer(gen_obj, "GenI"); 
    113111 
    114112        genicomXY(sld2i, npoints, qx, qy, I_out); 
    … …  
    128126        double *I_out; 
    129127        PyObject *gen_obj; 
    130         GenI *sld2i; 
    131128 
    132129        if (!PyArg_ParseTuple(args, "OiOO",  &gen_obj, &npoints, &q_obj, &I_out_obj)) return NULL; 
    … …  
    138135 
    139136        // Set the array pointers 
    140         sld2i = (GenI *)PyCapsule_GetPointer(gen_obj, "GenI"); 
     137        GenI *sld2i = (GenI *)PyCapsule_GetPointer(gen_obj, "GenI"); 
    141138 
    142139        genicom(sld2i, npoints, q, I_out); 
Note: See TracChangeset for help on using the changeset viewer.