schain Files · schainpy/model/proc/extensions.c

Update plot codes, add meta attribute to dataOut to send metadata to plots

José Chávez - - Load All Authors

File last commit:

r1004:3c7daf19daec


                r1139:ade57c0ecace

Download file

             extensions.c
        
                    139 lines
            
             | 4.3 KiB
            
                | text/x-c
            
             |
                CLexer
            
             / schainpy / model / proc / extensions.c
          
                    History
                
                 |
                  Source
                 | Raw
                 |Copy content
                 |Copy permalink

        José Chávez
    
cambios xmax xmin

              r1004
            
      #define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION

      #define NUM_CPY_THREADS 8

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
      #include <Python.h>

      #include <numpy/arrayobject.h>

      #include <math.h>

        José Chávez
    
cambios xmax xmin

              r1004
            
      #include <complex.h>

      #include <time.h>

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
        José Chávez
    
cambios xmax xmin

              r1004
            
      // void printArr(int *array);

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
      static PyObject *hildebrand_sekhon(PyObject *self, PyObject *args);

        José Chávez
    
cambios xmax xmin

              r1004
            
      static PyObject *correlateByBlock(PyObject *self, PyObject *args);

      #ifndef PyMODINIT_FUNC	/* declarations for DLL import/export */

      #define PyMODINIT_FUNC void

      #endif

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
      static PyMethodDef extensionsMethods[] = {

        José Chávez
    
cambios xmax xmin

              r1004
            
        { "correlateByBlock", (PyCFunction)correlateByBlock, METH_VARARGS, "get correlation by block" },

        { "hildebrand_sekhon", (PyCFunction)hildebrand_sekhon, METH_VARARGS, "get noise with hildebrand_sekhon" },

        { NULL, NULL, 0, NULL }

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
      };

      PyMODINIT_FUNC initcSchain() {

        Py_InitModule("cSchain", extensionsMethods);

        import_array();

      }

        José Chávez
    
cambios xmax xmin

              r1004
            
      static PyObject *correlateByBlock(PyObject *self, PyObject *args) {

        // int *x  = (int*) malloc(4000000 * 216 * sizeof(int));;

        // int a = 5;

        // x = &a;

        // int b = 6;

        // x = &b;

        // printf("Antes de imprimir x \n");

        // printf("%d \n", x[0]);

        PyObject *data_obj1, *data_obj2;

        PyArrayObject *data_array1, *data_array2, *correlateRow, *out, *dataRow, *codeRow; //, , 

        int mode;

        if (!PyArg_ParseTuple(args, "OOi", &data_obj1, &data_obj2, &mode)) return NULL;

        data_array1 = (PyArrayObject *) PyArray_FROM_OTF(data_obj1, NPY_COMPLEX128, NPY_ARRAY_DEFAULT);

        data_array2 = (PyArrayObject *) PyArray_FROM_OTF(data_obj2, NPY_FLOAT64, NPY_ARRAY_DEFAULT);

        npy_intp dims[1];

        dims[0] = 200;

        npy_intp dims_code[1];

        dims_code[0] = 16;

        double complex * dataRaw;

        double * codeRaw;

        dataRaw = (double complex*) PyArray_DATA(data_array1);

        codeRaw = (double *) PyArray_DATA(data_array2);

        double complex ** outC = malloc(40000*200*sizeof(double complex));

        int i;

        clock_t start = clock();

        for(i=0; i<40000; i++){

          // codeRow = PyArray_SimpleNewFromData(1, dims_code, NPY_FLOAT64, codeRaw + 16 * i);

          // dataRow = PyArray_SimpleNewFromData(1, dims, NPY_COMPLEX128, dataRaw + 200 * i);

          // Py_INCREF(codeRow);

          // Py_INCREF(dataRow);

          // PyArray_ENABLEFLAGS(codeRow, NPY_ARRAY_OWNDATA);

          // PyArray_ENABLEFLAGS(dataRow, NPY_ARRAY_OWNDATA);

          correlateRow = (PyArrayObject *) PyArray_Correlate2(PyArray_SimpleNewFromData(1, dims_code, NPY_FLOAT64, codeRaw + 16 * i), PyArray_SimpleNewFromData(1, dims, NPY_COMPLEX128, dataRaw + 200 * i), (npy_intp) 2);

          //Py_INCREF(correlateRow);

          // PyArray_ENABLEFLAGS(correlateRow, NPY_ARRAY_OWNDATA);

          memcpy(outC + 200*i, (double complex*) PyArray_DATA(correlateRow), 200 * sizeof(double complex));

          Py_DECREF(correlateRow);

          // Py_DECREF(codeRow);

          // Py_DECREF(dataRow);

        }

        clock_t end = clock();

        float seconds = (float)(end - start) / CLOCKS_PER_SEC;

        printf("%f", seconds);

        // 

        npy_intp dimsret[2];

        dimsret[0] = 40000;

        dimsret[1] = 200;

        out = PyArray_SimpleNewFromData(2, dimsret, NPY_COMPLEX128, outC);

        PyArray_ENABLEFLAGS(out, NPY_ARRAY_OWNDATA);

        //Py_INCREF(out);

        Py_DECREF(data_array1);

        Py_DECREF(data_array2);

        // PyArray_DebugPrint(out);

        // Py_DECREF(data_obj2);

        // Py_DECREF(data_obj1);

        // Py_DECREF(codeRow);

        // Py_DECREF(dataRow);

        // free(dataRaw);

        // free(codeRaw);

        return PyArray_Return(out);

      }

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
      static PyObject *hildebrand_sekhon(PyObject *self, PyObject *args) {

        double navg;

        PyObject *data_obj, *data_array;

        if (!PyArg_ParseTuple(args, "Od", &data_obj, &navg)) return NULL;

        José Chávez
    
cambios xmax xmin

              r1004
            
        data_array = PyArray_FROM_OTF(data_obj, NPY_FLOAT64, NPY_ARRAY_DEFAULT);

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
        if (data_array == NULL) {

            Py_XDECREF(data_array);

            Py_XDECREF(data_obj);

            return NULL;

        }

        double *sortdata = (double*)PyArray_DATA(data_array);

        int lenOfData = (int)PyArray_SIZE(data_array) ;

        double nums_min = lenOfData*0.2;

        if (nums_min <= 5) nums_min = 5;

        double sump = 0;

        double sumq = 0;

        int j = 0;

        int cont = 1;

        double rtest = 0;

        while ((cont == 1) && (j < lenOfData)) {

          sump = sump + sortdata[j];

          sumq = sumq + pow(sortdata[j], 2);

          if (j > nums_min) {

            rtest = (double)j/(j-1) + 1/navg;

            if ((sumq*j) > (rtest*pow(sump, 2))) {

              j = j - 1;

              sump = sump - sortdata[j];

              sumq = sumq - pow(sortdata[j],2);

              cont = 0;

            }

          }

          j = j + 1;

        }

        double lnoise = sump / j;

        Juan C. Valdez
    
Fix memory leak in c extension

              r880
            
        Py_DECREF(data_array);

        Juan C. Valdez
    
Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup

              r878
            
        return Py_BuildValue("d", lnoise);

      }

        José Chávez
    
cambios xmax xmin

              r1004

	Site-wide shortcuts
/	Use quick search box
g h	Goto home page
g g	Goto my private gists page
g G	Goto my public gists page
g 0-9	Goto bookmarked items from 0-9
n r	New repository page
n g	New gist page

	Repositories
g s	Goto summary page
g c	Goto changelog page
g f	Goto files page
g F	Goto files page with file search activated
g p	Goto pull requests page
g o	Goto repository settings
g O	Goto repository access permissions settings
t s	Toggle sidebar on some pages

José Chávez cambios xmax xmin	r1004	#define NPY_NO_DEPRECATED_API NPY_1_7_API_VERSION
		#define NUM_CPY_THREADS 8
Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878	#include <Python.h>
		#include <numpy/arrayobject.h>
		#include <math.h>
José Chávez cambios xmax xmin	r1004	#include <complex.h>
		#include <time.h>
Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878
José Chávez cambios xmax xmin	r1004	// void printArr(int *array);
Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878	static PyObject hildebrand_sekhon(PyObject self, PyObject *args);
José Chávez cambios xmax xmin	r1004	static PyObject correlateByBlock(PyObject self, PyObject *args);
		#ifndef PyMODINIT_FUNC /* declarations for DLL import/export */
		#define PyMODINIT_FUNC void
		#endif
Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878
		static PyMethodDef extensionsMethods[] = {
José Chávez cambios xmax xmin	r1004	{ "correlateByBlock", (PyCFunction)correlateByBlock, METH_VARARGS, "get correlation by block" },
		{ "hildebrand_sekhon", (PyCFunction)hildebrand_sekhon, METH_VARARGS, "get noise with hildebrand_sekhon" },
		{ NULL, NULL, 0, NULL }
Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878	};

		PyMODINIT_FUNC initcSchain() {
		Py_InitModule("cSchain", extensionsMethods);
		import_array();
		}

José Chávez cambios xmax xmin	r1004	static PyObject correlateByBlock(PyObject self, PyObject *args) {

		// int x = (int) malloc(4000000 * 216 * sizeof(int));;
		// int a = 5;
		// x = &a;
		// int b = 6;
		// x = &b;
		// printf("Antes de imprimir x \n");
		// printf("%d \n", x[0]);

		PyObject data_obj1, data_obj2;
		PyArrayObject data_array1, data_array2, correlateRow, out, dataRow, codeRow; //, ,
		int mode;

		if (!PyArg_ParseTuple(args, "OOi", &data_obj1, &data_obj2, &mode)) return NULL;

		data_array1 = (PyArrayObject *) PyArray_FROM_OTF(data_obj1, NPY_COMPLEX128, NPY_ARRAY_DEFAULT);
		data_array2 = (PyArrayObject *) PyArray_FROM_OTF(data_obj2, NPY_FLOAT64, NPY_ARRAY_DEFAULT);

		npy_intp dims[1];
		dims[0] = 200;
		npy_intp dims_code[1];
		dims_code[0] = 16;

		double complex * dataRaw;
		double * codeRaw;
		dataRaw = (double complex*) PyArray_DATA(data_array1);
		codeRaw = (double *) PyArray_DATA(data_array2);
		double complex ** outC = malloc(40000200sizeof(double complex));
		int i;

		clock_t start = clock();
		for(i=0; i<40000; i++){
		// codeRow = PyArray_SimpleNewFromData(1, dims_code, NPY_FLOAT64, codeRaw + 16 * i);
		// dataRow = PyArray_SimpleNewFromData(1, dims, NPY_COMPLEX128, dataRaw + 200 * i);
		// Py_INCREF(codeRow);
		// Py_INCREF(dataRow);
		// PyArray_ENABLEFLAGS(codeRow, NPY_ARRAY_OWNDATA);
		// PyArray_ENABLEFLAGS(dataRow, NPY_ARRAY_OWNDATA);
		correlateRow = (PyArrayObject ) PyArray_Correlate2(PyArray_SimpleNewFromData(1, dims_code, NPY_FLOAT64, codeRaw + 16 i), PyArray_SimpleNewFromData(1, dims, NPY_COMPLEX128, dataRaw + 200 * i), (npy_intp) 2);
		//Py_INCREF(correlateRow);
		// PyArray_ENABLEFLAGS(correlateRow, NPY_ARRAY_OWNDATA);
		memcpy(outC + 200i, (double complex) PyArray_DATA(correlateRow), 200 * sizeof(double complex));

		Py_DECREF(correlateRow);
		// Py_DECREF(codeRow);
		// Py_DECREF(dataRow);
		}
		clock_t end = clock();
		float seconds = (float)(end - start) / CLOCKS_PER_SEC;
		printf("%f", seconds);
		//
		npy_intp dimsret[2];
		dimsret[0] = 40000;
		dimsret[1] = 200;
		out = PyArray_SimpleNewFromData(2, dimsret, NPY_COMPLEX128, outC);
		PyArray_ENABLEFLAGS(out, NPY_ARRAY_OWNDATA);
		//Py_INCREF(out);
		Py_DECREF(data_array1);
		Py_DECREF(data_array2);
		// PyArray_DebugPrint(out);
		// Py_DECREF(data_obj2);
		// Py_DECREF(data_obj1);
		// Py_DECREF(codeRow);
		// Py_DECREF(dataRow);
		// free(dataRaw);
		// free(codeRaw);

		return PyArray_Return(out);
		}

Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878	static PyObject hildebrand_sekhon(PyObject self, PyObject *args) {
		double navg;
		PyObject data_obj, data_array;

		if (!PyArg_ParseTuple(args, "Od", &data_obj, &navg)) return NULL;
José Chávez cambios xmax xmin	r1004	data_array = PyArray_FROM_OTF(data_obj, NPY_FLOAT64, NPY_ARRAY_DEFAULT);
Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878	if (data_array == NULL) {
		Py_XDECREF(data_array);
		Py_XDECREF(data_obj);
		return NULL;
		}
		double sortdata = (double)PyArray_DATA(data_array);
		int lenOfData = (int)PyArray_SIZE(data_array) ;
		double nums_min = lenOfData*0.2;
		if (nums_min <= 5) nums_min = 5;
		double sump = 0;
		double sumq = 0;
		int j = 0;
		int cont = 1;
		double rtest = 0;
		while ((cont == 1) && (j < lenOfData)) {
		sump = sump + sortdata[j];
		sumq = sumq + pow(sortdata[j], 2);
		if (j > nums_min) {
		rtest = (double)j/(j-1) + 1/navg;
		if ((sumqj) > (rtestpow(sump, 2))) {
		j = j - 1;
		sump = sump - sortdata[j];
		sumq = sumq - pow(sortdata[j],2);
		cont = 0;
		}
		}
		j = j + 1;
		}

		double lnoise = sump / j;
Juan C. Valdez Fix memory leak in c extension	r880
		Py_DECREF(data_array);

Juan C. Valdez Add C extension for hildebrand_sekhon noise estimation and update requirements versions in setup	r878	return Py_BuildValue("d", lnoise);
		}
José Chávez cambios xmax xmin	r1004