Back to index

python-biopython  1.60
_pwm.c
Go to the documentation of this file.
00001 #include <Python.h>
00002 #include "numpy/arrayobject.h"
00003 
00004 
00005 
00006 static PyObject*
00007 calculate(const char sequence[], int s, PyObject* matrix, npy_intp m)
00008 {
00009     npy_intp n = s - m + 1;
00010     npy_intp i, j;
00011     char c;
00012     double score;
00013     int ok;
00014     PyObject* result;
00015     float* p;
00016     npy_intp shape = (npy_intp)n;
00017     float nan = 0.0;
00018     nan /= nan;
00019     if ((int)shape!=n)
00020     {
00021         PyErr_SetString(PyExc_ValueError, "integer overflow");
00022         return NULL;
00023     }
00024     result = PyArray_SimpleNew(1, &shape, NPY_FLOAT32);
00025     if (!result)
00026     {
00027         PyErr_SetString(PyExc_MemoryError, "failed to create output data");
00028         return NULL;
00029     }
00030     p = PyArray_DATA(result);
00031     for (i = 0; i < n; i++)
00032     {
00033         score = 0.0;
00034         ok = 1;
00035         for (j = 0; j < m; j++)
00036         {
00037             c = sequence[i+j];
00038             switch (c)
00039             {
00040                 case 'A':
00041                 case 'a':
00042                     score += *((double*)PyArray_GETPTR2(matrix, j, 0)); break;
00043                 case 'C':
00044                 case 'c':
00045                     score += *((double*)PyArray_GETPTR2(matrix, j, 1)); break;
00046                 case 'G':
00047                 case 'g':
00048                     score += *((double*)PyArray_GETPTR2(matrix, j, 2)); break;
00049                 case 'T':
00050                 case 't':
00051                     score += *((double*)PyArray_GETPTR2(matrix, j, 3)); break;
00052                 default:
00053                     ok = 0;
00054             }
00055         }
00056         if (ok) *p = (float)score;
00057         else *p = nan;
00058         p++;
00059     }
00060     return result;
00061 }
00062 
00063 static char calculate__doc__[] = 
00064 "    calculate(sequence, pwm) -> array of score values\n"
00065 "\n"
00066 "This function calculates the position-weight matrix scores for all\n"
00067 "positions along the sequence, and returns them as a Numerical Python\n"
00068 "array.\n";
00069 
00070 static PyObject*
00071 py_calculate(PyObject* self, PyObject* args, PyObject* keywords)
00072 {
00073     const char* sequence;
00074     PyObject* matrix = NULL;
00075     static char* kwlist[] = {"sequence", "matrix", NULL};
00076     npy_intp m;
00077     int s;
00078     PyObject* result;
00079     if(!PyArg_ParseTupleAndKeywords(args, keywords, "s#O&", kwlist,
00080                                     &sequence,
00081                                     &s,
00082                                     PyArray_Converter,
00083                                     &matrix)) return NULL;
00084 
00085     if (PyArray_TYPE(matrix) != NPY_DOUBLE)
00086     {
00087         PyErr_SetString(PyExc_ValueError,
00088             "position-weight matrix should contain floating-point values");
00089         result = NULL;
00090     }
00091     else if (PyArray_NDIM(matrix) != 2) /* Checking number of dimensions */
00092     {
00093         result = PyErr_Format(PyExc_ValueError, 
00094             "position-weight matrix has incorrect rank (%d expected 2)",
00095             PyArray_NDIM(matrix));
00096     }
00097     else if(PyArray_DIM(matrix, 1) != 4)
00098     {
00099         result = PyErr_Format(PyExc_ValueError,
00100             "position-weight matrix should have four columns (%" NPY_INTP_FMT
00101             " columns found)", PyArray_DIM(matrix, 1));
00102     }
00103     else
00104     {
00105         m = PyArray_DIM(matrix, 0);
00106         result = calculate(sequence, s, matrix, m);
00107     }
00108     Py_DECREF(matrix);
00109     return result;
00110 }
00111 
00112 
00113 static struct PyMethodDef methods[] = {
00114    {"calculate", (PyCFunction)py_calculate, METH_VARARGS | METH_KEYWORDS, calculate__doc__},
00115    {NULL,          NULL, 0, NULL} /* sentinel */
00116 };
00117 
00118 
00119 #if PY_MAJOR_VERSION >= 3
00120 
00121 static struct PyModuleDef moduledef = {
00122         PyModuleDef_HEAD_INIT,
00123         "_pwm",
00124         "Fast calculations involving position-weight matrices",
00125         -1,
00126         methods,
00127         NULL,
00128         NULL,
00129         NULL,
00130         NULL
00131 };
00132 
00133 PyObject*
00134 PyInit__pwm(void)
00135 
00136 #else
00137 
00138 void init_pwm(void)
00139 #endif
00140 {
00141   PyObject *m;
00142 
00143   import_array();
00144 
00145 #if PY_MAJOR_VERSION >= 3
00146   m = PyModule_Create(&moduledef);
00147   if (m==NULL) return NULL;
00148 #else
00149   m = Py_InitModule4("_pwm",
00150                      methods,
00151                      "Fast calculations involving position-weight matrices",
00152                      NULL,
00153                      PYTHON_API_VERSION);
00154   if (m==NULL) return;
00155 #endif
00156 
00157   if (PyErr_Occurred()) Py_FatalError("can't initialize module _pwm");
00158 #if PY_MAJOR_VERSION >= 3
00159     return m;
00160 #endif
00161 }