Changeset 0957bb3a in sasview

Timestamp:

Nov 7, 2017 10:00:25 AM (7 years ago)

Author:

Paul Kienzle <pkienzle@…>

Branches:

master, magnetic_scatt, release-4.2.2, ticket-1009, ticket-1094-headless, ticket-1242-2d-resolution, ticket-1243, ticket-1249, ticket885, unittest-saveload

Children:

54b0650

Parents:

0fc5a03

Message:

convert sldi from C++ to C

Files:

• setup.py (modified) (1 diff)
• src/sas/sascalc/calculator/c_extensions/sld2i.c (moved) (moved from src/sas/sascalc/calculator/c_extensions/sld2i.cpp) (13 diffs)
• src/sas/sascalc/calculator/c_extensions/sld2i.h (added)
• src/sas/sascalc/calculator/c_extensions/sld2i.hh (deleted)
• src/sas/sascalc/calculator/c_extensions/sld2i_module.c (moved) (moved from src/sas/sascalc/calculator/c_extensions/sld2i_module.cpp) (8 diffs)

setup.py

@@ -242 +242 @@
 ext_modules.append(Extension("sas.sascalc.calculator.core.sld2i",
                              sources=[
-                                 os.path.join(gen_dir, "sld2i_module.cpp"),
-                                 os.path.join(gen_dir, "sld2i.cpp"),
+                                 os.path.join(gen_dir, "sld2i_module.c"),
+                                 os.path.join(gen_dir, "sld2i.c"),
                                  os.path.join(gen_dir, "libfunc.c"),
                                  os.path.join(gen_dir, "librefl.c"),

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

@@ -2 +2 @@
 Computes the (magnetic) scattering form sld (n and m) profile
  */
-#include "sld2i.hh"
 #include <stdio.h>
 #include <math.h>
-using namespace std;
-extern "C" {
-        #include "libfunc.h"
-        #include "librefl.h"
-}
+#include "sld2i.h"
+#include "libfunc.h"
+#include "librefl.h"
 /**
  * Constructor for GenI
@@ -28 +25 @@
  * @param s_theta: angle (from x-axis) of the up spin in degree
  */
-GenI :: GenI(int npix, double* x, double* y, double* z, double* sldn,
+void initGenI(GenI* this, int npix, double* x, double* y, double* z, double* sldn,
                         double* mx, double* my, double* mz, double* voli,
                         double in_spin, double out_spin,
                         double s_theta) {
-        //this->qx_val = qx;
-        //this->qy_val = qy;
-        //this->qz_val = qz;
         this->n_pix = npix;
         this->x_val = x;
@@ -47 +41 @@
         this->outspin = out_spin;
         this->stheta = s_theta;
-};
+}
 
 /**
  * Compute 2D anisotropic
  */
-void GenI :: genicomXY(int npoints, double *qx, double *qy, double *I_out){
-        //npoints is given negative for angular averaging 
+void genicomXY(GenI* this, int npoints, double *qx, double *qy, double *I_out){
+        //npoints is given negative for angular averaging
         // Assumes that q doesn't have qz component and sld_n is all real
         //double q = 0.0;
@@ -76 +70 @@
         //int y_size = 0; //in Ang
         //int z_size = 0; //in Ang
-        
+
         // Loop over q-values and multiply apply matrix
-        
+
         for(int i=0; i<npoints; i++){
                 //I_out[i] = 0.0;
@@ -84 +78 @@
                 sumj_ud = cassign(0.0, 0.0);
                 sumj_du = cassign(0.0, 0.0);
-                sumj_dd = cassign(0.0, 0.0);            
+                sumj_dd = cassign(0.0, 0.0);
                 //printf ("%d ", i);
                 //q = sqrt(qx[i]*qx[i] + qy[i]*qy[i]); // + qz[i]*qz[i]);
 
-                for(int j=0; j<n_pix; j++){
-                        if (sldn_val[j]!=0.0||mx_val[j]!=0.0||my_val[j]!=0.0||mz_val[j]!=0.0)
-                        {       
+                for(int j=0; j<this->n_pix; j++){
+                        if (this->sldn_val[j]!=0.0
+                                ||this->mx_val[j]!=0.0
+                                ||this->my_val[j]!=0.0
+                                ||this->mz_val[j]!=0.0)
+                        {
                                 //anisotropic
                                 temp_fi = cassign(0.0, 0.0);
-                                b_sld = cal_msld(0, qx[i], qy[i], sldn_val[j],
-                                                         mx_val[j], my_val[j], mz_val[j],
-                                                         inspin, outspin, stheta);
-                                qr = (qx[i]*x_val[j] + qy[i]*y_val[j]);
+                                b_sld = cal_msld(0, qx[i], qy[i], this->sldn_val[j],
+                                                         this->mx_val[j], this->my_val[j], this->mz_val[j],
+                                                         this->inspin, this->outspin, this->stheta);
+                                qr = (qx[i]*this->x_val[j] + qy[i]*this->y_val[j]);
                                 iqr = cassign(0.0, qr);
                                 ephase = cplx_exp(iqr);
-                                
+
                                 //Let's multiply pixel(atomic) volume here
-                                ephase = rcmult(vol_pix[j], ephase);
+                                ephase = rcmult(this->vol_pix[j], ephase);
                                 //up_up
-                                if (inspin > 0.0 && outspin > 0.0){
+                                if (this->inspin > 0.0 && this->outspin > 0.0){
                                         comp_sld = cassign(b_sld.uu, 0.0);
                                         temp_fi = cplx_mult(comp_sld, ephase);
@@ -109 +106 @@
                                 }
                                 //down_down
-                                if (inspin < 1.0 && outspin < 1.0){
+                                if (this->inspin < 1.0 && this->outspin < 1.0){
                                         comp_sld = cassign(b_sld.dd, 0.0);
                                         temp_fi = cplx_mult(comp_sld, ephase);
@@ -115 +112 @@
                                 }
                                 //up_down
-                                if (inspin > 0.0 && outspin < 1.0){
+                                if (this->inspin > 0.0 && this->outspin < 1.0){
                                         comp_sld = cassign(b_sld.re_ud, b_sld.im_ud);
                                         temp_fi = cplx_mult(comp_sld, ephase);
@@ -121 +118 @@
                                 }
                                 //down_up
-                                if (inspin < 1.0 && outspin > 0.0){
+                                if (this->inspin < 1.0 && this->outspin > 0.0){
                                         comp_sld = cassign(b_sld.re_du, b_sld.im_du);
                                         temp_fi = cplx_mult(comp_sld, ephase);
@@ -129 +126 @@
 
                                 if (i == 0){
-                                        count += vol_pix[j];
+                                        count += this->vol_pix[j];
                                 }
                         }
@@ -149 +146 @@
  * Also assumes there is no polarization: No dependency on spin
  */
-void GenI :: genicom(int npoints, double *q, double *I_out){
-        //npoints is given negative for angular averaging 
+void genicom(GenI* this, int npoints, double *q, double *I_out){
+        //npoints is given negative for angular averaging
         // Assumes that q doesn't have qz component and sld_n is all real
         //double Pi = 4.0*atan(1.0);
-        int is_sym = 0;
+        int is_sym = this->n_pix < 0;
         double qr = 0.0;
         double sumj;
         double sld_j = 0.0;
         double count = 0.0;
-        if (n_pix < 0 ){
-                is_sym = 1;
-                n_pix = n_pix * -1;
-        }
+        int n_pix = is_sym ? -this->n_pix : this->n_pix;
         //Assume that pixel volumes are given in vol_pix in A^3 unit
         // Loop over q-values and multiply apply matrix
         for(int i=0; i<npoints; i++){
-                sumj =0.0;              
+                sumj =0.0;
                 for(int j=0; j<n_pix; j++){
                         //Isotropic: Assumes all slds are real (no magnetic)
@@ -171 +165 @@
                         if (is_sym == 1){
                                 // approximation for a spherical symmetric particle
-                                qr = sqrt(x_val[j]*x_val[j]+y_val[j]*y_val[j]+z_val[j]*z_val[j])*q[i];
+                                qr = sqrt(this->x_val[j]*this->x_val[j]+this->y_val[j]*this->y_val[j]+this->z_val[j]*this->z_val[j])*q[i];
                                 if (qr > 0.0){
                                         qr = sin(qr) / qr;
-                                        sumj += sldn_val[j] * vol_pix[j] * qr;
+                                        sumj += this->sldn_val[j] * this->vol_pix[j] * qr;
                                 }
                                 else{
-                                        sumj += sldn_val[j] * vol_pix[j];
+                                        sumj += this->sldn_val[j] * this->vol_pix[j];
                                 }
                         }
@@ -184 +178 @@
                                 //pragma omp parallel for
                                 for(int k=0; k<n_pix; k++){
-                                        sld_j =  sldn_val[j] * sldn_val[k] * vol_pix[j] * vol_pix[k];
-                                        qr = (x_val[j]-x_val[k])*(x_val[j]-x_val[k])+
-                                                      (y_val[j]-y_val[k])*(y_val[j]-y_val[k])+ 
-                                                      (z_val[j]-z_val[k])*(z_val[j]-z_val[k]);
+                                        sld_j =  this->sldn_val[j] * this->sldn_val[k] * this->vol_pix[j] * this->vol_pix[k];
+                                        qr = (this->x_val[j]-this->x_val[k])*(this->x_val[j]-this->x_val[k])+
+                                                      (this->y_val[j]-this->y_val[k])*(this->y_val[j]-this->y_val[k])+
+                                                      (this->z_val[j]-this->z_val[k])*(this->z_val[j]-this->z_val[k]);
                                         qr = sqrt(qr) * q[i];
                                         if (qr > 0.0){
@@ -198 +192 @@
                         }
                         if (i == 0){
-                                count += vol_pix[j];
+                                count += this->vol_pix[j];
                         }
                 }

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

@@ -4 +4 @@
 #include <Python.h>
 #include <stdio.h>
-#include <sld2i.hh>
+#include <sld2i.h>
 
 #if PY_MAJOR_VERSION < 3
@@ -35 +35 @@
 void
 del_sld2i(PyObject *obj){
-        GenI* sld2i = static_cast<GenI *>(PyCapsule_GetPointer(obj, "GenI"));
-        delete sld2i;
-        return;
+        GenI* sld2i = (GenI *)(PyCapsule_GetPointer(obj, "GenI"));
+        PyMem_Free((void *)sld2i);
 }
 
@@ -43 +42 @@
  * Create a GenI as a python object by supplying arrays
  */
-PyObject * new_GenI(PyObject *, PyObject *args) {
+PyObject * new_GenI(PyObject *self, PyObject *args) {
         PyObject *x_val_obj;
         PyObject *y_val_obj;
@@ -79 +78 @@
         OUTVECTOR(mz_val_obj, mz_val, n_x);
         OUTVECTOR(vol_pix_obj, vol_pix, n_x);
-        GenI* sld2i = new GenI(n_pix,x_val,y_val,z_val,sldn_val,mx_val,my_val,mz_val,vol_pix,inspin,outspin,stheta);
+        GenI* sld2i =  PyMem_Malloc(sizeof(GenI));
+        if (sld2i != NULL) {
+                initGenI(sld2i, n_pix,x_val,y_val,z_val,sldn_val,mx_val,my_val,mz_val,vol_pix,inspin,outspin,stheta);
+        }
         return PyCapsule_New(sld2i, "GenI", del_sld2i);
 }
@@ -86 +88 @@
  * GenI the given input (2D) according to a given object
  */
-PyObject * genicom_inputXY(PyObject *, PyObject *args) {
+PyObject * genicom_inputXY(PyObject *self, PyObject *args) {
         int npoints;
         PyObject *qx_obj;
@@ -106 +108 @@
 
         // Set the array pointers
-        void *temp = PyCapsule_GetPointer(gen_obj, "GenI");
-        GenI* s = static_cast<GenI *>(temp);
-
-        s->genicomXY(npoints, qx, qy, I_out);
+        GenI* sld2i = (GenI *)PyCapsule_GetPointer(gen_obj, "GenI");
+
+        genicomXY(sld2i, npoints, qx, qy, I_out);
         //return PyCObject_FromVoidPtr(s, del_genicom);
         return Py_BuildValue("i",1);
@@ -117 +118 @@
  * GenI the given 1D input according to a given object
  */
-PyObject * genicom_input(PyObject *, PyObject *args) {
+PyObject * genicom_input(PyObject *self, PyObject *args) {
         int npoints;
         PyObject *q_obj;
@@ -134 +135 @@
 
         // Set the array pointers
-        void *temp = PyCapsule_GetPointer(gen_obj, "GenI");
-        GenI* s = static_cast<GenI *>(temp);
-
-        s->genicom(npoints, q, I_out);
+        GenI *sld2i = (GenI *)PyCapsule_GetPointer(gen_obj, "GenI");
+
+        genicom(sld2i, npoints, q, I_out);
         //return PyCObject_FromVoidPtr(s, del_genicom);
         return Py_BuildValue("i",1);