diff --git a/src/ClusterReader.c b/src/ClusterReader.c index 3b53021..e0c864c 100644 --- a/src/ClusterReader.c +++ b/src/ClusterReader.c @@ -27,7 +27,6 @@ static int ClusterFileReader_init(ClusterFileReader *self, PyObject *args, static char *kwlist[] = {"fname", "chunk", NULL}; - if (!PyArg_ParseTupleAndKeywords(args, kwds, "O|n", kwlist, &fname_obj, &self->chunk)) { return -1; @@ -42,8 +41,7 @@ static int ClusterFileReader_init(ClusterFileReader *self, PyObject *args, printf("Opening: %s\n chunk: %lu\n", fname, self->chunk); #endif - self->fp = fopen((const char*)fname, "rb"); - + self->fp = fopen((const char *)fname, "rb"); // Keep the return code to not return before releasing buffer int rc = 0; @@ -116,14 +114,12 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, noise_map = (double *)(PyArray_DATA((PyArrayObject *)(noise_array))); - /* for (int i=0; i< ndim_noise; i++) { */ /* printf("Dimension %d size %d pointer \n",i,noise_shape[i], * noise_map); */ /* } */ - if (ndim_noise == 2) { nx = noise_shape[0]; @@ -143,7 +139,7 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, } // Create an uninitialized numpy array - PyObject *clusters = PyArray_SimpleNewFromDescr(ndim, dims, cluster_dt()); + PyObject *clusters = PyArray_SimpleNewFromDescr(ndim, dims, cluster_dt()); // Fill with zeros PyArray_FILLWBYTE((PyArrayObject *)clusters, 0); @@ -158,7 +154,7 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, n_read = read_clusters_with_cut(self->fp, size, buf, &self->n_left, noise_map, nx, ny); else - n_read = read_clusters(self->fp, size, buf, &self->n_left); + n_read = read_clusters(self->fp, size, buf, &self->n_left); if (n_read != size) { // resize the array to match the number of read photons @@ -180,10 +176,11 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, } /* // clusterize method */ -/* static PyObject *ClusterFileReader_clusterize(ClusterFileReader *self, PyObject *args) { */ +/* static PyObject *ClusterFileReader_clusterize(ClusterFileReader *self, + * PyObject *args) { */ /* const int ndim = 1; */ - + /* Py_ssize_t size = 0; */ /* PyObject *data_obj; */ /* if (!PyArg_ParseTuple(args, "nO", &size,&data_obj)) { */ @@ -195,36 +192,37 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, /* // */ -/* // Create two numpy arrays from the passed objects, if possible numpy will */ -/* // use the underlying buffer, otherwise it will create a copy, for example */ +/* // Create two numpy arrays from the passed objects, if possible numpy + * will */ +/* // use the underlying buffer, otherwise it will create a copy, for + * example */ /* // if data type is different or we pass in a list. The */ /* // NPY_ARRAY_C_CONTIGUOUS flag ensures that we have contiguous memory. */ -/* PyObject *data_array = PyArray_FROM_OTF(data_obj, NPY_INT32, NPY_ARRAY_C_CONTIGUOUS); */ +/* PyObject *data_array = PyArray_FROM_OTF(data_obj, NPY_INT32, + * NPY_ARRAY_C_CONTIGUOUS); */ /* int nx=0,ny=0; */ /* int32_t *data=NULL; */ - /* // If parsing of a or b fails we throw an exception in Python */ /* if (data_array ) { */ - + /* int ndim_data = PyArray_NDIM((PyArrayObject *)(data_array)); */ /* npy_intp *data_shape = PyArray_SHAPE((PyArrayObject *)(data_array)); */ - - -/* // For the C++ function call we need pointers (or another C++ type/data */ + +/* // For the C++ function call we need pointers (or another C++ type/data + */ /* // structure) */ - + /* data = (int32_t *)(PyArray_DATA((PyArrayObject *)(data_array))); */ - - /* /\* for (int i=0; i< ndim_noise; i++) { *\/ */ -/* /\* printf("Dimension %d size %d pointer \n",i,noise_shape[i], noise_map); *\/ */ - +/* /\* printf("Dimension %d size %d pointer \n",i,noise_shape[i], + * noise_map); *\/ */ + /* /\* } *\/ */ /* if (ndim_data==2) { */ - + /* nx=data_shape[0]; */ /* ny=data_shape[1]; */ /* if (ny!=9) { */ @@ -240,14 +238,15 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, /* "Wrong data type."); */ /* } */ - + /* } */ /* // Create an uninitialized numpy array */ /* //npy_intp dims[] = {nx}; */ /* // printf("%d %d\n",ndim,nx); */ /* npy_intp dims[] = {nx}; */ -/* PyObject *ca = PyArray_SimpleNewFromDescr(ndim, dims, cluster_analysis_dt()); */ +/* PyObject *ca = PyArray_SimpleNewFromDescr(ndim, dims, + * cluster_analysis_dt()); */ /* // printf("1\n"); */ @@ -260,7 +259,7 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, /* // Call the standalone C code to read clusters from file */ /* // Here goes the looping, removing frame numbers etc. */ - + /* // printf("3\n"); */ /* int n_read=analyze_clusters(nx,data,buf,size); */ /* if (n_read != nx) { */ @@ -280,31 +279,15 @@ static PyObject *ClusterFileReader_read(ClusterFileReader *self, /* } */ /* return ca; */ - + /* } */ - - - - - - - - - - - - - - - - - // List all methods in our ClusterFileReader class static PyMethodDef ClusterFileReader_methods[] = { {"read", (PyCFunction)ClusterFileReader_read, METH_VARARGS, "Read clusters"}, - /* {"clusterize", (PyCFunction)ClusterFileReader_clusterize, METH_VARARGS, */ + /* {"clusterize", (PyCFunction)ClusterFileReader_clusterize, METH_VARARGS, + */ /* "Analyze clusters"}, */ {NULL, NULL, 0, NULL} /* Sentinel */ }; diff --git a/src/RawFileReader.c b/src/RawFileReader.c index beaed49..5023470 100644 --- a/src/RawFileReader.c +++ b/src/RawFileReader.c @@ -48,7 +48,7 @@ static int RawFileReader_init(RawFileReader *self, PyObject *args, printf("fname: %s\n read_header: %d detector type: %d\n", fname, self->read_header, self->detector); #endif - self->fp = fopen((const char*)fname, "rb"); + self->fp = fopen((const char *)fname, "rb"); // Keep the return code to not return before releasing buffer int rc = 0; @@ -110,12 +110,13 @@ static PyObject *RawFileReader_read(RawFileReader *self, PyObject *args) { switch (self->detector) { case DT_MOENCH_03: - n_read = read_raw_m03(self->fp, n_frames, out_buf, (Header*)header_out); + n_read = + read_raw_m03(self->fp, n_frames, out_buf, (Header *)header_out); break; case DT_MOENCH_04_A: case DT_MOENCH_04_AD: - n_read = - read_raw_m04(self->fp, n_frames, out_buf, digital_out, (Header*)header_out); + n_read = read_raw_m04(self->fp, n_frames, out_buf, digital_out, + (Header *)header_out); break; default: break; @@ -136,8 +137,9 @@ static PyObject *RawFileReader_read(RawFileReader *self, PyObject *args) { // resize the array to match the number of clusters read PyArray_Resize((PyArrayObject *)frames, &new_shape, 1, NPY_ANYORDER); - if(digital_frames){ - PyArray_Resize((PyArrayObject *)digital_frames, &new_shape, 1, NPY_ANYORDER); + if (digital_frames) { + PyArray_Resize((PyArrayObject *)digital_frames, &new_shape, 1, + NPY_ANYORDER); } // if we also read header we need to reshape the header @@ -146,32 +148,29 @@ static PyObject *RawFileReader_read(RawFileReader *self, PyObject *args) { PyArray_Resize((PyArrayObject *)header, &new_shape, 1, NPY_ANYORDER); } - - } // Build up a tuple with the return values PyObject *ret = PyTuple_Pack(1, frames); - if (self->detector == DT_MOENCH_04_AD){ + if (self->detector == DT_MOENCH_04_AD) { _PyTuple_Resize(&ret, 2); PyTuple_SET_ITEM(ret, 1, digital_frames); } - if (self->read_header){ + if (self->read_header) { Py_ssize_t old_size = PyTuple_GET_SIZE(ret); - _PyTuple_Resize(&ret, old_size+1); + _PyTuple_Resize(&ret, old_size + 1); PyTuple_SET_ITEM(ret, old_size, header); } // if we only have one item in the tuple lets return it instead of the tuple - if(PyTuple_GET_SIZE(ret) == 1){ + if (PyTuple_GET_SIZE(ret) == 1) { Py_DECREF(ret); return frames; - }else{ + } else { Py_DECREF(frames); return ret; } - } // List all methods in our ClusterFileReader class diff --git a/src/cluster_reader.c b/src/cluster_reader.c index f4aeba7..fc953ff 100644 --- a/src/cluster_reader.c +++ b/src/cluster_reader.c @@ -1,14 +1,14 @@ #include "cluster_reader.h" #include -size_t read_clusters(FILE *fp, size_t n_clusters, Cluster *buf, uint32_t *n_left) { +size_t read_clusters(FILE *fp, size_t n_clusters, Cluster *buf, + uint32_t *n_left) { int32_t iframe = 0; // frame number needs to be 4 bytes! uint32_t nph = *n_left; // number of clusters in frame needs to be 4 bytes! size_t nph_read = 0; uint32_t nn = *n_left; - // if there are photons left from previous frame read them first if (nph) { if (nph > n_clusters) { @@ -19,7 +19,7 @@ size_t read_clusters(FILE *fp, size_t n_clusters, Cluster *buf, uint32_t *n_left nn = nph; } nph_read += fread((void *)(buf + nph_read), sizeof(Cluster), nn, fp); - *n_left = nph - nn; //write back the number of photons left + *n_left = nph - nn; // write back the number of photons left } if (nph_read < n_clusters) { // keep on reading frames and photons until reaching n_clusters @@ -31,7 +31,8 @@ size_t read_clusters(FILE *fp, size_t n_clusters, Cluster *buf, uint32_t *n_left else nn = nph; - nph_read += fread((void *)(buf + nph_read), sizeof(Cluster), nn, fp); + nph_read += + fread((void *)(buf + nph_read), sizeof(Cluster), nn, fp); *n_left = nph - nn; } if (nph_read >= n_clusters) @@ -43,7 +44,8 @@ size_t read_clusters(FILE *fp, size_t n_clusters, Cluster *buf, uint32_t *n_left } size_t read_clusters_with_cut(FILE *fp, size_t n_clusters, Cluster *buf, - uint32_t *n_left, double *noise_map, int nx, int ny) { + uint32_t *n_left, double *noise_map, int nx, + int ny) { int iframe = 0; uint32_t nph = *n_left; @@ -62,7 +64,7 @@ size_t read_clusters_with_cut(FILE *fp, size_t n_clusters, Cluster *buf, if (n_read != 1) { return nph_read; } - //TODO! error handling on read + // TODO! error handling on read good = 1; if (noise_map) { if (ptr->x >= 0 && ptr->x < nx && ptr->y >= 0 && ptr->y < ny) { @@ -92,105 +94,104 @@ size_t read_clusters_with_cut(FILE *fp, size_t n_clusters, Cluster *buf, if (nph_read < n_clusters) { // keep on reading frames and photons until reaching n_clusters while (fread(&iframe, sizeof(iframe), 1, fp)) { - //printf("%d\n",nph_read); + // printf("%d\n",nph_read); - if (fread(&nph, sizeof(nph), 1, fp)) { - //printf("** %d\n",nph); - *n_left = nph; - for (size_t iph=0; iphx >= 0 && ptr->x < nx && ptr->y >= 0 && + ptr->y < ny) { + tot1 = ptr->data[4]; + analyze_cluster(*ptr, &t2max, &tot3, NULL, NULL, + NULL, NULL, NULL); + noise = noise_map[ptr->y * nx + ptr->x]; + if (tot1 > noise && t2max > 2 * noise && + tot3 > 3 * noise) { + ; + } else + good = 0; + } else { + printf("Bad pixel number %d %d\n", ptr->x, ptr->y); + good = 0; + } + } + if (good) { + ptr++; + nph_read++; + (*n_left)--; + } + if (nph_read >= n_clusters) + break; + } + } + + if (nph_read >= n_clusters) + break; } - good=1; - if (noise_map) { - if (ptr->x>=0 && ptr->xy>=0 && ptr->ydata[4]; - analyze_cluster(*ptr, &t2max, &tot3, NULL, NULL, NULL, NULL,NULL); - noise=noise_map[ptr->y*nx+ptr->x]; - if (tot1>noise && t2max>2*noise && tot3>3*noise) { - ; - } else - good=0; - } else{ - printf("Bad pixel number %d %d\n",ptr->x,ptr->y); - good=0; - } - } - if (good) { - ptr++; - nph_read++; - (*n_left)--; - } - if (nph_read >= n_clusters) - break; - - } - } - - if (nph_read >= n_clusters) - break; - } } // printf("%d\n",nph_read); assert(nph_read <= n_clusters); // sanity check in debug mode return nph_read; } - - - - - - -int analyze_clusters(int64_t n_clusters, int32_t *cin, ClusterAnalysis *co, int csize) { +int analyze_clusters(int64_t n_clusters, int32_t *cin, ClusterAnalysis *co, + int csize) { char quad; int32_t tot; double etax, etay; - int nc=0; - //printf("csize is %d\n",csize); + int nc = 0; + // printf("csize is %d\n",csize); int ret; - for (int ic = 0; ic < n_clusters; ic++) { - switch (csize) { - case 2: - ret=analyze_data((cin+9*ic), &tot, NULL, &quad, &etax,&etay, NULL, NULL); - break; - default: - ret=analyze_data((cin+9*ic), NULL, &tot, &quad, NULL, NULL, &etax,&etay); - } - if (ret==0) { - printf("%d %d %d %f %f\n",ic,tot,quad,etax,etay); - - } - nc+=ret; - //printf("%d %d %d %d\n", ic , quad , t2max , tot3); - (co + ic)->c = quad; - (co + ic)->tot = tot; - (co + ic)->etax = etax; - (co + ic)->etay = etay; - //printf("%g %g\n",etax, etay); - /* if (tot<=0) */ - /* printf("%d %d %d %d %d %d\n",ic,(cin+ic)->x, (cin+ic)->y, */ - /* (cout+ic)->c, (cout+ic)->tot2, (cout+ic)->tot3); */ - + for (int ic = 0; ic < n_clusters; ic++) { + switch (csize) { + case 2: + ret = analyze_data((cin + 9 * ic), &tot, NULL, &quad, &etax, &etay, + NULL, NULL); + break; + default: + ret = analyze_data((cin + 9 * ic), NULL, &tot, &quad, NULL, NULL, + &etax, &etay); + } + if (ret == 0) { + printf("%d %d %d %f %f\n", ic, tot, quad, etax, etay); + } + nc += ret; + // printf("%d %d %d %d\n", ic , quad , t2max , tot3); + (co + ic)->c = quad; + (co + ic)->tot = tot; + (co + ic)->etax = etax; + (co + ic)->etay = etay; + // printf("%g %g\n",etax, etay); + /* if (tot<=0) */ + /* printf("%d %d %d %d %d %d\n",ic,(cin+ic)->x, (cin+ic)->y, */ + /* (cout+ic)->c, (cout+ic)->tot2, (cout+ic)->tot3); */ } return nc; } -int analyze_cluster(Cluster cl, int32_t *t2, int32_t *t3, char *quad, double *eta2x, double *eta2y, double *eta3x, double *eta3y) { - - return analyze_data(cl.data, t2, t3, quad, eta2x, eta2y, eta3x, eta3y); +int analyze_cluster(Cluster cl, int32_t *t2, int32_t *t3, char *quad, + double *eta2x, double *eta2y, double *eta3x, + double *eta3y) { + return analyze_data(cl.data, t2, t3, quad, eta2x, eta2y, eta3x, eta3y); } -int analyze_data(int32_t *data, int32_t *t2, int32_t *t3, char *quad, double *eta2x, double *eta2y, double *eta3x, double *eta3y) { +int analyze_data(int32_t *data, int32_t *t2, int32_t *t3, char *quad, + double *eta2x, double *eta2y, double *eta3x, double *eta3y) { - - int ok=1; + int ok = 1; int32_t tot2[4]; - int32_t t2max=0; + int32_t t2max = 0; char c; int32_t val, tot3; @@ -199,23 +200,22 @@ int analyze_data(int32_t *data, int32_t *t2, int32_t *t3, char *quad, double *et tot2[i] = 0; for (int ix = 0; ix < 3; ix++) { - for (int iy = 0; iy < 3; iy++) { - val = data[iy * 3 + ix]; - // printf ("%d ",data[iy * 3 + ix]); - tot3 += val; - if (ix <= 1 && iy <= 1) - tot2[cBottomLeft] += val; - if (ix >= 1 && iy <= 1) - tot2[cBottomRight] += val; - if (ix <= 1 && iy >= 1) - tot2[cTopLeft] += val; - if (ix >= 1 && iy >= 1) - tot2[cTopRight] += val; - } - // printf ("\n"); - + for (int iy = 0; iy < 3; iy++) { + val = data[iy * 3 + ix]; + // printf ("%d ",data[iy * 3 + ix]); + tot3 += val; + if (ix <= 1 && iy <= 1) + tot2[cBottomLeft] += val; + if (ix >= 1 && iy <= 1) + tot2[cBottomRight] += val; + if (ix <= 1 && iy >= 1) + tot2[cTopLeft] += val; + if (ix >= 1 && iy >= 1) + tot2[cTopRight] += val; + } + // printf ("\n"); } - //printf ("\n"); + // printf ("\n"); if (t2 || quad) { @@ -227,57 +227,57 @@ int analyze_data(int32_t *data, int32_t *t2, int32_t *t3, char *quad, double *et c = i; } } - - if (quad) - *quad = c; - if (t2) - *t2 = t2max; + + if (quad) + *quad = c; + if (t2) + *t2 = t2max; } if (t3) - *t3 = tot3; - + *t3 = tot3; + if (eta2x || eta2y) { - if (eta2x ) - *eta2x=0; - if (eta2y ) - *eta2y=0; - switch (c) { - case cBottomLeft: - if (eta2x && (data[3]+data[4])!=0) - *eta2x=(double)(data[4])/(data[3]+data[4]); - if (eta2y && (data[1]+data[4])!=0) - *eta2y=(double)(data[4])/(data[1]+data[4]); - break; - case cBottomRight: - if (eta2x && (data[2]+data[5])!=0) - *eta2x=(double)(data[5])/(data[4]+data[5]); - if (eta2y && (data[1]+data[4])!=0) - *eta2y=(double)(data[4])/(data[1]+data[4]); - break; - case cTopLeft: - if (eta2x && (data[7]+data[4])!=0) - *eta2x=(double)(data[4])/(data[3]+data[4]); - if (eta2y && (data[7]+data[4])!=0) - *eta2y=(double)(data[7])/(data[7]+data[4]); - break; - case cTopRight: - if (eta2x && t2max!=0) - *eta2x=(double)(data[5])/(data[5]+data[4]); - if (eta2y && t2max!=0) - *eta2y=(double)(data[7])/(data[7]+data[4]); - break; - default: - ; - } + if (eta2x) + *eta2x = 0; + if (eta2y) + *eta2y = 0; + switch (c) { + case cBottomLeft: + if (eta2x && (data[3] + data[4]) != 0) + *eta2x = (double)(data[4]) / (data[3] + data[4]); + if (eta2y && (data[1] + data[4]) != 0) + *eta2y = (double)(data[4]) / (data[1] + data[4]); + break; + case cBottomRight: + if (eta2x && (data[2] + data[5]) != 0) + *eta2x = (double)(data[5]) / (data[4] + data[5]); + if (eta2y && (data[1] + data[4]) != 0) + *eta2y = (double)(data[4]) / (data[1] + data[4]); + break; + case cTopLeft: + if (eta2x && (data[7] + data[4]) != 0) + *eta2x = (double)(data[4]) / (data[3] + data[4]); + if (eta2y && (data[7] + data[4]) != 0) + *eta2y = (double)(data[7]) / (data[7] + data[4]); + break; + case cTopRight: + if (eta2x && t2max != 0) + *eta2x = (double)(data[5]) / (data[5] + data[4]); + if (eta2y && t2max != 0) + *eta2y = (double)(data[7]) / (data[7] + data[4]); + break; + default:; + } } - + if (eta3x || eta3y) { - if (eta3x && (data[3]+data[4]+data[5])!=0) - *eta3x=(double)(-data[3]+data[3+2])/(data[3]+data[4]+data[5]); - if (eta3y && (data[1]+data[4]+data[7])!=0) - *eta3y=(double)(-data[1]+data[2*3+1])/(data[1]+data[4]+data[7]); + if (eta3x && (data[3] + data[4] + data[5]) != 0) + *eta3x = (double)(-data[3] + data[3 + 2]) / + (data[3] + data[4] + data[5]); + if (eta3y && (data[1] + data[4] + data[7]) != 0) + *eta3y = (double)(-data[1] + data[2 * 3 + 1]) / + (data[1] + data[4] + data[7]); } - return ok; } diff --git a/src/cluster_reader.h b/src/cluster_reader.h index 4892c59..32668d0 100644 --- a/src/cluster_reader.h +++ b/src/cluster_reader.h @@ -3,17 +3,20 @@ #include #include "data_types.h" -//Pure C implementation to read a cluster file +// Pure C implementation to read a cluster file -size_t read_clusters(FILE* fp, size_t n_clusters, Cluster* buf, uint32_t *n_left); +size_t read_clusters(FILE *fp, size_t n_clusters, Cluster *buf, + uint32_t *n_left); -size_t read_clusters_with_cut(FILE* fp, size_t n_clusters, Cluster* buf, uint32_t *n_left, double *noise_map, int nx, int ny); +size_t read_clusters_with_cut(FILE *fp, size_t n_clusters, Cluster *buf, + uint32_t *n_left, double *noise_map, int nx, + int ny); -int analyze_clusters(int64_t n_clusters, int32_t* cin, ClusterAnalysis *cout, int csize); +int analyze_clusters(int64_t n_clusters, int32_t *cin, ClusterAnalysis *cout, + int csize); +int analyze_data(int32_t *data, int32_t *t2, int32_t *t3, char *quad, + double *eta2x, double *eta2y, double *eta3x, double *eta3y); - - -int analyze_data(int32_t *data, int32_t *t2, int32_t *t3, char *quad, double *eta2x, double *eta2y, double *eta3x, double *eta3y); - -int analyze_cluster(Cluster data, int32_t *t2, int32_t *t3, char *quad, double *eta2x, double *eta2y, double *eta3x, double *eta3y); +int analyze_cluster(Cluster data, int32_t *t2, int32_t *t3, char *quad, + double *eta2x, double *eta2y, double *eta3x, double *eta3y); diff --git a/src/creader_module.c b/src/creader_module.c index 575b30a..3ffbefc 100644 --- a/src/creader_module.c +++ b/src/creader_module.c @@ -30,7 +30,8 @@ /\* (PyArrayObject *)cl_obj, cluster_dt(), NPY_ARRAY_C_CONTIGUOUS); *\/ /\* if (cl_array == NULL) { *\/ /\* PyErr_SetString(PyExc_TypeError, *\/ - /\* "Could not convert first argument to numpy array."); *\/ + /\* "Could not convert first argument to numpy array."); +*\/ /\* return NULL; *\/ /\* } *\/ @@ -56,92 +57,83 @@ */ - // clusterize method -//static PyObject *ClusterFileReader_clusterize(ClusterFileReader *self, PyObject *args) { +// static PyObject *ClusterFileReader_clusterize(ClusterFileReader *self, +// PyObject *args) { static PyObject *clusterize(PyObject *Py_UNUSED(self), PyObject *args) { const int ndim = 1; - + Py_ssize_t size = 0; PyObject *data_obj; if (!PyArg_ParseTuple(args, "nO", &size, &data_obj)) { - PyErr_SetString( - PyExc_TypeError, - "Could not parse args."); - return NULL; - } + PyErr_SetString(PyExc_TypeError, "Could not parse args."); + return NULL; + } - // + // // Create two numpy arrays from the passed objects, if possible numpy will // use the underlying buffer, otherwise it will create a copy, for example // if data type is different or we pass in a list. The // NPY_ARRAY_C_CONTIGUOUS flag ensures that we have contiguous memory. - PyObject *data_array = PyArray_FROM_OTF(data_obj, NPY_INT32, NPY_ARRAY_C_CONTIGUOUS); - int nx=0,ny=0; - int32_t *data=NULL; - + PyObject *data_array = + PyArray_FROM_OTF(data_obj, NPY_INT32, NPY_ARRAY_C_CONTIGUOUS); + int nx = 0, ny = 0; + int32_t *data = NULL; // If parsing of a or b fails we throw an exception in Python - if (data_array ) { - - int ndim_data = PyArray_NDIM((PyArrayObject *)(data_array)); - npy_intp *data_shape = PyArray_SHAPE((PyArrayObject *)(data_array)); - - - // For the C++ function call we need pointers (or another C++ type/data - // structure) - - data = (int32_t *)(PyArray_DATA((PyArrayObject *)(data_array))); + if (data_array) { + int ndim_data = PyArray_NDIM((PyArrayObject *)(data_array)); + npy_intp *data_shape = PyArray_SHAPE((PyArrayObject *)(data_array)); + // For the C++ function call we need pointers (or another C++ type/data + // structure) - /* for (int i=0; i< ndim_noise; i++) { */ - /* printf("Dimension %d size %d pointer \n",i,noise_shape[i], noise_map); */ - - /* } */ + data = (int32_t *)(PyArray_DATA((PyArrayObject *)(data_array))); - if (ndim_data==2) { - - nx=data_shape[0]; - ny=data_shape[1]; - if (ny!=9) { - PyErr_SetString( - PyExc_TypeError, - "Wrong data type."); - // printf("Data found size %d %d %d\n",nx,ny,ndim); - } + /* for (int i=0; i< ndim_noise; i++) { */ + /* printf("Dimension %d size %d pointer \n",i,noise_shape[i], + * noise_map); */ - } else { - PyErr_SetString( - PyExc_TypeError, - "Wrong data type."); + /* } */ - } - + if (ndim_data == 2) { + + nx = data_shape[0]; + ny = data_shape[1]; + if (ny != 9) { + PyErr_SetString(PyExc_TypeError, "Wrong data type."); + // printf("Data found size %d %d %d\n",nx,ny,ndim); + } + + } else { + PyErr_SetString(PyExc_TypeError, "Wrong data type."); + } } // Create an uninitialized numpy array - //npy_intp dims[] = {nx}; + // npy_intp dims[] = {nx}; // printf("%d %d\n",ndim,nx); npy_intp dims[] = {nx}; - PyObject *ca = PyArray_SimpleNewFromDescr(ndim, dims, cluster_analysis_dt()); + PyObject *ca = + PyArray_SimpleNewFromDescr(ndim, dims, cluster_analysis_dt()); - //printf("1\n"); + // printf("1\n"); // Fill with zeros PyArray_FILLWBYTE((PyArrayObject *)ca, 0); - //printf("2\n"); - // Get a pointer to the array memory + // printf("2\n"); + // Get a pointer to the array memory void *buf = PyArray_DATA((PyArrayObject *)ca); // Call the standalone C code to read clusters from file // Here goes the looping, removing frame numbers etc. - - //printf("3\n"); - int nc=analyze_clusters(nx,data,buf,size); - + + // printf("3\n"); + int nc = analyze_clusters(nx, data, buf, size); + // printf("aa %d %d\n",n_read, nx); /* if (nc != nx) { */ /* // resize the array to match the number of read photons */ @@ -159,16 +151,13 @@ static PyObject *clusterize(PyObject *Py_UNUSED(self), PyObject *args) { /* PyArray_Resize((PyArrayObject *)ca, &new_shape, 1, NPY_ANYORDER); */ /* } */ if (nc != nx) { - printf("%d %d\n",nx,nc); - PyErr_SetString(PyExc_TypeError, "Parsed wrong size array!"); + printf("%d %d\n", nx, nc); + PyErr_SetString(PyExc_TypeError, "Parsed wrong size array!"); } Py_DECREF(data_array); return ca; - } - - static PyObject *get_cluster_dt(PyObject *Py_UNUSED(self), PyObject *args) { if (!PyArg_ParseTuple(args, "")) return NULL; diff --git a/src/raw_reader.h b/src/raw_reader.h index e50e1a6..f3f96db 100644 --- a/src/raw_reader.h +++ b/src/raw_reader.h @@ -1,27 +1,16 @@ #pragma once +#include "data_types.h" #include #include #include -#include "data_types.h" -int64_t read_raw_m03( - FILE *fp, - int64_t n_frames, - char* frame_out, - Header* header_out -); +int64_t read_raw_m03(FILE *fp, int64_t n_frames, char *frame_out, + Header *header_out); -int64_t read_raw_m04( - FILE *fp, - int64_t n_frames, - char* frame_out, - char* digital_out, - Header* header_out -); +int64_t read_raw_m04(FILE *fp, int64_t n_frames, char *frame_out, + char *digital_out, Header *header_out); void decode_moench03(const uint16_t *buf, uint16_t *out_buf); -void decode_moench04(const uint16_t *analog_data, - const uint64_t *digital_data, - uint16_t *analog_frame, - uint8_t *digital_frame); \ No newline at end of file +void decode_moench04(const uint16_t *analog_data, const uint64_t *digital_data, + uint16_t *analog_frame, uint8_t *digital_frame); \ No newline at end of file