Anonymous avatar Anonymous committed fdcd1bc

Initial version. On the way to support numarray

Comments (0)

Files changed (2)

+"""
+Find out if Numeric of numarray is going to be used as its array package
+
+WARNING: PyGSL has used Numeric as its core. This is to test to see
+if we can support both, but sooner or later the favour has to be given
+to one of these packages.
+
+When imported this module:
+     1.) Looks on the command line if it can find a flag of type
+         --array-object=[Numeric|nummarray]
+
+     2.) Tries to import these modules.
+
+     3.) Tries to use the preferred one.
+
+     4.) Or goes with the other one found
+     
+         
+"""
+import re
+import sys
+import os
+from distutils.errors import DistutilsModuleError
+
+
+array_pattern = re.compile("--array-object=(.+)")
+pos=0
+array_preference = None
+while pos<len(sys.argv):
+	match=array_pattern.match(sys.argv[pos])
+	if match:
+            array_preference=match.group(1)
+            array_preference.strip()
+            sys.argv[pos:pos+1]=[]
+            break
+        pos+=1
+
+
+
+have_numeric = 0
+have_numarray = 0
+try:
+	import Numeric
+	have_numeric = 1
+except ImportError:
+	pass
+
+try:
+	import numarray
+	have_numarray = 1
+except ImportError:
+	pass
+
+
+use_numeric = 0
+use_numarray = 0
+
+if array_preference != None:
+    if array_preference == 'Numeric':
+        if have_numeric == 1:
+            use_numeric = 1
+        else:
+            print "Did not find the Numeric module you asked for"            
+    else:
+        if array_preference == 'numarray':
+            if have_numarray == 1:
+                use_numarray = 1
+            else:
+                print "Did not find the numarray module you asked for"
+
+if use_numeric == 0 and use_numarray == 0:            
+    if have_numeric == 1:
+        use_numeric = 1
+    elif have_numarray == 1:
+        use_numarray = 1
+    else:
+        raise  DistutilsModuleError, "I need either Numeric or nummarray!"
+
+if use_numeric == 1:
+	print "Using Numeric as array Package"
+	use_numarray = 0
+        nummodule = "Numeric"
+
+elif use_numarray == 1:    
+	print "Using nummarray as array Package"
+	use_numeric = 0
+        nummodule = "numarray"
+else:
+	raise  DistutilsModuleError, "I need either Numeric or nummarray!"
+
+# Write the chosen module to a file so it is automatically inculded when pygsl starts up
+file = open(os.path.join("pygsl", "_numobj.py"), "w")
+warnmsg = """
+       WARNING: File Generated during build. DO NOT MODIFY!!!
+"""
+file.write('"""\n')
+file.write(warnmsg)
+file.write('"""\n')
+file.write('\n')
+file.write('from %s import *' % nummodule)
+file.close()
+del file
+
+# Write the chosen module to a include header
+file = open(os.path.join("Include", "pygsl", "arrayobject.h"), "w")
+file.write('/*')
+file.write(warnmsg)
+file.write('*/\n')
+file.write('#include <%s/arrayobject.h>\n' % nummodule)
+file.close()
+del file
+
+# Write the chosen module to a include header
+file = open(os.path.join("pygsl", "_mlab.py"), "w")
+file.write('"""\n')
+file.write(warnmsg)
+file.write('"""\n')
+if use_numeric == 1:
+	file.write('from MLab import *\n')
+elif use_numarray == 1:
+	file.write('from numarray.linear_algebra.mlab import *\n')
+file.close()
+del file

src/statistics/_statmodule.c

+#include "statmodule.h"
+
+static void * __PyGSL_STATISTICS_API[8] = {NULL, NULL, NULL, NULL, NULL, NULL, NULL, NULL};
+
+
+static PyObject*
+PyGSL_statistics_d_A(PyObject *self, PyObject *args,
+		     double (*pointer)(const void *, size_t, size_t),
+		     int array_type, int basis_type_size)
+{
+    PyObject *input = NULL; 
+    PyArrayObject *data; 
+    double result;
+    int stride=1, n; 
+    
+    FUNC_MESS_BEGIN();
+    if(!(PyArg_ParseTuple(args, "O", &input))) 
+	return NULL;
+
+    data = PyGSL_PyArray_PREPARE_gsl_vector_view(input, array_type, 0, -1, 1, NULL);
+    if(data == NULL) 
+	return NULL;
+
+    
+    if(PyGSL_STRIDE_RECALC(data->strides[0], basis_type_size, &stride) != GSL_SUCCESS){
+	 Py_XDECREF(data); 
+	 return NULL;
+    }
+
+    n = data->dimensions[0];
+    result = pointer((void *)data->data, (size_t) stride, (size_t) n);
+    Py_DECREF(data);
+    FUNC_MESS_END();
+    return PyFloat_FromDouble(result);
+}
+
+
+static PyObject*
+PyGSL_statistics_l_A(PyObject *self, PyObject *args, 
+		     size_t (*pointer)(const void *, size_t, size_t),
+		     int array_type, int basis_type_size)
+{
+    PyObject *input = NULL; 
+    PyArrayObject *data; 
+    long result;
+    int stride=1, n; 
+    
+    if(!(PyArg_ParseTuple(args, "O", &input))) 
+	return NULL;
+
+    data = PyGSL_PyArray_PREPARE_gsl_vector_view(input, array_type, 0, -1, 1, NULL);
+    if(data == NULL) 
+	return NULL;
+
+    if(PyGSL_STRIDE_RECALC(data->strides[0], basis_type_size, &stride) != GSL_SUCCESS){
+	 Py_XDECREF(data); 
+	 return NULL;
+    }
+
+
+    n = data->dimensions[0]; 
+    result = pointer((void *)data->data, (size_t) stride, (size_t) n);
+    Py_DECREF(data);
+    return PyInt_FromLong(result);
+}
+
+
+static PyObject*
+PyGSL_statistics_d_Ad(PyObject *self, PyObject *args, 
+		      double (*pointer)(const void *, size_t, size_t, double),
+		      int array_type, int basis_type_size)
+{
+    PyObject *input = NULL; 
+    PyArrayObject *data; 
+    double result, mean;
+    size_t stride=1, n; 
+    
+    if(!(PyArg_ParseTuple(args, "Od", &input, &mean))) 
+	return NULL;
+
+    data = PyGSL_PyArray_PREPARE_gsl_vector_view(input, array_type, 0, -1, 1, NULL);
+    if(data == NULL) 
+	return NULL;
+
+    if(PyGSL_STRIDE_RECALC(data->strides[0], basis_type_size, &stride) != GSL_SUCCESS){
+	 Py_XDECREF(data); 
+	 return NULL;
+    }
+
+
+    n = data->dimensions[0];
+    result = pointer((void *)data->data, (size_t) stride, (size_t) n, mean);
+    Py_DECREF(data);
+    return PyFloat_FromDouble(result);
+}
+
+static PyObject*
+PyGSL_statistics_d_AA(PyObject *self, PyObject *args, 
+		      double (*pointer)(const void *, size_t,const void *, size_t, size_t),
+		      int array_type, int basis_type_size)
+{ 
+    PyObject *input1 = NULL, *input2 = NULL;
+    PyArrayObject *data1=NULL, *data2=NULL;
+    double result;
+    int stride1=1, stride2=1, n1;
+ 
+    if(!(PyArg_ParseTuple(args, "OO", &input1, &input2))) 
+	return NULL;
+
+    data1 = PyGSL_PyArray_PREPARE_gsl_vector_view(input1, array_type, 0, -1, 1, NULL);
+    if(data1 == NULL) 
+	return NULL;
+
+    if(PyGSL_STRIDE_RECALC(data1->strides[0], basis_type_size, &stride1) != GSL_SUCCESS){
+	 goto fail;
+    }
+
+    n1 = data1->dimensions[0];
+
+
+    data2 = PyGSL_PyArray_PREPARE_gsl_vector_view(input1, array_type, 0, n1, 1, NULL);
+    if(data2 == NULL){
+	 goto fail;
+    }
+
+    if(PyGSL_STRIDE_RECALC(data2->strides[0], basis_type_size, &stride2) != GSL_SUCCESS){
+	 goto fail;
+    }
+
+	goto fail; 
+    result = pointer((void *)data1->data, (size_t) stride1, (void *)data2->data, (size_t) stride2,  (size_t) n1);
+    Py_DECREF(data1); 
+    Py_DECREF(data2); 
+    return PyFloat_FromDouble(result); 
+ fail:
+    Py_XDECREF(data1);
+    Py_XDECREF(data2);
+    return NULL;
+}
+
+
+static PyObject*
+PyGSL_statistics_d_AAd(PyObject *self, PyObject *args, 
+		       double (*pointer)(const void *, size_t,const void *, size_t, size_t, double),
+		       int array_type, int basis_type_size)
+{ 
+    PyObject *input1 = NULL, *input2 = NULL;
+    PyArrayObject *data1=NULL, *data2=NULL;
+    double result, mean;
+    size_t stride1=1, stride2=1;
+    size_t n1;
+ 
+    if(!(PyArg_ParseTuple(args, "OOd", &input1, &input2, &mean))) 
+	return NULL;
+
+    data1 = PyGSL_PyArray_PREPARE_gsl_vector_view(input1, array_type, 0, -1, 1, NULL);
+    if(data1 == NULL) 
+	return NULL;
+
+    if(PyGSL_STRIDE_RECALC(data1->strides[0], basis_type_size, &stride1) != GSL_SUCCESS){
+	 goto fail;
+    }
+
+    n1 = data1->dimensions[0];
+
+    data2 = PyGSL_PyArray_PREPARE_gsl_vector_view(input1, array_type, 0, n1, 1, NULL);
+    if(data2 == NULL){
+	Py_XDECREF(data1); 
+	return NULL;
+    }
+
+    if(PyGSL_STRIDE_RECALC(data2->strides[0], basis_type_size, &stride2) != GSL_SUCCESS){
+	 goto fail;
+    }
+
+    result = pointer((void *)data1->data, (size_t) stride1, (void *)data2->data, (size_t) stride2, (size_t) n1, mean);
+    Py_DECREF(data1); 
+    Py_DECREF(data2); 
+    return PyFloat_FromDouble(result); 
+ fail:
+    Py_XDECREF(data1);
+    Py_XDECREF(data2);
+    return NULL;
+}
+
+
+
+static PyObject*
+PyGSL_statistics_d_AAdd(PyObject *self, PyObject *args, 
+			double (*pointer)(const void *, size_t,const void *, size_t, size_t, double, double),
+			int array_type, int basis_type_size)
+{ 
+    PyObject *input1 = NULL, *input2 = NULL;
+    PyArrayObject *data1=NULL, *data2=NULL;
+    double result, mean1, mean2;
+    size_t stride1=1, stride2=1;
+    size_t n1;
+ 
+    if(!(PyArg_ParseTuple(args, "OOdd", &input1, &input2, &mean1, &mean2))) 
+	return NULL;
+
+    data1 = PyGSL_PyArray_PREPARE_gsl_vector_view(input1, array_type, 0, -1, 1, NULL);
+    if(data1 == NULL) 
+	return NULL;
+
+    if(PyGSL_STRIDE_RECALC(data1->strides[0], basis_type_size, &stride1) != GSL_SUCCESS){
+	 goto fail;
+    }
+
+    n1 = data1->dimensions[0];
+
+    data2 = PyGSL_PyArray_PREPARE_gsl_vector_view(input1, array_type, 0, n1, 1, NULL);
+    if(data2 == NULL){
+	Py_XDECREF(data1); 
+	return NULL;
+    }
+
+    if(PyGSL_STRIDE_RECALC(data2->strides[0], basis_type_size, &stride2) != GSL_SUCCESS){
+	 goto fail;
+    }
+
+    result = pointer((void *)data1->data, (size_t) stride1, (void *)data2->data, (size_t) stride2, (size_t) n1, mean1, mean2);
+    Py_DECREF(data1); 
+    Py_DECREF(data2); 
+    return PyFloat_FromDouble(result); 
+ fail:
+    Py_XDECREF(data1);
+    Py_XDECREF(data2);
+    return NULL;
+}
+
+
+
+static PyObject*
+PyGSL_statistics_d_Add(PyObject *self, PyObject *args, 
+		       double (*pointer)(const void *, size_t, size_t, double, double),
+		       int array_type, int basis_type_size)
+{
+    PyObject *input = NULL; 
+    PyArrayObject *data; 
+    double result, mean1, mean2;
+    size_t stride=1, n; 
+    
+    if(!(PyArg_ParseTuple(args, "Odd", &input, &mean1, &mean2))) 
+	return NULL;
+
+    data = PyGSL_PyArray_PREPARE_gsl_vector_view(input, array_type, 0, -1, 1, NULL);
+    if(data == NULL) 
+	return NULL;
+
+    if(PyGSL_STRIDE_RECALC(data->strides[0], basis_type_size, &stride) != GSL_SUCCESS){
+	 Py_XDECREF(data); 
+	 return NULL;
+    }
+
+
+    n = data->dimensions[0];
+    result = pointer((void *)data->data, (size_t) stride, (size_t) n, mean1, mean2);
+    Py_DECREF(data);
+    return PyFloat_FromDouble(result);
+}
+
+
+static PyObject*
+PyGSL_statistics_ll_A(PyObject *self, PyObject *args, 
+		      void (*pointer)(size_t *, size_t *, const void *, size_t, size_t),
+		      int array_type, int basis_type_size)
+{
+    PyObject *input = NULL; 
+    PyArrayObject *data; 
+    size_t result1, result2;
+    size_t stride=1, n; 
+    
+    if(!(PyArg_ParseTuple(args, "O", &input))) 
+	return NULL;
+
+    data = PyGSL_PyArray_PREPARE_gsl_vector_view(input, array_type, 0, -1, 1, NULL);
+    if(data == NULL) 
+	return NULL;
+
+    if(PyGSL_STRIDE_RECALC(data->strides[0], basis_type_size, &stride) != GSL_SUCCESS){
+	 Py_XDECREF(data); 
+	 return NULL;
+    }
+
+
+    n = data->dimensions[0];
+    pointer(&result1, &result2, (void *)data->data, (size_t) stride, (size_t) n);
+    Py_DECREF(data);
+    return Py_BuildValue("ll", (long) result1, (long) result2);
+}
+
+void 
+set_api_pointer(void)
+{
+
+     FUNC_MESS_BEGIN();
+
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_d_A_NUM]    = (void *) PyGSL_statistics_d_A; 
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_l_A_NUM]    = (void *) PyGSL_statistics_l_A;
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_d_Ad_NUM]   = (void *) PyGSL_statistics_d_Ad;
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_d_AA_NUM]   = (void *) PyGSL_statistics_d_AA;
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_d_AAd_NUM]  = (void *) PyGSL_statistics_d_AAd;
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_d_AAdd_NUM] = (void *) PyGSL_statistics_d_AAdd;
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_d_Add_NUM]  = (void *) PyGSL_statistics_d_Add;
+    __PyGSL_STATISTICS_API[PyGSL_STATISTICS_ll_A_NUM]   = (void *) PyGSL_statistics_ll_A;
+     PyGSL_STATISTICS_API = __PyGSL_STATISTICS_API;
+     /* fprintf(stderr, "__PyGSL_RNG_API @ %p\n", (void *) __PyGSL_RNG_API); */
+     FUNC_MESS_END();
+}
+
+/* initialization */
+
+DL_EXPORT(void) init_stat(void)
+{
+     PyObject *api, *dict, *m;
+
+     import_array();
+     init_pygsl();
+     m = Py_InitModule("_stat", NULL);
+
+     if(m == NULL)
+	  goto fail;
+     dict = PyModule_GetDict(m);
+     if(dict == NULL)
+	  goto fail;
+
+     set_api_pointer();
+     api = PyCObject_FromVoidPtr((void *) PyGSL_STATISTICS_API, NULL);
+     assert(api);
+     if (PyDict_SetItemString(dict, "_PYGSL_STATISTICS_API", api) != 0){
+	  PyErr_SetString(PyExc_ImportError, 
+			  "I could not add  _PYGSL_RNG_API!");
+	  goto fail;
+     }
+    return;
+ fail:
+     if(!PyErr_Occurred()){
+	  PyErr_SetString(PyExc_ImportError, "I could not init statistics._stat module!");
+     }
+
+}
+
Tip: Filter by directory path e.g. /media app.js to search for public/media/app.js.
Tip: Use camelCasing e.g. ProjME to search for ProjectModifiedEvent.java.
Tip: Filter by extension type e.g. /repo .js to search for all .js files in the /repo directory.
Tip: Separate your search with spaces e.g. /ssh pom.xml to search for src/ssh/pom.xml.
Tip: Use ↑ and ↓ arrow keys to navigate and return to view the file.
Tip: You can also navigate files with Ctrl+j (next) and Ctrl+k (previous) and view the file with Ctrl+o.
Tip: You can also navigate files with Alt+j (next) and Alt+k (previous) and view the file with Alt+o.