Changes in src/sas/sascalc/calculator/c_extensions/sld2i.c [e6f2009:f54e82cf] in sasview
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
src/sas/sascalc/calculator/c_extensions/sld2i.c
re6f2009 rf54e82cf 53 53 polar_sld b_sld; 54 54 double qr = 0.0; 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; 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; 64 66 65 67 double count = 0.0; 66 68 //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); 67 73 68 74 //Assume that pixel volumes are given in vol_pix in A^3 unit … … 72 78 73 79 // Loop over q-values and multiply apply matrix 74 int i;75 80 for(i=0; i<npoints; i++){ 76 81 //I_out[i] = 0.0; 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);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); 82 87 //q = sqrt(qx[i]*qx[i] + qy[i]*qy[i]); // + qz[i]*qz[i]); 83 int j;84 88 for(j=0; j<this->n_pix; j++){ 85 89 if (this->sldn_val[j]!=0.0 … … 88 92 ||this->mz_val[j]!=0.0) 89 93 { 94 // printf("i,j: %d,%d\n", i,j); 90 95 //anisotropic 91 temp_fi = cassign(0.0, 0.0);92 b_sld = cal_msld(0, qx[i], qy[i], this->sldn_val[j],96 cassign(&temp_fi, 0.0, 0.0); 97 cal_msld(&b_sld, 0, qx[i], qy[i], this->sldn_val[j], 93 98 this->mx_val[j], this->my_val[j], this->mz_val[j], 94 99 this->inspin, this->outspin, this->stheta); 95 100 qr = (qx[i]*this->x_val[j] + qy[i]*this->y_val[j]); 96 iqr = cassign(0.0, qr);97 ephase = cplx_exp(iqr);101 cassign(&iqr, 0.0, qr); 102 cplx_exp(&ephase, iqr); 98 103 99 104 //Let's multiply pixel(atomic) volume here 100 ephase = rcmult(this->vol_pix[j], ephase);105 rcmult(&ephase, this->vol_pix[j], ephase); 101 106 //up_up 102 107 if (this->inspin > 0.0 && this->outspin > 0.0){ 103 c omp_sld = cassign(b_sld.uu, 0.0);104 temp_fi = cplx_mult(comp_sld, ephase);105 sumj_uu = cplx_add(sumj_uu, temp_fi);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); 106 111 } 107 112 //down_down 108 113 if (this->inspin < 1.0 && this->outspin < 1.0){ 109 c omp_sld = cassign(b_sld.dd, 0.0);110 temp_fi = cplx_mult(comp_sld, ephase);111 sumj_dd = cplx_add(sumj_dd, temp_fi);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); 112 117 } 113 118 //up_down 114 119 if (this->inspin > 0.0 && this->outspin < 1.0){ 115 c omp_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);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); 118 123 } 119 124 //down_up 120 125 if (this->inspin < 1.0 && this->outspin > 0.0){ 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 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 } 126 130 127 131 if (i == 0){ … … 158 162 //Assume that pixel volumes are given in vol_pix in A^3 unit 159 163 // Loop over q-values and multiply apply matrix 160 int i;164 int i, j, k; 161 165 for(i=0; i<npoints; i++){ 162 166 sumj =0.0; 163 int j;164 167 for(j=0; j<n_pix; j++){ 165 168 //Isotropic: Assumes all slds are real (no magnetic) … … 179 182 //full calculation 180 183 //pragma omp parallel for 181 int k;182 184 for(k=0; k<n_pix; k++){ 183 185 sld_j = this->sldn_val[j] * this->sldn_val[k] * this->vol_pix[j] * this->vol_pix[k];
Note: See TracChangeset
for help on using the changeset viewer.