From 627719cae59bbcc47149e9390da9de196eb0f0c2 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Mon, 28 Jul 2014 09:08:10 +0200 Subject: [PATCH 01/60] Vectorize all d2dtf() imputs. For later automation, it is easier to have all inputs vectorized. --- cython_numpy/erfa.pyx | 24 ++++++++++++++---------- 1 file changed, 14 insertions(+), 10 deletions(-) diff --git a/cython_numpy/erfa.pyx b/cython_numpy/erfa.pyx index 6504f38..01abb66 100644 --- a/cython_numpy/erfa.pyx +++ b/cython_numpy/erfa.pyx @@ -74,14 +74,16 @@ def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, php def d2dtf(scale, ndp, d1, d2): - shape = np.broadcast(d1, d2).shape + shape = np.broadcast(scale, ndp, d1, d2).shape iy = np.empty(shape, dtype=np.int) im = np.empty(shape, dtype=np.int) id = np.empty(shape, dtype=np.int) ihmsf = np.empty(shape, dtype=[('h','i'),('m','i'),('s','i'),('f','i')]) - cdef np.broadcast it = np.broadcast(d1, d2, iy, im, id, ihmsf) + cdef np.broadcast it = np.broadcast(scale, ndp, d1, d2, iy, im, id, ihmsf) + cdef char *_scale + cdef int _ndp cdef int _iy cdef int _im cdef int _id @@ -89,16 +91,18 @@ def d2dtf(scale, ndp, d1, d2): while np.PyArray_MultiIter_NOTDONE(it): - _d1 = (np.PyArray_MultiIter_DATA(it, 0))[0] - _d2 = (np.PyArray_MultiIter_DATA(it, 1))[0] + _scale = ( np.PyArray_MultiIter_DATA(it, 0)) + _ndp = ( np.PyArray_MultiIter_DATA(it, 1))[0] + _d1 = (np.PyArray_MultiIter_DATA(it, 2))[0] + _d2 = (np.PyArray_MultiIter_DATA(it, 3))[0] - ret = eraD2dtf(scale, ndp, _d1, _d2, &_iy, &_im, &_id, _ihmsf) + ret = eraD2dtf(_scale, _ndp, _d1, _d2, &_iy, &_im, &_id, _ihmsf) - (np.PyArray_MultiIter_DATA(it, 2))[0] = _iy - (np.PyArray_MultiIter_DATA(it, 3))[0] = _im - (np.PyArray_MultiIter_DATA(it, 4))[0] = _id - (np.PyArray_MultiIter_DATA(it, 5))[0:4] = _ihmsf + (np.PyArray_MultiIter_DATA(it, 4))[0] = _iy + (np.PyArray_MultiIter_DATA(it, 5))[0] = _im + (np.PyArray_MultiIter_DATA(it, 6))[0] = _id + (np.PyArray_MultiIter_DATA(it, 7))[0:4] = _ihmsf np.PyArray_MultiIter_NEXT(it) - return iy, im, id, ihmsf + return iy, im, id, ihmsf \ No newline at end of file From b815ba73db86f13006e9176c1a3c3b3590476074 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Mon, 28 Jul 2014 09:09:24 +0200 Subject: [PATCH 02/60] Added aper() Just because one of its parameters is an input and an output at the same time. And it also uses the complex eraASTROM structure. --- cython_numpy/erfa.pyx | 71 +++++++++++++++++++++++++++++++++++++++++-- 1 file changed, 69 insertions(+), 2 deletions(-) diff --git a/cython_numpy/erfa.pyx b/cython_numpy/erfa.pyx index 01abb66..62a6ca5 100644 --- a/cython_numpy/erfa.pyx +++ b/cython_numpy/erfa.pyx @@ -3,7 +3,8 @@ cimport numpy as np np.import_array() -__all__ = ['atco13', 'd2dtf'] +__all__ = ['atco13', 'd2dtf', 'aper'] + cdef extern from "erfa.h": int eraAtco13(double rc, double dc, @@ -15,10 +16,50 @@ cdef extern from "erfa.h": double *dob, double *rob, double *eo) int eraD2dtf(const char *scale, int ndp, double d1, double d2, int *iy, int *im, int *id, int ihmsf[4]) + void eraAper(double theta, eraASTROM *astrom) + +cdef struct eraASTROM: + double pmt + double eb[3] + double eh[3] + double em + double v[3] + double bm1 + double bpn[3][3] + double along + double phi + double xpl + double ypl + double sphi + double cphi + double diurab + double eral + double refa + double refb + +dt_eraASTROM = np.dtype([('pmt','d'), + ('eb','d',(3,)), + ('eh','d',(3,)), + ('em','d'), + ('v','d',(3,)), + ('bm1 ','d'), + ('bpn','d',(3,3)), + ('along','d'), + ('phi','d'), + ('xpl','d'), + ('ypl','d'), + ('sphi','d'), + ('cphi','d'), + ('diurab','d'), + ('eral','d'), + ('refa','d'), + ('refb','d')], align=True) + #Note: the pattern used here follows https://github.com/cython/cython/wiki/tutorials-numpy#dimensionally-simple-functions + def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl): shape = np.broadcast(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl).shape @@ -72,6 +113,8 @@ def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, php return aob, zob, hob, dob, rob, eo + + def d2dtf(scale, ndp, d1, d2): shape = np.broadcast(scale, ndp, d1, d2).shape @@ -105,4 +148,28 @@ def d2dtf(scale, ndp, d1, d2): np.PyArray_MultiIter_NEXT(it) - return iy, im, id, ihmsf \ No newline at end of file + return iy, im, id, ihmsf + + + +def aper(theta, astrom): + + shape = np.broadcast(theta, astrom).shape + astrom_out = np.empty(shape, dtype=dt_eraASTROM) + np.copyto(astrom_out, astrom) + + cdef np.broadcast it = np.broadcast(theta, astrom_out) + + cdef double _theta + cdef eraASTROM *_astrom + + while np.PyArray_MultiIter_NOTDONE(it): + + _theta = ( np.PyArray_MultiIter_DATA(it, 0))[0] + _astrom = (np.PyArray_MultiIter_DATA(it, 1)) + + eraAper(_theta, _astrom) + + np.PyArray_MultiIter_NEXT(it) + + return astrom_out \ No newline at end of file From f1b7ce9e41dd01388c605891c540eabcc6a977e7 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Mon, 28 Jul 2014 09:09:45 +0200 Subject: [PATCH 03/60] Added a test for aper() --- tests/test_erfa.py | 6 ++++++ 1 file changed, 6 insertions(+) diff --git a/tests/test_erfa.py b/tests/test_erfa.py index 80653bc..9bdd66d 100644 --- a/tests/test_erfa.py +++ b/tests/test_erfa.py @@ -22,3 +22,9 @@ print(iy.shape, ihmsf.shape, ihmsf.dtype) print(ihmsf) +astrom = np.zeros([2],dtype=erfa.dt_eraASTROM) +theta = np.arange(0,10.0) +print(theta.shape) +print(astrom.shape) +astrom = erfa.aper(theta[:,None], astrom[None,:]) +print(astrom) From 32097053d5906264d1aa28d0aa46bb127fb53629 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Mon, 28 Jul 2014 09:26:09 +0200 Subject: [PATCH 04/60] Postfixed all output arguments with '_out' --- cython_numpy/erfa.pyx | 28 ++++++++++++++-------------- 1 file changed, 14 insertions(+), 14 deletions(-) diff --git a/cython_numpy/erfa.pyx b/cython_numpy/erfa.pyx index 62a6ca5..b0d8030 100644 --- a/cython_numpy/erfa.pyx +++ b/cython_numpy/erfa.pyx @@ -63,14 +63,14 @@ dt_eraASTROM = np.dtype([('pmt','d'), def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl): shape = np.broadcast(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl).shape - aob = np.empty(shape, dtype=np.double) - zob = np.empty(shape, dtype=np.double) - hob = np.empty(shape, dtype=np.double) - dob = np.empty(shape, dtype=np.double) - rob = np.empty(shape, dtype=np.double) - eo = np.empty(shape, dtype=np.double) + aob_out = np.empty(shape, dtype=np.double) + zob_out = np.empty(shape, dtype=np.double) + hob_out = np.empty(shape, dtype=np.double) + dob_out = np.empty(shape, dtype=np.double) + rob_out = np.empty(shape, dtype=np.double) + eo_out = np.empty(shape, dtype=np.double) - cdef np.broadcast it = np.broadcast(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl, aob, zob, hob, dob, rob, eo) + cdef np.broadcast it = np.broadcast(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out) cdef double _aob cdef double _zob @@ -111,19 +111,19 @@ def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, php np.PyArray_MultiIter_NEXT(it) - return aob, zob, hob, dob, rob, eo + return aob_out, zob_out, hob_out, dob_out, rob_out, eo_out def d2dtf(scale, ndp, d1, d2): shape = np.broadcast(scale, ndp, d1, d2).shape - iy = np.empty(shape, dtype=np.int) - im = np.empty(shape, dtype=np.int) - id = np.empty(shape, dtype=np.int) - ihmsf = np.empty(shape, dtype=[('h','i'),('m','i'),('s','i'),('f','i')]) + iy_out = np.empty(shape, dtype=np.int) + im_out = np.empty(shape, dtype=np.int) + id_out = np.empty(shape, dtype=np.int) + ihmsf_out = np.empty(shape, dtype=[('h','i'),('m','i'),('s','i'),('f','i')]) - cdef np.broadcast it = np.broadcast(scale, ndp, d1, d2, iy, im, id, ihmsf) + cdef np.broadcast it = np.broadcast(scale, ndp, d1, d2, iy_out, im_out, id_out, ihmsf_out) cdef char *_scale cdef int _ndp @@ -148,7 +148,7 @@ def d2dtf(scale, ndp, d1, d2): np.PyArray_MultiIter_NEXT(it) - return iy, im, id, ihmsf + return iy_out, im_out, id_out, ihmsf_out From 64e30685b9147943224ad227d65fae2609896ed9 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Mon, 28 Jul 2014 10:05:36 +0200 Subject: [PATCH 05/60] cdef and assign all C arguments before use --- cython_numpy/erfa.pyx | 66 +++++++++++++++++++++++++++---------------- 1 file changed, 42 insertions(+), 24 deletions(-) diff --git a/cython_numpy/erfa.pyx b/cython_numpy/erfa.pyx index b0d8030..eb2ef87 100644 --- a/cython_numpy/erfa.pyx +++ b/cython_numpy/erfa.pyx @@ -72,12 +72,30 @@ def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, php cdef np.broadcast it = np.broadcast(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out) - cdef double _aob - cdef double _zob - cdef double _hob - cdef double _dob - cdef double _rob - cdef double _eo + cdef double _rc + cdef double _dc + cdef double _pr + cdef double _pd + cdef double _px + cdef double _rv + cdef double _utc1 + cdef double _utc2 + cdef double _dut1 + cdef double _elong + cdef double _phi + cdef double _hm + cdef double _xp + cdef double _yp + cdef double _phpa + cdef double _tk + cdef double _rh + cdef double _wl + cdef double *_aob + cdef double *_zob + cdef double *_hob + cdef double *_dob + cdef double *_rob + cdef double *_eo while np.PyArray_MultiIter_NOTDONE(it): @@ -99,15 +117,14 @@ def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, php _tk = (np.PyArray_MultiIter_DATA(it, 15))[0] _rh = (np.PyArray_MultiIter_DATA(it, 16))[0] _wl = (np.PyArray_MultiIter_DATA(it, 17))[0] + _aob = (np.PyArray_MultiIter_DATA(it, 18)) + _zob = (np.PyArray_MultiIter_DATA(it, 19)) + _hob = (np.PyArray_MultiIter_DATA(it, 20)) + _dob = (np.PyArray_MultiIter_DATA(it, 21)) + _rob = (np.PyArray_MultiIter_DATA(it, 22)) + _eo = (np.PyArray_MultiIter_DATA(it, 23)) - ret = eraAtco13(_rc, _dc, _pr, _pd, _px, _rv, _utc1, _utc2, _dut1, _elong, _phi, _hm, _xp, _yp, _phpa, _tk, _rh, _wl, &_aob, &_zob, &_hob, &_dob, &_rob, &_eo) - - (np.PyArray_MultiIter_DATA(it, 18))[0] = _aob - (np.PyArray_MultiIter_DATA(it, 19))[0] = _zob - (np.PyArray_MultiIter_DATA(it, 20))[0] = _hob - (np.PyArray_MultiIter_DATA(it, 21))[0] = _dob - (np.PyArray_MultiIter_DATA(it, 22))[0] = _rob - (np.PyArray_MultiIter_DATA(it, 23))[0] = _eo + ret = eraAtco13(_rc, _dc, _pr, _pd, _px, _rv, _utc1, _utc2, _dut1, _elong, _phi, _hm, _xp, _yp, _phpa, _tk, _rh, _wl, _aob, _zob, _hob, _dob, _rob, _eo) np.PyArray_MultiIter_NEXT(it) @@ -127,10 +144,12 @@ def d2dtf(scale, ndp, d1, d2): cdef char *_scale cdef int _ndp - cdef int _iy - cdef int _im - cdef int _id - cdef int _ihmsf[4] + cdef double _d1 + cdef double _d2 + cdef int *_iy + cdef int *_im + cdef int *_id + cdef int *_ihmsf while np.PyArray_MultiIter_NOTDONE(it): @@ -138,13 +157,12 @@ def d2dtf(scale, ndp, d1, d2): _ndp = ( np.PyArray_MultiIter_DATA(it, 1))[0] _d1 = (np.PyArray_MultiIter_DATA(it, 2))[0] _d2 = (np.PyArray_MultiIter_DATA(it, 3))[0] + _iy = ( np.PyArray_MultiIter_DATA(it, 4)) + _im = ( np.PyArray_MultiIter_DATA(it, 5)) + _id = ( np.PyArray_MultiIter_DATA(it, 6)) + _ihmsf = ( np.PyArray_MultiIter_DATA(it, 7)) - ret = eraD2dtf(_scale, _ndp, _d1, _d2, &_iy, &_im, &_id, _ihmsf) - - (np.PyArray_MultiIter_DATA(it, 4))[0] = _iy - (np.PyArray_MultiIter_DATA(it, 5))[0] = _im - (np.PyArray_MultiIter_DATA(it, 6))[0] = _id - (np.PyArray_MultiIter_DATA(it, 7))[0:4] = _ihmsf + ret = eraD2dtf(_scale, _ndp, _d1, _d2, _iy, _im, _id, _ihmsf) np.PyArray_MultiIter_NEXT(it) From 0cdfa62a819ea5dc5b4884189f555299aa059563 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Mon, 28 Jul 2014 10:06:11 +0200 Subject: [PATCH 06/60] Use copy constructor array() and broadcast_arrays()... MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit …instead of empty() and copyto() --- cython_numpy/erfa.pyx | 3 +-- 1 file changed, 1 insertion(+), 2 deletions(-) diff --git a/cython_numpy/erfa.pyx b/cython_numpy/erfa.pyx index eb2ef87..f3943e8 100644 --- a/cython_numpy/erfa.pyx +++ b/cython_numpy/erfa.pyx @@ -173,8 +173,7 @@ def d2dtf(scale, ndp, d1, d2): def aper(theta, astrom): shape = np.broadcast(theta, astrom).shape - astrom_out = np.empty(shape, dtype=dt_eraASTROM) - np.copyto(astrom_out, astrom) + astrom_out = np.array(np.broadcast_arrays(theta, astrom)[1], dtype=dt_eraASTROM) cdef np.broadcast it = np.broadcast(theta, astrom_out) From 9c9d296895373e950afafb363feddc1b223be496 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 29 Jul 2014 11:29:31 +0200 Subject: [PATCH 07/60] Added *.py[cod] to ignore list --- .gitignore | 3 ++- 1 file changed, 2 insertions(+), 1 deletion(-) diff --git a/.gitignore b/.gitignore index 1bb95c0..6a9fce5 100644 --- a/.gitignore +++ b/.gitignore @@ -1,5 +1,6 @@ /build .project - .pydevproject + +*.py[cod] From 3dbc8f8e5bcec274ecb75749682ec9b038289b4e Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 29 Jul 2014 11:30:49 +0200 Subject: [PATCH 08/60] Initial commit of code analyzer --- cython_numpy_auto/code_analyzer.py | 152 +++++++++++++++++++++++++++++ 1 file changed, 152 insertions(+) create mode 100644 cython_numpy_auto/code_analyzer.py diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py new file mode 100644 index 0000000..3489372 --- /dev/null +++ b/cython_numpy_auto/code_analyzer.py @@ -0,0 +1,152 @@ +import os.path +import re + + +__all__ = ['Function'] + + +ctype_to_dtype = {'double' : "np.double", + 'double *' : "np.double", + 'int' : "np.int", + 'int *' : "np.int", + 'int[4]' : "np.dtype([('', 'i')]*4)", + 'eraASTROM *': "dt_eraASTROM", + } + + +class FunctionDoc: + + def __init__(self, doc): + self.doc = doc.replace("**"," ").replace("/*"," ").replace("*/"," ") + self.__input = None + self.__output = None + + @property + def input(self): + if self.__input is None: + self.__input = [] + __input = re.search("Given:\n(.+?) \n", self.doc, re.DOTALL).group(1) + for i in __input.split("\n"): + arg_doc = ArgumentDoc(i) + if arg_doc.name is not None: + self.__input.append(arg_doc) + return self.__input + + @property + def output(self): + if self.__output is None: + self.__output = [] + __output = re.search("Returned:\n(.+?) \n", self.doc, re.DOTALL).group(1) + for i in __output.split("\n"): + arg_doc = ArgumentDoc(i) + if arg_doc.name is not None: + self.__output.append(arg_doc) + return self.__output + + def __repr__(self): + return self.doc + +class ArgumentDoc: + + def __init__(self, doc): + match = re.search("^ ([^ ]+)[ ]+([^ ]+)[ ]+(.+)", doc) + if match is not None: + self.name = match.group(1) + self.type = match.group(2) + self.doc = match.group(3) + else: + self.name = None + self.type = None + self.doc = None + + def __repr__(self): + return " {0:15} {1:15} {2}".format(self.name, self.type, self.doc) + +class Argument: + + def __init__(self, definition, doc): + self.__doc = doc + self.__inout_state = None + self.definition = definition.strip() + if self.definition[-1] == "]": + self.ctype, self.name = self.definition.split(" ",1) + self.name, arr = self.name.split("[") + self.ctype += ("["+arr) + elif "*" in self.definition: + self.ctype, self.name = self.definition.split("*", 1) + self.ctype += "*" + else: + self.ctype, self.name = self.definition.split(" ", 1) + + @property + def inout_state(self): + if self.__inout_state is None: + self.__inout_state = '' + for i in self.__doc.input: + if self.name in i.name.split(','): + self.__inout_state = 'in' + for o in self.__doc.output: + if self.name in o.name.split(','): + if self.__inout_state == 'in': + self.__inout_state = 'inout' + else: + self.__inout_state = 'out' + return self.__inout_state + + @property + def is_in(self): + return self.inout_state == 'in' + + @property + def is_out(self): + return self.inout_state == 'out' + + @property + def is_inout(self): + return self.inout_state == 'inout' + + @property + def ctype_ptr(self): + if self.ctype[-1] == ']': + return self.ctype.split('[')[0]+" *" + elif self.ctype[:6] == 'const ': + return self.ctype[6:] + else: + return self.ctype + + @property + def dtype(self): + return ctype_to_dtype[self.ctype] + + def __repr__(self): + return "Argument('{0}', name='{1}', ctype='{2}', inout_state='{3}')".format(self.definition, self.name, self.ctype, self.inout_state) + +class Function: + + def __init__(self, name, source_path): + self.name = name + self.pyname = name.split('era')[-1].lower() + self.filename = name.split("era")[-1].lower()+".c" + self.filepath = os.path.join(os.path.normpath(source_path), self.filename) + pattern = "\n([^\n]+{0}\([^)]+\)).+?(/\*.+?\*/)".format(name) + p = re.compile(pattern, flags=re.DOTALL|re.MULTILINE) + with open(self.filepath) as f: + search = p.search(f.read()) + self.cfunc = search.group(1) + self.__doc = FunctionDoc(search.group(2)) + self.args = [] + for arg in re.search("\(([^)]+)\)", self.cfunc, flags=re.MULTILINE|re.DOTALL).group(1).split(','): + self.args.append(Argument(arg, self.__doc)) + + def args_by_inout(self, inout_filter, prop=None, join=None): + result = [] + for arg in self.args: + if arg.inout_state in inout_filter.split('|'): + if prop is None: + result.append(arg) + else: + result.append(getattr(arg, prop)) + if join is not None: + return join.join(result) + else: + return result From 7523521c2643ab655b59027b1a8eb5c68bd305b6 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 29 Jul 2014 11:31:22 +0200 Subject: [PATCH 09/60] Initial commit of cython_generator --- cython_numpy_auto/cython_generator.py | 62 +++++++++++++++++++++++++++ 1 file changed, 62 insertions(+) create mode 100644 cython_numpy_auto/cython_generator.py diff --git a/cython_numpy_auto/cython_generator.py b/cython_numpy_auto/cython_generator.py new file mode 100644 index 0000000..a0d84a8 --- /dev/null +++ b/cython_numpy_auto/cython_generator.py @@ -0,0 +1,62 @@ +from code_analyzer import * + + +def generate(era_func_names, source_path): + + era_funcs = [Function(name, source_path) for name in era_func_names] + + wrapper = [] + wrapper.append("import numpy as np") + wrapper.append("cimport numpy as np") + wrapper.append("") + wrapper.append("np.import_array()") + wrapper.append("") + wrapper.append("__all__ = [{0}]".format(", ".join(["'{0}'".format(func.pyname) for func in era_funcs]+["'dt_eraASTROM'"]))) + wrapper.append("") + wrapper.append("cdef extern from 'erfa.h':") + wrapper.append(" struct eraASTROM:") + wrapper.append(" pass") + for func in era_funcs: + wrapper.append(" {0}".format(func.cfunc)) + wrapper.append("") + wrapper.append("dt_eraASTROM = np.dtype([('pmt','d'),('eb','d',(3,)),('eh','d',(3,)),('em','d'),('v','d',(3,)),('bm1 ','d'),('bpn','d',(3,3)),('along','d'),('phi','d'),('xpl','d'),('ypl','d'),('sphi','d'),('cphi','d'),('diurab','d'),('eral','d'),('refa','d'),('refb','d')], align=True)") + wrapper.append("") + for func in era_funcs: + wrapper.append("#=== {0} ===".format(func.name)) + wrapper.append("") + wrapper.append("def {0}({1}):".format(func.pyname, func.args_by_inout('in|inout','name',', '))) + wrapper.append(" ") + wrapper.append(" shape = np.broadcast({0}).shape".format(func.args_by_inout('in|inout','name',', '))) + for arg in func.args_by_inout('out'): + wrapper.append(" {0}_out = np.empty(shape, dtype={1})".format(arg.name, arg.dtype)) + for arg in func.args_by_inout('inout'): + wrapper.append(" {0}_out = np.array(np.broadcast_arrays({1})[{2}], dtype={3})".format(arg.name, func.args_by_inout('in|inout','name',', '), func.args.index(arg), arg.dtype)) + wrapper.append(" ") + wrapper.append(" cdef np.broadcast it = np.broadcast({0})".format(', '.join([arg.name if arg.is_in else (arg.name+"_out") for arg in func.args_by_inout('in|inout|out')]))) + wrapper.append(" ") + for arg in func.args_by_inout('in|out|inout'): + wrapper.append(" cdef {0} _{1}".format(arg.ctype_ptr, arg.name)) + wrapper.append(" ") + wrapper.append(" while np.PyArray_MultiIter_NOTDONE(it):") + wrapper.append(" ") + for arg in func.args_by_inout('in|out|inout'): + if arg.ctype_ptr[-1] == '*': + wrapper.append(" _{0} = (<{1}>np.PyArray_MultiIter_DATA(it, {2}))".format(arg.name, arg.ctype_ptr, func.args.index(arg))) + else: + wrapper.append(" _{0} = (<{1}*>np.PyArray_MultiIter_DATA(it, {2}))[0]".format(arg.name, arg.ctype_ptr, func.args.index(arg))) + wrapper.append(" ") + wrapper.append(" {0}({1})".format(func.name, ", ".join(["_"+name for name in func.args_by_inout('in|out|inout','name')]))) + wrapper.append(" ") + wrapper.append(" np.PyArray_MultiIter_NEXT(it)") + wrapper.append(" ") + wrapper.append(" return {0}".format(', '.join([arg+"_out" for arg in func.args_by_inout('out|inout','name')]))) + wrapper.append("") + + with open("./erfa.pyx","w") as f: + f.write("\n".join(wrapper)) + + +if __name__ == '__main__': + era_func_names = ["eraAtco13", "eraD2dtf", "eraAper"] + source_path = "../../erfa/src" + generate(era_func_names, source_path) From 2a6093340b1c693d860741da77fc9a48626aeac9 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 29 Jul 2014 11:32:47 +0200 Subject: [PATCH 10/60] Ignore erfa.pyx and erfa.c in cython_numpy_auto --- .gitignore | 4 ++++ 1 file changed, 4 insertions(+) diff --git a/.gitignore b/.gitignore index 6a9fce5..157f002 100644 --- a/.gitignore +++ b/.gitignore @@ -4,3 +4,7 @@ .pydevproject *.py[cod] + +cython_numpy_auto/erfa.c + +cython_numpy_auto/erfa.pyx From 6e653c169c8d485b08ad7d8a81608e4479d9ce4d Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 29 Jul 2014 11:33:08 +0200 Subject: [PATCH 11/60] Added cython_numpy_auto extension to setup.py --- setup.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/setup.py b/setup.py index acfd3b0..4b7a75a 100644 --- a/setup.py +++ b/setup.py @@ -10,8 +10,14 @@ include_dirs = [np.get_include(), '/opt/local/include'], library_dirs = ['/opt/local/lib']) +cython_numpy_auto_ext = Extension("erfa.cython_numpy_auto", + ["cython_numpy_auto/erfa.pyx"], + libraries=['erfa'], + include_dirs = [np.get_include(), '/opt/local/include'], + library_dirs = ['/opt/local/lib']) + setup(name = "erfa", cmdclass = { "build_ext": build_ext }, packages = ["erfa"], package_dir = {"erfa":"."}, - ext_modules = [cython_numpy_ext]) + ext_modules = [cython_numpy_ext, cython_numpy_auto_ext]) From bcbe7091ecf49cead69b36f25558c47f61e9e89f Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 29 Jul 2014 11:39:17 +0200 Subject: [PATCH 12/60] Run test against both cython_numpy and cython_numpy_auto --- tests/test_erfa.py | 63 +++++++++++++++++++++++++--------------------- 1 file changed, 35 insertions(+), 28 deletions(-) diff --git a/tests/test_erfa.py b/tests/test_erfa.py index 9bdd66d..7126c05 100644 --- a/tests/test_erfa.py +++ b/tests/test_erfa.py @@ -1,30 +1,37 @@ +import importlib import numpy as np -import erfa.cython_numpy as erfa -jd = np.linspace(2456855.5, 2456855.5+1.0/24.0/60.0, 60*2+1) -ra = np.linspace(0.0,np.pi*2.0,5) -dec = np.linspace(-np.pi/2.0,+np.pi/2.0,4) - -aob, zob, hob, dob, rob, eo = erfa.atco13(0.0,0.0,0.0,0.0,0.0,0.0,jd,0.0,0.0,0.0,np.pi/4.0,0.0,0.0,0.0,1014.0,0.0,0.0,0.5) -print(aob.shape) - -aob, zob, hob, dob, rob, eo = erfa.atco13(0.0,0.0,0.0,0.0,0.0,0.0,jd[0],0.0,0.0,0.0,np.pi/4.0,0.0,0.0,0.0,1014.0,0.0,0.0,0.5) -print(aob.shape) - -aob, zob, hob, dob, rob, eo = erfa.atco13(ra[:,None,None],dec[None,:,None],0.0,0.0,0.0,0.0,jd[None,None,:],0.0,0.0,0.0,np.pi/4.0,0.0,0.0,0.0,1014.0,0.0,0.0,0.5) -print(aob.shape) - -iy, im, id, ihmsf = erfa.d2dtf("UTC", 3, jd, 0.0) -print(iy.shape, ihmsf.shape, ihmsf.dtype) -print(ihmsf) - -iy, im, id, ihmsf = erfa.d2dtf("UTC", 3, jd[0], 0.0) -print(iy.shape, ihmsf.shape, ihmsf.dtype) -print(ihmsf) - -astrom = np.zeros([2],dtype=erfa.dt_eraASTROM) -theta = np.arange(0,10.0) -print(theta.shape) -print(astrom.shape) -astrom = erfa.aper(theta[:,None], astrom[None,:]) -print(astrom) +for erfa_wrapper_name in ['cython_numpy', 'cython_numpy_auto']: + + erfa_module_name = '.'.join(['erfa', erfa_wrapper_name]) + print("Testing with {0}...".format(erfa_module_name)) + + erfa = importlib.import_module(erfa_module_name) + + jd = np.linspace(2456855.5, 2456855.5+1.0/24.0/60.0, 60*2+1) + ra = np.linspace(0.0,np.pi*2.0,5) + dec = np.linspace(-np.pi/2.0,+np.pi/2.0,4) + + aob, zob, hob, dob, rob, eo = erfa.atco13(0.0,0.0,0.0,0.0,0.0,0.0,jd,0.0,0.0,0.0,np.pi/4.0,0.0,0.0,0.0,1014.0,0.0,0.0,0.5) + print(aob.shape) + + aob, zob, hob, dob, rob, eo = erfa.atco13(0.0,0.0,0.0,0.0,0.0,0.0,jd[0],0.0,0.0,0.0,np.pi/4.0,0.0,0.0,0.0,1014.0,0.0,0.0,0.5) + print(aob.shape) + + aob, zob, hob, dob, rob, eo = erfa.atco13(ra[:,None,None],dec[None,:,None],0.0,0.0,0.0,0.0,jd[None,None,:],0.0,0.0,0.0,np.pi/4.0,0.0,0.0,0.0,1014.0,0.0,0.0,0.5) + print(aob.shape) + + iy, im, id, ihmsf = erfa.d2dtf("UTC", 3, jd, 0.0) + print(iy.shape, ihmsf.shape, ihmsf.dtype) + print(ihmsf) + + iy, im, id, ihmsf = erfa.d2dtf("UTC", 3, jd[0], 0.0) + print(iy.shape, ihmsf.shape, ihmsf.dtype) + print(ihmsf) + + astrom = np.zeros([2],dtype=erfa.dt_eraASTROM) + theta = np.arange(0,10.0) + print(theta.shape) + print(astrom.shape) + astrom = erfa.aper(theta[:,None], astrom[None,:]) + print(astrom) From d1dcc3b513e3767dd973c582c96d044f17ec123e Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 30 Jul 2014 09:53:42 +0200 Subject: [PATCH 13/60] Base all classes on `object` for attribute access. --- cython_numpy_auto/code_analyzer.py | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 3489372..9d376a4 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -14,7 +14,7 @@ } -class FunctionDoc: +class FunctionDoc(object): def __init__(self, doc): self.doc = doc.replace("**"," ").replace("/*"," ").replace("*/"," ") @@ -46,7 +46,7 @@ def output(self): def __repr__(self): return self.doc -class ArgumentDoc: +class ArgumentDoc(object): def __init__(self, doc): match = re.search("^ ([^ ]+)[ ]+([^ ]+)[ ]+(.+)", doc) @@ -62,7 +62,7 @@ def __init__(self, doc): def __repr__(self): return " {0:15} {1:15} {2}".format(self.name, self.type, self.doc) -class Argument: +class Argument(object): def __init__(self, definition, doc): self.__doc = doc @@ -121,7 +121,7 @@ def dtype(self): def __repr__(self): return "Argument('{0}', name='{1}', ctype='{2}', inout_state='{3}')".format(self.definition, self.name, self.ctype, self.inout_state) -class Function: +class Function(object): def __init__(self, name, source_path): self.name = name From 729d2e0e3f65deda215a7e78a50cfbebdb0c7100 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 30 Jul 2014 09:54:12 +0200 Subject: [PATCH 14/60] Turn `cython_numpy_auto` into a package --- cython_numpy_auto/__init__.py | 0 1 file changed, 0 insertions(+), 0 deletions(-) create mode 100644 cython_numpy_auto/__init__.py diff --git a/cython_numpy_auto/__init__.py b/cython_numpy_auto/__init__.py new file mode 100644 index 0000000..e69de29 From b6da8b8b53aea480f0bc94d4eeff3b4a90ccfb01 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 30 Jul 2014 09:55:08 +0200 Subject: [PATCH 15/60] Use jinja2 for code generation --- cython_numpy_auto/cython_generator.py | 67 ++++----------------------- cython_numpy_auto/erfa.pyx.in | 67 +++++++++++++++++++++++++++ 2 files changed, 76 insertions(+), 58 deletions(-) create mode 100644 cython_numpy_auto/erfa.pyx.in diff --git a/cython_numpy_auto/cython_generator.py b/cython_numpy_auto/cython_generator.py index a0d84a8..d8b4f43 100644 --- a/cython_numpy_auto/cython_generator.py +++ b/cython_numpy_auto/cython_generator.py @@ -1,62 +1,13 @@ -from code_analyzer import * +from jinja2 import Environment, PackageLoader +from cython_numpy_auto.code_analyzer import Function +env = Environment(loader=PackageLoader('cython_numpy_auto', '.')) +erfa_pyx_in = env.get_template('erfa.pyx.in') -def generate(era_func_names, source_path): +era_func_names = ["eraAtco13", "eraD2dtf", "eraAper"] +source_path = "../../erfa/src" - era_funcs = [Function(name, source_path) for name in era_func_names] - - wrapper = [] - wrapper.append("import numpy as np") - wrapper.append("cimport numpy as np") - wrapper.append("") - wrapper.append("np.import_array()") - wrapper.append("") - wrapper.append("__all__ = [{0}]".format(", ".join(["'{0}'".format(func.pyname) for func in era_funcs]+["'dt_eraASTROM'"]))) - wrapper.append("") - wrapper.append("cdef extern from 'erfa.h':") - wrapper.append(" struct eraASTROM:") - wrapper.append(" pass") - for func in era_funcs: - wrapper.append(" {0}".format(func.cfunc)) - wrapper.append("") - wrapper.append("dt_eraASTROM = np.dtype([('pmt','d'),('eb','d',(3,)),('eh','d',(3,)),('em','d'),('v','d',(3,)),('bm1 ','d'),('bpn','d',(3,3)),('along','d'),('phi','d'),('xpl','d'),('ypl','d'),('sphi','d'),('cphi','d'),('diurab','d'),('eral','d'),('refa','d'),('refb','d')], align=True)") - wrapper.append("") - for func in era_funcs: - wrapper.append("#=== {0} ===".format(func.name)) - wrapper.append("") - wrapper.append("def {0}({1}):".format(func.pyname, func.args_by_inout('in|inout','name',', '))) - wrapper.append(" ") - wrapper.append(" shape = np.broadcast({0}).shape".format(func.args_by_inout('in|inout','name',', '))) - for arg in func.args_by_inout('out'): - wrapper.append(" {0}_out = np.empty(shape, dtype={1})".format(arg.name, arg.dtype)) - for arg in func.args_by_inout('inout'): - wrapper.append(" {0}_out = np.array(np.broadcast_arrays({1})[{2}], dtype={3})".format(arg.name, func.args_by_inout('in|inout','name',', '), func.args.index(arg), arg.dtype)) - wrapper.append(" ") - wrapper.append(" cdef np.broadcast it = np.broadcast({0})".format(', '.join([arg.name if arg.is_in else (arg.name+"_out") for arg in func.args_by_inout('in|inout|out')]))) - wrapper.append(" ") - for arg in func.args_by_inout('in|out|inout'): - wrapper.append(" cdef {0} _{1}".format(arg.ctype_ptr, arg.name)) - wrapper.append(" ") - wrapper.append(" while np.PyArray_MultiIter_NOTDONE(it):") - wrapper.append(" ") - for arg in func.args_by_inout('in|out|inout'): - if arg.ctype_ptr[-1] == '*': - wrapper.append(" _{0} = (<{1}>np.PyArray_MultiIter_DATA(it, {2}))".format(arg.name, arg.ctype_ptr, func.args.index(arg))) - else: - wrapper.append(" _{0} = (<{1}*>np.PyArray_MultiIter_DATA(it, {2}))[0]".format(arg.name, arg.ctype_ptr, func.args.index(arg))) - wrapper.append(" ") - wrapper.append(" {0}({1})".format(func.name, ", ".join(["_"+name for name in func.args_by_inout('in|out|inout','name')]))) - wrapper.append(" ") - wrapper.append(" np.PyArray_MultiIter_NEXT(it)") - wrapper.append(" ") - wrapper.append(" return {0}".format(', '.join([arg+"_out" for arg in func.args_by_inout('out|inout','name')]))) - wrapper.append("") - - with open("./erfa.pyx","w") as f: - f.write("\n".join(wrapper)) +erfa_pyx = erfa_pyx_in.render(funcs=[Function(name, source_path) for name in era_func_names]) - -if __name__ == '__main__': - era_func_names = ["eraAtco13", "eraD2dtf", "eraAper"] - source_path = "../../erfa/src" - generate(era_func_names, source_path) +with open("erfa.pyx", "w") as f: + f.write(erfa_pyx) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in new file mode 100644 index 0000000..3673a32 --- /dev/null +++ b/cython_numpy_auto/erfa.pyx.in @@ -0,0 +1,67 @@ +import numpy as np +cimport numpy as np + +np.import_array() + +__all__ = ['{{ funcs|map(attribute='name')|join("', '") }}'] + + +cdef extern from "erfa.h": + struct eraASTROM: + pass +{%- for func in funcs %} + {{ func.cfunc }} +{%- endfor %} + +dt_eraASTROM = np.dtype([('pmt','d'), + ('eb','d',(3,)), + ('eh','d',(3,)), + ('em','d'), + ('v','d',(3,)), + ('bm1 ','d'), + ('bpn','d',(3,3)), + ('along','d'), + ('phi','d'), + ('xpl','d'), + ('ypl','d'), + ('sphi','d'), + ('cphi','d'), + ('diurab','d'), + ('eral','d'), + ('refa','d'), + ('refb','d')], align=True) + + +{% for func in funcs %} +def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}): + + shape = np.broadcast({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}).shape + {%- for arg in func.args_by_inout('out') %} + {{ arg.name }}_out = np.empty(shape, dtype={{ arg.dtype }}) + {%- endfor %} + {%- for arg in func.args_by_inout('inout') %} + {{ arg.name }}_out = np.array(np.broadcast_arrays({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }})[{{ func.args.index(arg) }}], dtype={{ arg.dtype }}) + {%- endfor %} + + cdef np.broadcast it = np.broadcast({{ func.args_by_inout('in')|map(attribute='name')|join(', ') }}, {{func.args_by_inout('inout|out')|map(attribute='name')|join('_out, ') }}_out) + {%- for arg in func.args_by_inout('in|inout|out') %} + cdef {{ arg.ctype_ptr }} _{{ arg.name }} + {%- endfor %} + + while np.PyArray_MultiIter_NOTDONE(it): + + {%- for arg in func.args_by_inout('in|inout|out') %} + {%- if arg.ctype_ptr[-1] == '*' %} + _{{ arg.name }} = (<{{ arg.ctype_ptr }}>np.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }})) + {%- else %} + _{{ arg.name }} = (<{{ arg.ctype_ptr }}*>np.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }}))[0] + {%- endif %} + {%- endfor %} + + {{ func.name }}(_{{ func.args_by_inout('in|out|inout')|map(attribute='name')|join(', _') }}) + + np.PyArray_MultiIter_NEXT(it) + + return {{ func.args_by_inout('inout|out')|map(attribute='name')|join('_out, ') }}_out + +{% endfor %} From 1ec74bcd06a25eebe340633a2cd4b3f14e9bae73 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:31:43 +0200 Subject: [PATCH 16/60] Add more ctype to dtype conversions --- cython_numpy_auto/code_analyzer.py | 17 +++++++++++------ 1 file changed, 11 insertions(+), 6 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 9d376a4..4344790 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -5,12 +5,17 @@ __all__ = ['Function'] -ctype_to_dtype = {'double' : "np.double", - 'double *' : "np.double", - 'int' : "np.int", - 'int *' : "np.int", - 'int[4]' : "np.dtype([('', 'i')]*4)", - 'eraASTROM *': "dt_eraASTROM", +ctype_to_dtype = {'double' : "numpy.double", + 'double *' : "numpy.double", + 'int' : "numpy.int", + 'int *' : "numpy.int", + 'int[4]' : "numpy.dtype([('', 'i', (4,))])", + 'double[2]' : "numpy.dtype([('', 'd', (2,))])", + 'double[3]' : "numpy.dtype([('p', 'd', (3,))])", + 'double[2][3]' : "numpy.dtype([('pv', 'd', (2,3))])", + 'double[3][3]' : "numpy.dtype([('r', 'd', (3,3))])", + 'eraASTROM *' : "dt_eraASTROM", + 'char *' : "numpy.dtype('S1')", } From 940753a2d0db12e38726bd308c9a8c4073059b95 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:32:30 +0200 Subject: [PATCH 17/60] Deal with empty input/output sections in doc --- cython_numpy_auto/code_analyzer.py | 24 ++++++++++++++---------- 1 file changed, 14 insertions(+), 10 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 4344790..9fa218d 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -30,22 +30,26 @@ def __init__(self, doc): def input(self): if self.__input is None: self.__input = [] - __input = re.search("Given:\n(.+?) \n", self.doc, re.DOTALL).group(1) - for i in __input.split("\n"): - arg_doc = ArgumentDoc(i) - if arg_doc.name is not None: - self.__input.append(arg_doc) + result = re.search("Given[^\n]*:\n(.+?) \n", self.doc, re.DOTALL) + if result is not None: + __input = result.group(1) + for i in __input.split("\n"): + arg_doc = ArgumentDoc(i) + if arg_doc.name is not None: + self.__input.append(arg_doc) return self.__input @property def output(self): if self.__output is None: self.__output = [] - __output = re.search("Returned:\n(.+?) \n", self.doc, re.DOTALL).group(1) - for i in __output.split("\n"): - arg_doc = ArgumentDoc(i) - if arg_doc.name is not None: - self.__output.append(arg_doc) + result = re.search("Returned:\n(.+?) \n", self.doc, re.DOTALL) + if result is not None: + __output = result.group(1) + for i in __output.split("\n"): + arg_doc = ArgumentDoc(i) + if arg_doc.name is not None: + self.__output.append(arg_doc) return self.__output def __repr__(self): From 41cbfe12e081e47d6d9487b0a36f5f6953058e0d Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:33:21 +0200 Subject: [PATCH 18/60] Argument documentations have variable leading spaces --- cython_numpy_auto/code_analyzer.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 9fa218d..69ade73 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -58,7 +58,7 @@ def __repr__(self): class ArgumentDoc(object): def __init__(self, doc): - match = re.search("^ ([^ ]+)[ ]+([^ ]+)[ ]+(.+)", doc) + match = re.search("^ +([^ ]+)[ ]+([^ ]+)[ ]+(.+)", doc) if match is not None: self.name = match.group(1) self.type = match.group(2) From 7506ff08d0c3a3053b3bcc8c07a92fdced7fd522 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:34:00 +0200 Subject: [PATCH 19/60] Deal with multidimensional C arrays --- cython_numpy_auto/code_analyzer.py | 11 +++++------ 1 file changed, 5 insertions(+), 6 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 69ade73..da1a49e 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -77,15 +77,14 @@ def __init__(self, definition, doc): self.__doc = doc self.__inout_state = None self.definition = definition.strip() - if self.definition[-1] == "]": - self.ctype, self.name = self.definition.split(" ",1) - self.name, arr = self.name.split("[") - self.ctype += ("["+arr) - elif "*" in self.definition: + if "*" in self.definition: self.ctype, self.name = self.definition.split("*", 1) self.ctype += "*" else: - self.ctype, self.name = self.definition.split(" ", 1) + self.ctype, self.name = self.definition.rsplit(" ", 1) + if "[" in self.name: + self.name, arr = self.name.split("[", 1) + self.ctype += ("["+arr) @property def inout_state(self): From 7bbcabb4e44f36000c343855d49217726fedf7d2 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:34:31 +0200 Subject: [PATCH 20/60] Pointer arguments are always of input type --- cython_numpy_auto/code_analyzer.py | 2 ++ 1 file changed, 2 insertions(+) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index da1a49e..701985a 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -93,6 +93,8 @@ def inout_state(self): for i in self.__doc.input: if self.name in i.name.split(','): self.__inout_state = 'in' + if "*" in self.ctype_ptr: + self.__inout_state = 'in' for o in self.__doc.output: if self.name in o.name.split(','): if self.__inout_state == 'in': From 922cba1bab5bd2e0b7443888ae7361d346d57f56 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:34:58 +0200 Subject: [PATCH 21/60] Remove dead code: is_in(), is_out(), is_inout() --- cython_numpy_auto/code_analyzer.py | 12 ------------ 1 file changed, 12 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 701985a..3c3986b 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -103,18 +103,6 @@ def inout_state(self): self.__inout_state = 'out' return self.__inout_state - @property - def is_in(self): - return self.inout_state == 'in' - - @property - def is_out(self): - return self.inout_state == 'out' - - @property - def is_inout(self): - return self.inout_state == 'inout' - @property def ctype_ptr(self): if self.ctype[-1] == ']': From fa3612b9a9d6128e1fd0690f2adc2b893b53b680 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:35:50 +0200 Subject: [PATCH 22/60] Implement Return() class to use along regular arguments --- cython_numpy_auto/code_analyzer.py | 12 ++++++++++++ 1 file changed, 12 insertions(+) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 3c3986b..e00364d 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -119,6 +119,18 @@ def dtype(self): def __repr__(self): return "Argument('{0}', name='{1}', ctype='{2}', inout_state='{3}')".format(self.definition, self.name, self.ctype, self.inout_state) +class Return(object): + + def __init__(self, ctype, doc): + self.name = 'ret' + self.ctype = ctype + self.inout_state = 'ret' + self.ctype_ptr = ctype + + @property + def dtype(self): + return ctype_to_dtype[self.ctype] + class Function(object): def __init__(self, name, source_path): From 1ebb796c9ad46e28416d03be8c847d19cd3204f8 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:37:12 +0200 Subject: [PATCH 23/60] Spaces between function name and opening parenthesis --- cython_numpy_auto/code_analyzer.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index e00364d..2e7cc18 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -138,7 +138,7 @@ def __init__(self, name, source_path): self.pyname = name.split('era')[-1].lower() self.filename = name.split("era")[-1].lower()+".c" self.filepath = os.path.join(os.path.normpath(source_path), self.filename) - pattern = "\n([^\n]+{0}\([^)]+\)).+?(/\*.+?\*/)".format(name) + pattern = "\n([^\n]+{0} ?\([^)]+\)).+?(/\*.+?\*/)".format(name) p = re.compile(pattern, flags=re.DOTALL|re.MULTILINE) with open(self.filepath) as f: search = p.search(f.read()) From 3d3521823a013b795ecca51ec84546db32abad3e Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:37:54 +0200 Subject: [PATCH 24/60] Add returned parameter at end of arguments list --- cython_numpy_auto/code_analyzer.py | 3 +++ 1 file changed, 3 insertions(+) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 2e7cc18..5306e96 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -147,6 +147,9 @@ def __init__(self, name, source_path): self.args = [] for arg in re.search("\(([^)]+)\)", self.cfunc, flags=re.MULTILINE|re.DOTALL).group(1).split(','): self.args.append(Argument(arg, self.__doc)) + self.ret = re.search("^(.*){0}".format(name), self.cfunc).group(1).strip() + if self.ret == 'double': + self.args.append(Return(self.ret, self.__doc)) def args_by_inout(self, inout_filter, prop=None, join=None): result = [] From 0de22b98a30c1ad5eac4e0e75b2158cdbcd47a66 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:41:07 +0200 Subject: [PATCH 25/60] np instead of numpy collides with c code arguments --- cython_numpy_auto/erfa.pyx.in | 24 ++++++++++++------------ 1 file changed, 12 insertions(+), 12 deletions(-) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in index 3673a32..9b037e0 100644 --- a/cython_numpy_auto/erfa.pyx.in +++ b/cython_numpy_auto/erfa.pyx.in @@ -1,7 +1,7 @@ -import numpy as np -cimport numpy as np +import numpy +cimport numpy -np.import_array() +numpy.import_array() __all__ = ['{{ funcs|map(attribute='name')|join("', '") }}'] @@ -13,7 +13,7 @@ cdef extern from "erfa.h": {{ func.cfunc }} {%- endfor %} -dt_eraASTROM = np.dtype([('pmt','d'), +dt_eraASTROM = numpy.dtype([('pmt','d'), ('eb','d',(3,)), ('eh','d',(3,)), ('em','d'), @@ -35,32 +35,32 @@ dt_eraASTROM = np.dtype([('pmt','d'), {% for func in funcs %} def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}): - shape = np.broadcast({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}).shape + shape = numpy.broadcast({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}).shape {%- for arg in func.args_by_inout('out') %} - {{ arg.name }}_out = np.empty(shape, dtype={{ arg.dtype }}) + {{ arg.name }}_out = numpy.empty(shape, dtype={{ arg.dtype }}) {%- endfor %} {%- for arg in func.args_by_inout('inout') %} - {{ arg.name }}_out = np.array(np.broadcast_arrays({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }})[{{ func.args.index(arg) }}], dtype={{ arg.dtype }}) + {{ arg.name }}_out = numpy.array(numpy.broadcast_arrays({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }})[{{ func.args.index(arg) }}], dtype={{ arg.dtype }}) {%- endfor %} - cdef np.broadcast it = np.broadcast({{ func.args_by_inout('in')|map(attribute='name')|join(', ') }}, {{func.args_by_inout('inout|out')|map(attribute='name')|join('_out, ') }}_out) + cdef numpy.broadcast it = numpy.broadcast({{ func.args_by_inout('in')|map(attribute='name')|join(', ') }}, {{func.args_by_inout('inout|out')|map(attribute='name')|join('_out, ') }}_out) {%- for arg in func.args_by_inout('in|inout|out') %} cdef {{ arg.ctype_ptr }} _{{ arg.name }} {%- endfor %} - while np.PyArray_MultiIter_NOTDONE(it): + while numpy.PyArray_MultiIter_NOTDONE(it): {%- for arg in func.args_by_inout('in|inout|out') %} {%- if arg.ctype_ptr[-1] == '*' %} - _{{ arg.name }} = (<{{ arg.ctype_ptr }}>np.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }})) + _{{ arg.name }} = (<{{ arg.ctype_ptr }}>numpy.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }})) {%- else %} - _{{ arg.name }} = (<{{ arg.ctype_ptr }}*>np.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }}))[0] + _{{ arg.name }} = (<{{ arg.ctype_ptr }}*>numpy.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }}))[0] {%- endif %} {%- endfor %} {{ func.name }}(_{{ func.args_by_inout('in|out|inout')|map(attribute='name')|join(', _') }}) - np.PyArray_MultiIter_NEXT(it) + numpy.PyArray_MultiIter_NEXT(it) return {{ func.args_by_inout('inout|out')|map(attribute='name')|join('_out, ') }}_out From 763c08c49a156a30664943d9513564061493f253 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:44:41 +0200 Subject: [PATCH 26/60] Three new jinja filters: prefix, postfix, surround --- cython_numpy_auto/cython_generator.py | 11 +++++++++++ 1 file changed, 11 insertions(+) diff --git a/cython_numpy_auto/cython_generator.py b/cython_numpy_auto/cython_generator.py index d8b4f43..faf36b4 100644 --- a/cython_numpy_auto/cython_generator.py +++ b/cython_numpy_auto/cython_generator.py @@ -2,6 +2,17 @@ from cython_numpy_auto.code_analyzer import Function env = Environment(loader=PackageLoader('cython_numpy_auto', '.')) + +def prefix(a_list, pre): + return [pre+'{0}'.format(an_element) for an_element in a_list] +def postfix(a_list, post): + return ['{0}'.format(an_element)+post for an_element in a_list] +def surround(a_list, pre, post): + return [pre+'{0}'.format(an_element)+post for an_element in a_list] +env.filters['prefix'] = prefix +env.filters['postfix'] = postfix +env.filters['surround'] = surround + erfa_pyx_in = env.get_template('erfa.pyx.in') era_func_names = ["eraAtco13", "eraD2dtf", "eraAper"] From d03c0c1e34beb4ab393b6cead6e8739d243d4307 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:45:50 +0200 Subject: [PATCH 27/60] Process all ERFA functions --- cython_numpy_auto/cython_generator.py | 22 ++++++++++++++++------ 1 file changed, 16 insertions(+), 6 deletions(-) diff --git a/cython_numpy_auto/cython_generator.py b/cython_numpy_auto/cython_generator.py index faf36b4..8809d1b 100644 --- a/cython_numpy_auto/cython_generator.py +++ b/cython_numpy_auto/cython_generator.py @@ -1,6 +1,10 @@ +import re from jinja2 import Environment, PackageLoader from cython_numpy_auto.code_analyzer import Function +ERFA_SOURCES = "../../erfa/src" + +#Prepare the jinja2 templating environment env = Environment(loader=PackageLoader('cython_numpy_auto', '.')) def prefix(a_list, pre): @@ -15,10 +19,16 @@ def surround(a_list, pre, post): erfa_pyx_in = env.get_template('erfa.pyx.in') -era_func_names = ["eraAtco13", "eraD2dtf", "eraAper"] -source_path = "../../erfa/src" - -erfa_pyx = erfa_pyx_in.render(funcs=[Function(name, source_path) for name in era_func_names]) -with open("erfa.pyx", "w") as f: - f.write(erfa_pyx) +#Extract all the ERFA function names from erfa.h +with open(ERFA_SOURCES+"/erfa.h", "r") as f: + func_names = re.findall(' (\w+)\(.*?\);', f.read(), flags=re.DOTALL) + funcs = [] + for name in func_names: + print("Parsing {0}...".format(name)) + funcs.append(Function(name, ERFA_SOURCES)) + print("Done!") + #Render the template and save + erfa_pyx = erfa_pyx_in.render(funcs=funcs) + with open("erfa.pyx", "w") as f: + f.write(erfa_pyx) From c9be43ef2567f78c507b7bff590d4d012e999899 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:47:32 +0200 Subject: [PATCH 28/60] Add eraLDBODY structure support --- cython_numpy_auto/erfa.pyx.in | 2 ++ 1 file changed, 2 insertions(+) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in index 9b037e0..66b068c 100644 --- a/cython_numpy_auto/erfa.pyx.in +++ b/cython_numpy_auto/erfa.pyx.in @@ -9,6 +9,8 @@ __all__ = ['{{ funcs|map(attribute='name')|join("', '") }}'] cdef extern from "erfa.h": struct eraASTROM: pass + struct eraLDBODY: + pass {%- for func in funcs %} {{ func.cfunc }} {%- endfor %} From 56d80cad668c14fd091efd4b96da6973d86642dd Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:51:35 +0200 Subject: [PATCH 29/60] Use prefix, postfix, and surround in template. --- cython_numpy_auto/erfa.pyx.in | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in index 66b068c..be89fc6 100644 --- a/cython_numpy_auto/erfa.pyx.in +++ b/cython_numpy_auto/erfa.pyx.in @@ -3,7 +3,7 @@ cimport numpy numpy.import_array() -__all__ = ['{{ funcs|map(attribute='name')|join("', '") }}'] +__all__ = [{{ funcs|map(attribute='name')|surround("'","'")|join(", ") }}] cdef extern from "erfa.h": @@ -45,7 +45,7 @@ def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|jo {{ arg.name }}_out = numpy.array(numpy.broadcast_arrays({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }})[{{ func.args.index(arg) }}], dtype={{ arg.dtype }}) {%- endfor %} - cdef numpy.broadcast it = numpy.broadcast({{ func.args_by_inout('in')|map(attribute='name')|join(', ') }}, {{func.args_by_inout('inout|out')|map(attribute='name')|join('_out, ') }}_out) + cdef numpy.broadcast it = numpy.broadcast({{ (func.args_by_inout('in')|map(attribute='name')|list + func.args_by_inout('inout|out')|map(attribute='name')|postfix('_out')|list) |join(', ') }}) {%- for arg in func.args_by_inout('in|inout|out') %} cdef {{ arg.ctype_ptr }} _{{ arg.name }} {%- endfor %} @@ -60,10 +60,10 @@ def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|jo {%- endif %} {%- endfor %} - {{ func.name }}(_{{ func.args_by_inout('in|out|inout')|map(attribute='name')|join(', _') }}) + {{ func.name }}({{ func.args_by_inout('in|out|inout')|map(attribute='name')|prefix('_')|join(', ') }}) numpy.PyArray_MultiIter_NEXT(it) - return {{ func.args_by_inout('inout|out')|map(attribute='name')|join('_out, ') }}_out + return {{ func.args_by_inout('inout|out')|map(attribute='name')|postfix('_out')|join(', ') }} {% endfor %} From e3ccefdd5d7628fd29cf19c6b82dbc56bada6ea2 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:52:08 +0200 Subject: [PATCH 30/60] Re-generate extern function signature with pointers --- cython_numpy_auto/erfa.pyx.in | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in index be89fc6..7eddf13 100644 --- a/cython_numpy_auto/erfa.pyx.in +++ b/cython_numpy_auto/erfa.pyx.in @@ -12,7 +12,7 @@ cdef extern from "erfa.h": struct eraLDBODY: pass {%- for func in funcs %} - {{ func.cfunc }} + {{ func.ret }} {{ func.name }}({{ func.args_by_inout('in|inout|out')|map(attribute='ctype_ptr')|join(', ') }}) {%- endfor %} dt_eraASTROM = numpy.dtype([('pmt','d'), From e641fd4fedf275e77e97ab6c6b4c558772188f3a Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:55:14 +0200 Subject: [PATCH 31/60] Add returned values --- cython_numpy_auto/erfa.pyx.in | 14 +++++++++----- 1 file changed, 9 insertions(+), 5 deletions(-) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in index 7eddf13..44348e0 100644 --- a/cython_numpy_auto/erfa.pyx.in +++ b/cython_numpy_auto/erfa.pyx.in @@ -38,15 +38,15 @@ dt_eraASTROM = numpy.dtype([('pmt','d'), def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}): shape = numpy.broadcast({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}).shape - {%- for arg in func.args_by_inout('out') %} + {%- for arg in func.args_by_inout('out|ret') %} {{ arg.name }}_out = numpy.empty(shape, dtype={{ arg.dtype }}) {%- endfor %} {%- for arg in func.args_by_inout('inout') %} {{ arg.name }}_out = numpy.array(numpy.broadcast_arrays({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }})[{{ func.args.index(arg) }}], dtype={{ arg.dtype }}) {%- endfor %} - cdef numpy.broadcast it = numpy.broadcast({{ (func.args_by_inout('in')|map(attribute='name')|list + func.args_by_inout('inout|out')|map(attribute='name')|postfix('_out')|list) |join(', ') }}) - {%- for arg in func.args_by_inout('in|inout|out') %} + cdef numpy.broadcast it = numpy.broadcast({{ (func.args_by_inout('in')|map(attribute='name')|list + func.args_by_inout('inout|out|ret')|map(attribute='name')|postfix('_out')|list) |join(', ') }}) + {%- for arg in func.args_by_inout('in|inout|out|ret') %} cdef {{ arg.ctype_ptr }} _{{ arg.name }} {%- endfor %} @@ -60,10 +60,14 @@ def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|jo {%- endif %} {%- endfor %} - {{ func.name }}({{ func.args_by_inout('in|out|inout')|map(attribute='name')|prefix('_')|join(', ') }}) + {{ func.args_by_inout('ret')|map(attribute='name')|surround('_',' = ')|join }}{{ func.name }}({{ func.args_by_inout('in|out|inout')|map(attribute='name')|prefix('_')|join(', ') }}) + + {%- for arg in func.args_by_inout('ret') %} + (<{{ arg.ctype_ptr }}*>numpy.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }}))[0] = _{{ arg.name }} + {%- endfor %} numpy.PyArray_MultiIter_NEXT(it) - return {{ func.args_by_inout('inout|out')|map(attribute='name')|postfix('_out')|join(', ') }} + return {{ func.args_by_inout('inout|out|ret')|map(attribute='name')|postfix('_out')|join(', ') }} {% endfor %} From dade3e0d1e537f8f949c372b6fbe95ce57544a5f Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 08:55:50 +0200 Subject: [PATCH 32/60] Functions with no inputs do not need vectorisation --- cython_numpy_auto/erfa.pyx.in | 8 ++++++++ 1 file changed, 8 insertions(+) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in index 44348e0..c9831e4 100644 --- a/cython_numpy_auto/erfa.pyx.in +++ b/cython_numpy_auto/erfa.pyx.in @@ -37,6 +37,7 @@ dt_eraASTROM = numpy.dtype([('pmt','d'), {% for func in funcs %} def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}): + {%- if func.args_by_inout('in|inout') %} shape = numpy.broadcast({{ func.args_by_inout('in|inout')|map(attribute='name')|join(', ') }}).shape {%- for arg in func.args_by_inout('out|ret') %} {{ arg.name }}_out = numpy.empty(shape, dtype={{ arg.dtype }}) @@ -67,6 +68,13 @@ def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|jo {%- endfor %} numpy.PyArray_MultiIter_NEXT(it) + {%- else %} + {%- for arg in func.args_by_inout('out') %} + {{ arg.name }}_out = numpy.empty(shape, dtype={{ arg.dtype }}) + {%- endfor %} + + {{ func.name }}({{ func.args_by_inout('in|out|inout')|map(attribute='name')|prefix('_')|join(', ') }}) + {%- endif %} return {{ func.args_by_inout('inout|out|ret')|map(attribute='name')|postfix('_out')|join(', ') }} From e16ff9aa70640518f32bf9ae62bbd30776ba0e8d Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 09:35:43 +0200 Subject: [PATCH 33/60] Improvements for no input functions --- cython_numpy_auto/erfa.pyx.in | 20 +++++++++++++++++++- 1 file changed, 19 insertions(+), 1 deletion(-) diff --git a/cython_numpy_auto/erfa.pyx.in b/cython_numpy_auto/erfa.pyx.in index c9831e4..3b65992 100644 --- a/cython_numpy_auto/erfa.pyx.in +++ b/cython_numpy_auto/erfa.pyx.in @@ -69,11 +69,29 @@ def {{ func.pyname }}({{ func.args_by_inout('in|inout')|map(attribute='name')|jo numpy.PyArray_MultiIter_NEXT(it) {%- else %} + shape = [] {%- for arg in func.args_by_inout('out') %} {{ arg.name }}_out = numpy.empty(shape, dtype={{ arg.dtype }}) {%- endfor %} - {{ func.name }}({{ func.args_by_inout('in|out|inout')|map(attribute='name')|prefix('_')|join(', ') }}) + cdef numpy.broadcast it = numpy.broadcast({{ (func.args_by_inout('in')|map(attribute='name')|list + func.args_by_inout('inout|out|ret')|map(attribute='name')|postfix('_out')|list) |join(', ') }}) + {%- for arg in func.args_by_inout('out') %} + cdef {{ arg.ctype_ptr }} _{{ arg.name }} + {%- endfor %} + + while numpy.PyArray_MultiIter_NOTDONE(it): + + {%- for arg in func.args_by_inout('in|inout|out') %} + {%- if arg.ctype_ptr[-1] == '*' %} + _{{ arg.name }} = (<{{ arg.ctype_ptr }}>numpy.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }})) + {%- else %} + _{{ arg.name }} = (<{{ arg.ctype_ptr }}*>numpy.PyArray_MultiIter_DATA(it, {{ func.args.index(arg) }}))[0] + {%- endif %} + {%- endfor %} + + {{ func.name }}({{ func.args_by_inout('in|out|inout')|map(attribute='name')|prefix('_')|join(', ') }}) + + numpy.PyArray_MultiIter_NEXT(it) {%- endif %} return {{ func.args_by_inout('inout|out|ret')|map(attribute='name')|postfix('_out')|join(', ') }} From 080325d11df48524249ddc6ae4405566b221e1a0 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 09:36:16 +0200 Subject: [PATCH 34/60] Pointer arguments are actually not necessarily inputs --- cython_numpy_auto/code_analyzer.py | 2 -- 1 file changed, 2 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index 5306e96..fb4f405 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -93,8 +93,6 @@ def inout_state(self): for i in self.__doc.input: if self.name in i.name.split(','): self.__inout_state = 'in' - if "*" in self.ctype_ptr: - self.__inout_state = 'in' for o in self.__doc.output: if self.name in o.name.split(','): if self.__inout_state == 'in': From 9981fb059ed90af7785111cde9e56f9e41349db6 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 10:01:18 +0200 Subject: [PATCH 35/60] Deal with "Given and returned" documentation sections --- cython_numpy_auto/code_analyzer.py | 26 ++++++++++++++++++++++---- 1 file changed, 22 insertions(+), 4 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index fb4f405..e2555ce 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -30,9 +30,18 @@ def __init__(self, doc): def input(self): if self.__input is None: self.__input = [] - result = re.search("Given[^\n]*:\n(.+?) \n", self.doc, re.DOTALL) + result = re.search("Given([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) if result is not None: - __input = result.group(1) + print(result.group(1)) + __input = result.group(2) + for i in __input.split("\n"): + arg_doc = ArgumentDoc(i) + if arg_doc.name is not None: + self.__input.append(arg_doc) + result = re.search("Given and returned([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) + if result is not None: + print(result.group(1)) + __input = result.group(2) for i in __input.split("\n"): arg_doc = ArgumentDoc(i) if arg_doc.name is not None: @@ -43,9 +52,18 @@ def input(self): def output(self): if self.__output is None: self.__output = [] - result = re.search("Returned:\n(.+?) \n", self.doc, re.DOTALL) + result = re.search("Returned([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) + if result is not None: + print(result.group(1)) + __output = result.group(2) + for i in __output.split("\n"): + arg_doc = ArgumentDoc(i) + if arg_doc.name is not None: + self.__output.append(arg_doc) + result = re.search("Given and returned([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) if result is not None: - __output = result.group(1) + print(result.group(1)) + __output = result.group(2) for i in __output.split("\n"): arg_doc = ArgumentDoc(i) if arg_doc.name is not None: From a85fcf9bcd5cd494dbc2ce52caf8ef03af57dfa3 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 10:39:31 +0200 Subject: [PATCH 36/60] Remove unwanted print statements --- cython_numpy_auto/code_analyzer.py | 4 ---- 1 file changed, 4 deletions(-) diff --git a/cython_numpy_auto/code_analyzer.py b/cython_numpy_auto/code_analyzer.py index e2555ce..bde0348 100644 --- a/cython_numpy_auto/code_analyzer.py +++ b/cython_numpy_auto/code_analyzer.py @@ -32,7 +32,6 @@ def input(self): self.__input = [] result = re.search("Given([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) if result is not None: - print(result.group(1)) __input = result.group(2) for i in __input.split("\n"): arg_doc = ArgumentDoc(i) @@ -40,7 +39,6 @@ def input(self): self.__input.append(arg_doc) result = re.search("Given and returned([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) if result is not None: - print(result.group(1)) __input = result.group(2) for i in __input.split("\n"): arg_doc = ArgumentDoc(i) @@ -54,7 +52,6 @@ def output(self): self.__output = [] result = re.search("Returned([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) if result is not None: - print(result.group(1)) __output = result.group(2) for i in __output.split("\n"): arg_doc = ArgumentDoc(i) @@ -62,7 +59,6 @@ def output(self): self.__output.append(arg_doc) result = re.search("Given and returned([^\n]*):\n(.+?) \n", self.doc, re.DOTALL) if result is not None: - print(result.group(1)) __output = result.group(2) for i in __output.split("\n"): arg_doc = ArgumentDoc(i) From c43c83cd3f1e37d680cda0714a509320be6c6260 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 31 Jul 2014 10:41:02 +0200 Subject: [PATCH 37/60] Identify sections/subsections and process only Astronomy section --- cython_numpy_auto/cython_generator.py | 17 ++++++++++++----- 1 file changed, 12 insertions(+), 5 deletions(-) diff --git a/cython_numpy_auto/cython_generator.py b/cython_numpy_auto/cython_generator.py index 8809d1b..3b227fe 100644 --- a/cython_numpy_auto/cython_generator.py +++ b/cython_numpy_auto/cython_generator.py @@ -22,12 +22,19 @@ def surround(a_list, pre, post): #Extract all the ERFA function names from erfa.h with open(ERFA_SOURCES+"/erfa.h", "r") as f: - func_names = re.findall(' (\w+)\(.*?\);', f.read(), flags=re.DOTALL) + + erfa_h = f.read() + funcs = [] - for name in func_names: - print("Parsing {0}...".format(name)) - funcs.append(Function(name, ERFA_SOURCES)) - print("Done!") + section_subsection_functions = re.findall('/\* (\w*)/(\w*) \*/\n(.*?)\n\n', erfa_h, flags=re.DOTALL|re.MULTILINE) + for section, subsection, functions in section_subsection_functions: + print("{0}.{1}".format(section, subsection)) + if section == "Astronomy": + func_names = re.findall(' (\w+)\(.*?\);', functions, flags=re.DOTALL) + for name in func_names: + print("{0}.{1}.{2}...".format(section, subsection, name)) + funcs.append(Function(name, ERFA_SOURCES)) + print("Done!") #Render the template and save erfa_pyx = erfa_pyx_in.render(funcs=funcs) with open("erfa.pyx", "w") as f: From 2dd86f006d9507efbb27e703ff98d2ac7ab3ee3f Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 23 Oct 2014 14:46:13 +0200 Subject: [PATCH 38/60] First commit using nditer --- nditer/erfa.pyx | 68 +++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 68 insertions(+) create mode 100644 nditer/erfa.pyx diff --git a/nditer/erfa.pyx b/nditer/erfa.pyx new file mode 100644 index 0000000..bdd3198 --- /dev/null +++ b/nditer/erfa.pyx @@ -0,0 +1,68 @@ +import numpy as np +cimport numpy as np +from cpython.ref cimport PyObject + +np.import_array() + +__all__ = ['ab'] + + +cdef extern from "erfa.h": + void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) + +ctypedef void NpyIter + +cdef extern from "numpy/arrayobject.h": + char** GetDataPtrArray "NpyIter_GetDataPtrArray" (NpyIter* iter) + +ctypedef struct NewNpyArrayIterObject: + PyObject base + NpyIter *iter + +cdef inline NpyIter* GetNpyIter(object iter): + return (iter).iter + + +#void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) +def ab(pnat, v, s, bm1): + + #Turn all inputs into arrays of at least dimension 1 + pnat_in = np.array(pnat) + v_in = np.array(v) + s_in = np.array(s) + bm1_in = np.array(bm1) + s_in = s_in.reshape(s_in.shape+(1,)) + bm1_in = bm1_in.reshape(bm1_in.shape+(1,)) + + #Create the iterator, broadcasting on all but the last dimension + ndim = max([arr.ndim for arr in [pnat_in, v_in, s_in, bm1_in]]) + op_axes = [[-1]*(ndim-arr.ndim)+range(arr.ndim-1) for arr in [pnat_in, v_in, s_in, bm1_in]] + it = np.nditer([pnat_in, v_in, s_in, bm1_in], op_axes=op_axes, flags=["multi_index"]) + + #Create the output array, based on the iterator shape, adding the last dimension + ppr_out = np.empty(it.shape[::-1]+(3,), dtype=np.double) + print(ppr_out.shape) + + cdef double* pnat_ + cdef double* v_ + cdef double s_ + cdef double bm1_ + cdef double ppr_[3] + + #Iterate + cdef char** dataptrarray + cdef char* data + dataptrarray = GetDataPtrArray(GetNpyIter(it)) + while not it.finished: + #it.debug_print() + pnat_ = ((dataptrarray[0])) + v_ = ((dataptrarray[1])) + s_ = ((dataptrarray[2]))[0] + bm1_ = ((dataptrarray[3]))[0] + eraAb(pnat_, v_, s_, bm1_, ppr_) + print(it.multi_index) + ppr_out[it.multi_index[::-1]+(0,)] = ppr_[0] + ppr_out[it.multi_index[::-1]+(1,)] = ppr_[1] + ppr_out[it.multi_index[::-1]+(2,)] = ppr_[2] + it.iternext() + return ppr_out From e555bcff193dc91beb8df650f776b0e00330e2bb Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Thu, 23 Oct 2014 14:49:12 +0200 Subject: [PATCH 39/60] Updated setup.py --- setup.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/setup.py b/setup.py index 4b7a75a..ec52c41 100644 --- a/setup.py +++ b/setup.py @@ -16,8 +16,14 @@ include_dirs = [np.get_include(), '/opt/local/include'], library_dirs = ['/opt/local/lib']) +nditer_ext = Extension("erfa.nditer", + ["nditer/erfa.pyx"], + libraries=['erfa'], + include_dirs = [np.get_include(), '/opt/local/include'], + library_dirs = ['/opt/local/lib']) + setup(name = "erfa", cmdclass = { "build_ext": build_ext }, packages = ["erfa"], package_dir = {"erfa":"."}, - ext_modules = [cython_numpy_ext, cython_numpy_auto_ext]) + ext_modules = [cython_numpy_ext, cython_numpy_auto_ext, nditer_ext]) From 2769e3fa8ac0ae551fe7f73f20f4611d486b6c75 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Mon, 3 Nov 2014 09:57:25 +0100 Subject: [PATCH 40/60] Make `ppr` part of the iterator --- nditer/erfa.pyx | 31 ++++++++++++------------------- 1 file changed, 12 insertions(+), 19 deletions(-) diff --git a/nditer/erfa.pyx b/nditer/erfa.pyx index bdd3198..f1bbd01 100644 --- a/nditer/erfa.pyx +++ b/nditer/erfa.pyx @@ -26,43 +26,36 @@ cdef inline NpyIter* GetNpyIter(object iter): #void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) def ab(pnat, v, s, bm1): - #Turn all inputs into arrays of at least dimension 1 + #Turn all inputs into arrays pnat_in = np.array(pnat) v_in = np.array(v) s_in = np.array(s) bm1_in = np.array(bm1) - s_in = s_in.reshape(s_in.shape+(1,)) - bm1_in = bm1_in.reshape(bm1_in.shape+(1,)) - #Create the iterator, broadcasting on all but the last dimension - ndim = max([arr.ndim for arr in [pnat_in, v_in, s_in, bm1_in]]) - op_axes = [[-1]*(ndim-arr.ndim)+range(arr.ndim-1) for arr in [pnat_in, v_in, s_in, bm1_in]] - it = np.nditer([pnat_in, v_in, s_in, bm1_in], op_axes=op_axes, flags=["multi_index"]) + #Create the output array, based on the broadcasted shape, adding the last dimension if needed + shape = np.broadcast(pnat_in[...,0], v_in[...,0], s_in, bm1_in).shape + ppr_out = np.empty(shape+(3,), dtype=np.double) - #Create the output array, based on the iterator shape, adding the last dimension - ppr_out = np.empty(it.shape[::-1]+(3,), dtype=np.double) - print(ppr_out.shape) + #Create the iterator, broadcasting on all but the last dimension + ndim = max([arr.ndim for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in]]) + op_axes = [[-1]*(ndim-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] + it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes, flags=["multi_index"]) cdef double* pnat_ cdef double* v_ cdef double s_ cdef double bm1_ - cdef double ppr_[3] + cdef double* ppr_ #Iterate - cdef char** dataptrarray - cdef char* data - dataptrarray = GetDataPtrArray(GetNpyIter(it)) + cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) while not it.finished: - #it.debug_print() pnat_ = ((dataptrarray[0])) v_ = ((dataptrarray[1])) s_ = ((dataptrarray[2]))[0] bm1_ = ((dataptrarray[3]))[0] + ppr_ = ((dataptrarray[4])) eraAb(pnat_, v_, s_, bm1_, ppr_) - print(it.multi_index) - ppr_out[it.multi_index[::-1]+(0,)] = ppr_[0] - ppr_out[it.multi_index[::-1]+(1,)] = ppr_[1] - ppr_out[it.multi_index[::-1]+(2,)] = ppr_[2] it.iternext() + return ppr_out From 87c3d80971644964b7b09d26638e062feb731237 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 4 Nov 2014 17:11:50 +0100 Subject: [PATCH 41/60] Iterate with C rather than python --- nditer/erfa.pyx | 8 ++++++-- 1 file changed, 6 insertions(+), 2 deletions(-) diff --git a/nditer/erfa.pyx b/nditer/erfa.pyx index f1bbd01..ab9a673 100644 --- a/nditer/erfa.pyx +++ b/nditer/erfa.pyx @@ -11,8 +11,10 @@ cdef extern from "erfa.h": void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) ctypedef void NpyIter +ctypedef int (*IterNextFunc)(NpyIter * iter) nogil cdef extern from "numpy/arrayobject.h": + IterNextFunc GetIterNext "NpyIter_GetIterNext" (NpyIter *iter, char **) char** GetDataPtrArray "NpyIter_GetDataPtrArray" (NpyIter* iter) ctypedef struct NewNpyArrayIterObject: @@ -49,13 +51,15 @@ def ab(pnat, v, s, bm1): #Iterate cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) - while not it.finished: + cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) + cdef int status = 1 + while status: pnat_ = ((dataptrarray[0])) v_ = ((dataptrarray[1])) s_ = ((dataptrarray[2]))[0] bm1_ = ((dataptrarray[3]))[0] ppr_ = ((dataptrarray[4])) eraAb(pnat_, v_, s_, bm1_, ppr_) - it.iternext() + status = iternext(GetNpyIter(it)) return ppr_out From e7d414d81f271d0d1dd1a09aca52974c8fb52f3b Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 4 Nov 2014 17:12:09 +0100 Subject: [PATCH 42/60] The 'multi_index' flag is not needed anymore --- nditer/erfa.pyx | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/nditer/erfa.pyx b/nditer/erfa.pyx index ab9a673..4c21fbb 100644 --- a/nditer/erfa.pyx +++ b/nditer/erfa.pyx @@ -41,7 +41,7 @@ def ab(pnat, v, s, bm1): #Create the iterator, broadcasting on all but the last dimension ndim = max([arr.ndim for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in]]) op_axes = [[-1]*(ndim-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] - it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes, flags=["multi_index"]) + it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes) cdef double* pnat_ cdef double* v_ From d12f79f260595592b3fb683cf11b303c0e359ac0 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 4 Nov 2014 20:19:10 +0100 Subject: [PATCH 43/60] mdim and shape from the broadcasted array --- nditer/erfa.pyx | 7 +++---- 1 file changed, 3 insertions(+), 4 deletions(-) diff --git a/nditer/erfa.pyx b/nditer/erfa.pyx index 4c21fbb..48722fd 100644 --- a/nditer/erfa.pyx +++ b/nditer/erfa.pyx @@ -35,12 +35,11 @@ def ab(pnat, v, s, bm1): bm1_in = np.array(bm1) #Create the output array, based on the broadcasted shape, adding the last dimension if needed - shape = np.broadcast(pnat_in[...,0], v_in[...,0], s_in, bm1_in).shape - ppr_out = np.empty(shape+(3,), dtype=np.double) + broadcast = np.broadcast(pnat_in[...,0], v_in[...,0], s_in, bm1_in) + ppr_out = np.empty(broadcast.shape+(3,), dtype=np.double) #Create the iterator, broadcasting on all but the last dimension - ndim = max([arr.ndim for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in]]) - op_axes = [[-1]*(ndim-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes) cdef double* pnat_ From bd8c2e8f9c6510351c8ab564d9e61c41b23109b2 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 4 Nov 2014 20:19:29 +0100 Subject: [PATCH 44/60] Comments tweak --- nditer/erfa.pyx | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/nditer/erfa.pyx b/nditer/erfa.pyx index 48722fd..06599de 100644 --- a/nditer/erfa.pyx +++ b/nditer/erfa.pyx @@ -34,11 +34,11 @@ def ab(pnat, v, s, bm1): s_in = np.array(s) bm1_in = np.array(bm1) - #Create the output array, based on the broadcasted shape, adding the last dimension if needed + #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed broadcast = np.broadcast(pnat_in[...,0], v_in[...,0], s_in, bm1_in) ppr_out = np.empty(broadcast.shape+(3,), dtype=np.double) - #Create the iterator, broadcasting on all but the last dimension + #Create the iterator, broadcasting on all but the consumed dimensions op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes) From 4649370f9f90255ef79b6d2c736e8ffd03444a2a Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Tue, 4 Nov 2014 20:37:10 +0100 Subject: [PATCH 45/60] Do not necessarily copy input arrays --- nditer/erfa.pyx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/nditer/erfa.pyx b/nditer/erfa.pyx index 06599de..e5d2d91 100644 --- a/nditer/erfa.pyx +++ b/nditer/erfa.pyx @@ -29,10 +29,10 @@ cdef inline NpyIter* GetNpyIter(object iter): def ab(pnat, v, s, bm1): #Turn all inputs into arrays - pnat_in = np.array(pnat) - v_in = np.array(v) - s_in = np.array(s) - bm1_in = np.array(bm1) + pnat_in = np.array(pnat, dtype=np.double, copy=False, subok=True) + v_in = np.array(v, dtype=np.double, copy=False, subok=True) + s_in = np.array(s, dtype=np.double, copy=False, subok=True) + bm1_in = np.array(bm1, dtype=np.double, copy=False, subok=True) #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed broadcast = np.broadcast(pnat_in[...,0], v_in[...,0], s_in, bm1_in) From bc5bfebdbb854fec8953114c0eebac5d386f9537 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 08:36:13 +0100 Subject: [PATCH 46/60] Renamed nditer to cython_npyiter --- {nditer => cython_npyiter}/erfa.pyx | 0 setup.py | 6 +++--- 2 files changed, 3 insertions(+), 3 deletions(-) rename {nditer => cython_npyiter}/erfa.pyx (100%) diff --git a/nditer/erfa.pyx b/cython_npyiter/erfa.pyx similarity index 100% rename from nditer/erfa.pyx rename to cython_npyiter/erfa.pyx diff --git a/setup.py b/setup.py index ec52c41..0d84852 100644 --- a/setup.py +++ b/setup.py @@ -16,8 +16,8 @@ include_dirs = [np.get_include(), '/opt/local/include'], library_dirs = ['/opt/local/lib']) -nditer_ext = Extension("erfa.nditer", - ["nditer/erfa.pyx"], +cython_npyiter_ext = Extension("erfa.cython_npyiter", + ["cython_npyiter/erfa.pyx"], libraries=['erfa'], include_dirs = [np.get_include(), '/opt/local/include'], library_dirs = ['/opt/local/lib']) @@ -26,4 +26,4 @@ cmdclass = { "build_ext": build_ext }, packages = ["erfa"], package_dir = {"erfa":"."}, - ext_modules = [cython_numpy_ext, cython_numpy_auto_ext, nditer_ext]) + ext_modules = [cython_numpy_ext, cython_numpy_auto_ext, cython_npyiter_ext]) From ba4465feaebb70f4f4f861f7ecbcd1e542eb9a6c Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 08:36:36 +0100 Subject: [PATCH 47/60] Ignore cython generated cython_npyiter/erfa.c --- cython_npyiter/.gitignore | 1 + 1 file changed, 1 insertion(+) create mode 100644 cython_npyiter/.gitignore diff --git a/cython_npyiter/.gitignore b/cython_npyiter/.gitignore new file mode 100644 index 0000000..a72b5f9 --- /dev/null +++ b/cython_npyiter/.gitignore @@ -0,0 +1 @@ +/erfa.c From c49cda53dcfcfe245cf8cfb9ce4754338c3153a8 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 10:35:27 +0100 Subject: [PATCH 48/60] In cython_numpy, C int are numpy.int32 --- cython_numpy/erfa.pyx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/cython_numpy/erfa.pyx b/cython_numpy/erfa.pyx index f3943e8..87a3493 100644 --- a/cython_numpy/erfa.pyx +++ b/cython_numpy/erfa.pyx @@ -135,10 +135,10 @@ def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, php def d2dtf(scale, ndp, d1, d2): shape = np.broadcast(scale, ndp, d1, d2).shape - iy_out = np.empty(shape, dtype=np.int) - im_out = np.empty(shape, dtype=np.int) - id_out = np.empty(shape, dtype=np.int) - ihmsf_out = np.empty(shape, dtype=[('h','i'),('m','i'),('s','i'),('f','i')]) + iy_out = np.empty(shape, dtype=np.int32) + im_out = np.empty(shape, dtype=np.int32) + id_out = np.empty(shape, dtype=np.int32) + ihmsf_out = np.empty(shape, dtype=[('h','i4'),('m','i4'),('s','i4'),('f','i4')]) cdef np.broadcast it = np.broadcast(scale, ndp, d1, d2, iy_out, im_out, id_out, ihmsf_out) From 42f5a92580a78e46e3cab9f54fb1d2a6d7cfcc0a Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 10:37:00 +0100 Subject: [PATCH 49/60] Implemented d2dtf() --- cython_npyiter/erfa.pyx | 48 +++++++++++++++++++++++++++++++++++++++++ 1 file changed, 48 insertions(+) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index e5d2d91..86be085 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -9,6 +9,7 @@ __all__ = ['ab'] cdef extern from "erfa.h": void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) + void eraAper(double theta, eraASTROM *astrom) ctypedef void NpyIter ctypedef int (*IterNextFunc)(NpyIter * iter) nogil @@ -25,6 +26,53 @@ cdef inline NpyIter* GetNpyIter(object iter): return (iter).iter +#int eraD2dtf(const char *scale, int ndp, double d1, double d2, int *iy, int *im, int *id, int ihmsf[4]) +def d2dtf(scale, ndp, d1, d2): + + #Turn all inputs into arrays + ndp_in = np.array(ndp, dtype=np.int32, copy=False, subok=True) + d1_in = np.array(d1, dtype=np.double, copy=False, subok=True) + d2_in = np.array(d2, dtype=np.double, copy=False, subok=True) + + #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed + broadcast = np.broadcast(ndp_in[...], d1_in[...], d2_in[...]) + iy_out = np.zeros(broadcast.shape, dtype=np.int32) + im_out = np.zeros(broadcast.shape, dtype=np.int32) + id_out = np.zeros(broadcast.shape, dtype=np.int32) + ihmsf_out = np.zeros(broadcast.shape+(4,), dtype=np.int32) + + #Create the iterator, broadcasting on all but the consumed dimensions + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [ndp_in[...], d1_in[...], d2_in[...], iy_out[...], im_out[...], id_out[...], ihmsf_out[...,0]]] + op_flags = [['readonly'], ['readonly'], ['readonly'], ['readwrite'], ['readwrite'], ['readwrite'], ['readwrite']] + it = np.nditer([ndp_in, d1_in, d2_in, iy_out, im_out, id_out, ihmsf_out], op_axes=op_axes, op_flags=op_flags) + + cdef char* scale_ = scale + cdef int ndp_ + cdef double d1_ + cdef double d2_ + cdef int* iy_ + cdef int* im_ + cdef int* id_ + cdef int* ihmsf_ + + #Iterate + cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) + cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) + cdef int status = 1 + while status: + ndp_ = (< int *>(dataptrarray[0]))[0] + d1_ = ((dataptrarray[1]))[0] + d2_ = ((dataptrarray[2]))[0] + iy_ = ( (dataptrarray[3])) + im_ = ( (dataptrarray[4])) + id_ = ( (dataptrarray[5])) + ihmsf_ = ( (dataptrarray[6])) + ret = eraD2dtf(scale_, ndp_, d1_, d2_, iy_, im_, id_, ihmsf_) + status = iternext(GetNpyIter(it)) + + return iy_out, im_out, id_out, ihmsf_out + + #void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) def ab(pnat, v, s, bm1): From 738cc64e4f453936e0540be49a88129d57c79d38 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 10:37:11 +0100 Subject: [PATCH 50/60] Implemented aper() --- cython_npyiter/erfa.pyx | 56 +++++++++++++++++++++++++++++++++++++++++ 1 file changed, 56 insertions(+) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index 86be085..837898b 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -9,8 +9,32 @@ __all__ = ['ab'] cdef extern from "erfa.h": void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) + int eraD2dtf(const char *scale, int ndp, double d1, double d2, + int *iy, int *im, int *id, int ihmsf[4]) void eraAper(double theta, eraASTROM *astrom) +cdef extern from "erfam.h": + cdef struct eraASTROM: + pass + +dt_eraASTROM = np.dtype([('pmt','d'), + ('eb','d',(3,)), + ('eh','d',(3,)), + ('em','d'), + ('v','d',(3,)), + ('bm1 ','d'), + ('bpn','d',(3,3)), + ('along','d'), + ('phi','d'), + ('xpl','d'), + ('ypl','d'), + ('sphi','d'), + ('cphi','d'), + ('diurab','d'), + ('eral','d'), + ('refa','d'), + ('refb','d')], align=True) + ctypedef void NpyIter ctypedef int (*IterNextFunc)(NpyIter * iter) nogil @@ -73,6 +97,38 @@ def d2dtf(scale, ndp, d1, d2): return iy_out, im_out, id_out, ihmsf_out +#void eraAper(double theta, eraASTROM *astrom) +def aper(theta, astrom): + + #Turn all inputs into arrays + theta_in = np.array(theta, dtype=np.double, copy=False, subok=True) + astrom_in = np.array(astrom, dtype=dt_eraASTROM, copy=False, subok=True) + + #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed + broadcast = np.broadcast(theta_in[...], astrom_in[...]) + astrom_out = np.empty(broadcast.shape, dtype=dt_eraASTROM) + np.copyto(astrom_out, astrom_in) + + #Create the iterator, broadcasting on all but the consumed dimensions + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [theta_in[...], astrom_out[...]]] + it = np.nditer([theta_in, astrom_out], op_axes=op_axes) + + cdef double theta_ + cdef eraASTROM* _astrom + + #Iterate + cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) + cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) + cdef int status = 1 + while status: + theta_ = ( (dataptrarray[0]))[0] + astrom_ = ((dataptrarray[1])) + eraAper(theta_, astrom_) + status = iternext(GetNpyIter(it)) + + return astrom_out + + #void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) def ab(pnat, v, s, bm1): From a7a1e509efb491f8f4f735c0259da1a5b7b20f45 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 10:57:47 +0100 Subject: [PATCH 51/60] Implement atco31() --- cython_npyiter/erfa.pyx | 110 ++++++++++++++++++++++++++++++++++++++++ 1 file changed, 110 insertions(+) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index 837898b..cc7ba42 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -8,6 +8,13 @@ __all__ = ['ab'] cdef extern from "erfa.h": + int eraAtco13(double rc, double dc, + double pr, double pd, double px, double rv, + double utc1, double utc2, double dut1, + double elong, double phi, double hm, double xp, double yp, + double phpa, double tk, double rh, double wl, + double *aob, double *zob, double *hob, + double *dob, double *rob, double *eo) void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) int eraD2dtf(const char *scale, int ndp, double d1, double d2, int *iy, int *im, int *id, int ihmsf[4]) @@ -50,6 +57,109 @@ cdef inline NpyIter* GetNpyIter(object iter): return (iter).iter +#int eraAtco13(double rc, double dc, +# double pr, double pd, double px, double rv, +# double utc1, double utc2, double dut1, +# double elong, double phi, double hm, double xp, double yp, +# double phpa, double tk, double rh, double wl, +# double *aob, double *zob, double *hob, +# double *dob, double *rob, double *eo) +def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl): + + #Turn all inputs into arrays + rc_in = np.array(rc, dtype=np.double, copy=False, subok=True) + dc_in = np.array(dc, dtype=np.double, copy=False, subok=True) + pr_in = np.array(pr, dtype=np.double, copy=False, subok=True) + pd_in = np.array(pd, dtype=np.double, copy=False, subok=True) + px_in = np.array(px, dtype=np.double, copy=False, subok=True) + rv_in = np.array(rv, dtype=np.double, copy=False, subok=True) + utc1_in = np.array(utc1, dtype=np.double, copy=False, subok=True) + utc2_in = np.array(utc2, dtype=np.double, copy=False, subok=True) + dut1_in = np.array(dut1, dtype=np.double, copy=False, subok=True) + elong_in = np.array(elong, dtype=np.double, copy=False, subok=True) + phi_in = np.array(phi, dtype=np.double, copy=False, subok=True) + hm_in = np.array(hm, dtype=np.double, copy=False, subok=True) + xp_in = np.array(xp, dtype=np.double, copy=False, subok=True) + yp_in = np.array(yp, dtype=np.double, copy=False, subok=True) + phpa_in = np.array(phpa, dtype=np.double, copy=False, subok=True) + tk_in = np.array(tk, dtype=np.double, copy=False, subok=True) + rh_in = np.array(rh, dtype=np.double, copy=False, subok=True) + wl_in = np.array(wl, dtype=np.double, copy=False, subok=True) + + #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed + broadcast = np.broadcast(rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tk_in, rh_in, wl_in) + aob_out = np.empty(broadcast.shape, dtype=np.double) + zob_out = np.empty(broadcast.shape, dtype=np.double) + hob_out = np.empty(broadcast.shape, dtype=np.double) + dob_out = np.empty(broadcast.shape, dtype=np.double) + rob_out = np.empty(broadcast.shape, dtype=np.double) + eo_out = np.empty(broadcast.shape, dtype=np.double) + + #Create the iterator, broadcasting on all but the consumed dimensions + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tk_in, rh_in, wl_in, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out]] + op_flags = [['readonly']]*18+[['readwrite']]*6 + it = np.nditer([rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tk_in, rh_in, wl_in, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out], op_axes=op_axes, op_flags=op_flags) + + cdef double _rc + cdef double _dc + cdef double _pr + cdef double _pd + cdef double _px + cdef double _rv + cdef double _utc1 + cdef double _utc2 + cdef double _dut1 + cdef double _elong + cdef double _phi + cdef double _hm + cdef double _xp + cdef double _yp + cdef double _phpa + cdef double _tk + cdef double _rh + cdef double _wl + cdef double *_aob + cdef double *_zob + cdef double *_hob + cdef double *_dob + cdef double *_rob + cdef double *_eo + + #Iterate + cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) + cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) + cdef int status = 1 + while status: + _rc = ((dataptrarray[ 0]))[0] + _dc = ((dataptrarray[ 1]))[0] + _pr = ((dataptrarray[ 2]))[0] + _pd = ((dataptrarray[ 3]))[0] + _px = ((dataptrarray[ 4]))[0] + _rv = ((dataptrarray[ 5]))[0] + _utc1 = ((dataptrarray[ 6]))[0] + _utc2 = ((dataptrarray[ 7]))[0] + _dut1 = ((dataptrarray[ 8]))[0] + _elong = ((dataptrarray[ 9]))[0] + _phi = ((dataptrarray[10]))[0] + _hm = ((dataptrarray[11]))[0] + _xp = ((dataptrarray[12]))[0] + _yp = ((dataptrarray[13]))[0] + _phpa = ((dataptrarray[14]))[0] + _tk = ((dataptrarray[15]))[0] + _rh = ((dataptrarray[16]))[0] + _wl = ((dataptrarray[17]))[0] + _aob = ((dataptrarray[18])) + _zob = ((dataptrarray[19])) + _hob = ((dataptrarray[20])) + _dob = ((dataptrarray[21])) + _rob = ((dataptrarray[22])) + _eo = ((dataptrarray[23])) + ret = eraAtco13(_rc, _dc, _pr, _pd, _px, _rv, _utc1, _utc2, _dut1, _elong, _phi, _hm, _xp, _yp, _phpa, _tk, _rh, _wl, _aob, _zob, _hob, _dob, _rob, _eo) + status = iternext(GetNpyIter(it)) + + return aob_out, zob_out, hob_out, dob_out, rob_out, eo_out + + #int eraD2dtf(const char *scale, int ndp, double d1, double d2, int *iy, int *im, int *id, int ihmsf[4]) def d2dtf(scale, ndp, d1, d2): From f236e950e46cdeae16f52d72caecd36b4e0eb345 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 10:58:06 +0100 Subject: [PATCH 52/60] np.zeros() -> np.empty() (faster) --- cython_npyiter/erfa.pyx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index cc7ba42..dda0edc 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -170,10 +170,10 @@ def d2dtf(scale, ndp, d1, d2): #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed broadcast = np.broadcast(ndp_in[...], d1_in[...], d2_in[...]) - iy_out = np.zeros(broadcast.shape, dtype=np.int32) - im_out = np.zeros(broadcast.shape, dtype=np.int32) - id_out = np.zeros(broadcast.shape, dtype=np.int32) - ihmsf_out = np.zeros(broadcast.shape+(4,), dtype=np.int32) + iy_out = np.empty(broadcast.shape, dtype=np.int32) + im_out = np.empty(broadcast.shape, dtype=np.int32) + id_out = np.empty(broadcast.shape, dtype=np.int32) + ihmsf_out = np.empty(broadcast.shape+(4,), dtype=np.int32) #Create the iterator, broadcasting on all but the consumed dimensions op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [ndp_in[...], d1_in[...], d2_in[...], iy_out[...], im_out[...], id_out[...], ihmsf_out[...,0]]] From d1e1a9f5c206733f05eb61a70cce95dc01899ba2 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 10:58:28 +0100 Subject: [PATCH 53/60] Also test cython_npyiter --- tests/test_erfa.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/tests/test_erfa.py b/tests/test_erfa.py index 7126c05..cd4b52a 100644 --- a/tests/test_erfa.py +++ b/tests/test_erfa.py @@ -1,7 +1,7 @@ import importlib import numpy as np -for erfa_wrapper_name in ['cython_numpy', 'cython_numpy_auto']: +for erfa_wrapper_name in ['cython_numpy', 'cython_numpy_auto', 'cython_npyiter']: erfa_module_name = '.'.join(['erfa', erfa_wrapper_name]) print("Testing with {0}...".format(erfa_module_name)) From 9be207a877aaf9dcde01fd5d164248407b68f6e3 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 10:59:31 +0100 Subject: [PATCH 54/60] Remove bare ellipsis --- cython_npyiter/erfa.pyx | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index dda0edc..c110d78 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -169,14 +169,14 @@ def d2dtf(scale, ndp, d1, d2): d2_in = np.array(d2, dtype=np.double, copy=False, subok=True) #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed - broadcast = np.broadcast(ndp_in[...], d1_in[...], d2_in[...]) + broadcast = np.broadcast(ndp_in, d1_in, d2_in) iy_out = np.empty(broadcast.shape, dtype=np.int32) im_out = np.empty(broadcast.shape, dtype=np.int32) id_out = np.empty(broadcast.shape, dtype=np.int32) ihmsf_out = np.empty(broadcast.shape+(4,), dtype=np.int32) #Create the iterator, broadcasting on all but the consumed dimensions - op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [ndp_in[...], d1_in[...], d2_in[...], iy_out[...], im_out[...], id_out[...], ihmsf_out[...,0]]] + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [ndp_in, d1_in, d2_in, iy_out, im_out, id_out, ihmsf_out[...,0]]] op_flags = [['readonly'], ['readonly'], ['readonly'], ['readwrite'], ['readwrite'], ['readwrite'], ['readwrite']] it = np.nditer([ndp_in, d1_in, d2_in, iy_out, im_out, id_out, ihmsf_out], op_axes=op_axes, op_flags=op_flags) @@ -215,12 +215,12 @@ def aper(theta, astrom): astrom_in = np.array(astrom, dtype=dt_eraASTROM, copy=False, subok=True) #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed - broadcast = np.broadcast(theta_in[...], astrom_in[...]) + broadcast = np.broadcast(theta_in, astrom_in) astrom_out = np.empty(broadcast.shape, dtype=dt_eraASTROM) np.copyto(astrom_out, astrom_in) #Create the iterator, broadcasting on all but the consumed dimensions - op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [theta_in[...], astrom_out[...]]] + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [theta_in, astrom_out]] it = np.nditer([theta_in, astrom_out], op_axes=op_axes) cdef double theta_ From 5fc60f3dbddf6bc08953fea8e3cfa51f1b12d340 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 11:02:39 +0100 Subject: [PATCH 55/60] op_flags generalisation and cleanup --- cython_npyiter/erfa.pyx | 8 +++++--- 1 file changed, 5 insertions(+), 3 deletions(-) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index c110d78..1428544 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -177,7 +177,7 @@ def d2dtf(scale, ndp, d1, d2): #Create the iterator, broadcasting on all but the consumed dimensions op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [ndp_in, d1_in, d2_in, iy_out, im_out, id_out, ihmsf_out[...,0]]] - op_flags = [['readonly'], ['readonly'], ['readonly'], ['readwrite'], ['readwrite'], ['readwrite'], ['readwrite']] + op_flags = [['readonly']]*3+[['readwrite']]*4 it = np.nditer([ndp_in, d1_in, d2_in, iy_out, im_out, id_out, ihmsf_out], op_axes=op_axes, op_flags=op_flags) cdef char* scale_ = scale @@ -221,7 +221,8 @@ def aper(theta, astrom): #Create the iterator, broadcasting on all but the consumed dimensions op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [theta_in, astrom_out]] - it = np.nditer([theta_in, astrom_out], op_axes=op_axes) + op_flags = [['readonly']]*1+[['readwrite']]*1 + it = np.nditer([theta_in, astrom_out], op_axes=op_axes, op_flags=op_flags) cdef double theta_ cdef eraASTROM* _astrom @@ -254,7 +255,8 @@ def ab(pnat, v, s, bm1): #Create the iterator, broadcasting on all but the consumed dimensions op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] - it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes) + op_flags = [['readonly']]*4+[['readwrite']]*1 + it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes, op_flags=op_flags) cdef double* pnat_ cdef double* v_ From d86de31d5d736a58ba788d6aefc7c7349f9c3ec1 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 11:10:16 +0100 Subject: [PATCH 56/60] Cleanup of C variable naming --- cython_npyiter/erfa.pyx | 62 ++++++++++++++++++++--------------------- 1 file changed, 31 insertions(+), 31 deletions(-) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index 1428544..1b7d406 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -180,28 +180,28 @@ def d2dtf(scale, ndp, d1, d2): op_flags = [['readonly']]*3+[['readwrite']]*4 it = np.nditer([ndp_in, d1_in, d2_in, iy_out, im_out, id_out, ihmsf_out], op_axes=op_axes, op_flags=op_flags) - cdef char* scale_ = scale - cdef int ndp_ - cdef double d1_ - cdef double d2_ - cdef int* iy_ - cdef int* im_ - cdef int* id_ - cdef int* ihmsf_ + cdef char* _scale = scale + cdef int _ndp + cdef double _d1 + cdef double _d2 + cdef int* _iy + cdef int* _im + cdef int* _id + cdef int* _ihmsf #Iterate cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 while status: - ndp_ = (< int *>(dataptrarray[0]))[0] - d1_ = ((dataptrarray[1]))[0] - d2_ = ((dataptrarray[2]))[0] - iy_ = ( (dataptrarray[3])) - im_ = ( (dataptrarray[4])) - id_ = ( (dataptrarray[5])) - ihmsf_ = ( (dataptrarray[6])) - ret = eraD2dtf(scale_, ndp_, d1_, d2_, iy_, im_, id_, ihmsf_) + _ndp = (< int *>(dataptrarray[0]))[0] + _d1 = ((dataptrarray[1]))[0] + _d2 = ((dataptrarray[2]))[0] + _iy = ( (dataptrarray[3])) + _im = ( (dataptrarray[4])) + _id = ( (dataptrarray[5])) + _ihmsf = ( (dataptrarray[6])) + ret = eraD2dtf(_scale, _ndp, _d1, _d2, _iy, _im, _id, _ihmsf) status = iternext(GetNpyIter(it)) return iy_out, im_out, id_out, ihmsf_out @@ -224,7 +224,7 @@ def aper(theta, astrom): op_flags = [['readonly']]*1+[['readwrite']]*1 it = np.nditer([theta_in, astrom_out], op_axes=op_axes, op_flags=op_flags) - cdef double theta_ + cdef double _theta cdef eraASTROM* _astrom #Iterate @@ -232,9 +232,9 @@ def aper(theta, astrom): cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 while status: - theta_ = ( (dataptrarray[0]))[0] - astrom_ = ((dataptrarray[1])) - eraAper(theta_, astrom_) + _theta = ( (dataptrarray[0]))[0] + _astrom = ((dataptrarray[1])) + eraAper(_theta, _astrom) status = iternext(GetNpyIter(it)) return astrom_out @@ -258,23 +258,23 @@ def ab(pnat, v, s, bm1): op_flags = [['readonly']]*4+[['readwrite']]*1 it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes, op_flags=op_flags) - cdef double* pnat_ - cdef double* v_ - cdef double s_ - cdef double bm1_ - cdef double* ppr_ + cdef double* _pnat + cdef double* _v + cdef double _s + cdef double _bm1 + cdef double* _ppr #Iterate cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 while status: - pnat_ = ((dataptrarray[0])) - v_ = ((dataptrarray[1])) - s_ = ((dataptrarray[2]))[0] - bm1_ = ((dataptrarray[3]))[0] - ppr_ = ((dataptrarray[4])) - eraAb(pnat_, v_, s_, bm1_, ppr_) + _pnat = ((dataptrarray[0])) + _v = ((dataptrarray[1])) + _s = ((dataptrarray[2]))[0] + _bm1 = ((dataptrarray[3]))[0] + _ppr = ((dataptrarray[4])) + eraAb(_pnat, _v, _s, _bm1, _ppr) status = iternext(GetNpyIter(it)) return ppr_out From 24e472117a9fb4f36dadedca3e92e2d6bdee4c81 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Wed, 5 Nov 2014 12:10:50 +0100 Subject: [PATCH 57/60] Implemented rx() --- cython_npyiter/erfa.pyx | 34 ++++++++++++++++++++++++++++++++++ 1 file changed, 34 insertions(+) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index 1b7d406..b538527 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -8,6 +8,7 @@ __all__ = ['ab'] cdef extern from "erfa.h": + void eraRx(double phi, double r[9]) int eraAtco13(double rc, double dc, double pr, double pd, double px, double rv, double utc1, double utc2, double dut1, @@ -57,6 +58,39 @@ cdef inline NpyIter* GetNpyIter(object iter): return (iter).iter +#void eraRx(double phi, double r[3][3]) +def rx(phi, r): + + #Turn all inputs into arrays + phi_in = np.array(phi, dtype=np.double, copy=False, subok=True) + r_in = np.array(r, dtype=np.double, copy=False, subok=True) + + #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed + broadcast = np.broadcast(phi_in, r_in[...,0,0]) + r_out = np.empty(broadcast.shape+(3,3), dtype=np.double) + np.copyto(r_out, r_in) + + #Create the iterator, broadcasting on all but the consumed dimensions + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [phi_in, r_out[...,0,0]]] + op_flags = [['readonly']]*1+[['readwrite']]*1 + it = np.nditer([phi_in, r_out], op_axes=op_axes, op_flags=op_flags) + + cdef double _phi + cdef double* _r + + #Iterate + cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) + cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) + cdef int status = 1 + while status: + _phi = ((dataptrarray[0]))[0] + _r = ((dataptrarray[1])) + eraRx(_phi, _r) + status = iternext(GetNpyIter(it)) + + return r_out + + #int eraAtco13(double rc, double dc, # double pr, double pd, double px, double rv, # double utc1, double utc2, double dut1, From d35e5cc317e3948edbded58414fd92f5d4d15fc6 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Fri, 7 Nov 2014 11:11:45 +0100 Subject: [PATCH 58/60] Initial commit of cython_npyiter_auto This is based on the work in cython_numpy_auto and the non-automated tests in cython_npyiter --- .gitignore | 4 + cython_npyiter_auto/code_analyzer.py | 132 ++++++++++++++++++++++++ cython_npyiter_auto/config.py | 26 +++++ cython_npyiter_auto/cython_generator.py | 54 ++++++++++ cython_npyiter_auto/erfa.pyx.jinja2 | 90 ++++++++++++++++ 5 files changed, 306 insertions(+) create mode 100644 cython_npyiter_auto/code_analyzer.py create mode 100644 cython_npyiter_auto/config.py create mode 100644 cython_npyiter_auto/cython_generator.py create mode 100644 cython_npyiter_auto/erfa.pyx.jinja2 diff --git a/.gitignore b/.gitignore index 157f002..0041ebd 100644 --- a/.gitignore +++ b/.gitignore @@ -8,3 +8,7 @@ cython_numpy_auto/erfa.c cython_numpy_auto/erfa.pyx + +cython_npyiter_auto/erfa.c + +cython_npyiter_auto/erfa.pyx diff --git a/cython_npyiter_auto/code_analyzer.py b/cython_npyiter_auto/code_analyzer.py new file mode 100644 index 0000000..a53f65e --- /dev/null +++ b/cython_npyiter_auto/code_analyzer.py @@ -0,0 +1,132 @@ +import os.path +import re + + +__all__ = ['Function'] + + +ctype_to_dtype = {'double' : "np.double", + 'double *' : "np.double", + 'int' : "np.int32", + 'int *' : "np.int32", + 'eraASTROM' : "dt_eraASTROM", + 'char' : "np.dtype('S1')", + } + + +class Argument(object): + + def __init__(self, definition, inout_state): + self.inout_state = inout_state + self.definition = definition.strip() + if self.definition[:5] == "const": + self.definition = self.definition[6:].strip() + if "*" in self.definition: + self.ctype, self.name = self.definition.split("*", 1) + self.ptr = "*" + self.shape = () + elif "[]" in self.definition: + self.ctype, self.name = self.definition.split("[]", 1) + self.ptr = "*" + self.shape = () + else: + self.ctype, self.name = self.definition.rsplit(" ", 1) + if "[" in self.name: + self.name, arr = self.name.split("[", 1) + self.ptr = "*" + self.shape = tuple([int(s) for s in arr[:-1].split("][")]) + else: + self.ptr = "" + self.shape = () + self.ctype = self.ctype.strip() + + @property + def name_in(self): + return self.name+"_in" + + @property + def name_out(self): + return self.name+"_out" + + @property + def name_in_broadcast(self): + if len(self.shape)>0: + return self.name_in+"[..."+",0"*len(self.shape)+"]" + else: + return self.name_in + + @property + def name_out_broadcast(self): + if len(self.shape)>0: + return self.name_out+"[..."+",0"*len(self.shape)+"]" + else: + return self.name_out + + @property + def ctype_ptr(self): + return self.ctype+self.ptr + + @property + def dtype(self): + return ctype_to_dtype[self.ctype] + + def __repr__(self): + return "Argument('{0}', name='{1}', ctype='{2}', inout_state='{3}')".format(self.definition, self.name, self.ctype, self.inout_state) + + +class Return(object): + + def __init__(self, ctype): + self.name = 'ret' + self.ctype = ctype + self.inout_state = '=' + self.ctype_ptr = ctype + + @property + def dtype(self): + return ctype_to_dtype[self.ctype] + + +class Function(object): + """ + The following attributes are derived from the supplied ERFA function + name: C name of the function + pyname: Python name of the function + signature: Signature of the function, as extracted from the source file + documentation: Documentation block for the function, as extracted from the source file + args: List of arguments of function + ret: Return parameter of function + """ + + def __init__(self, name, inout_state, source_path): + self.name = name + self.pyname = name.split('era')[-1].lower() + #From source get: signature and documentation + pattern = "\n([^\n]+{0} ?\([^)]+\)).+?(/\*.+?\*/)".format(name) + p = re.compile(pattern, flags=re.DOTALL|re.MULTILINE) + filepath = os.path.join(os.path.normpath(source_path), self.pyname+".c") + with open(filepath) as f: + search = p.search(f.read()) + self.signature = " ".join(search.group(1).split()) #Single line with single spaces + self.documentation = search.group(2).replace("**"," ").replace("/*"," ").replace("*/"," ") + #Extract arguments from signature + self.args = [] + arg_list = self.signature.split(')')[0].split('(')[1].split(', ') + for arg_index in range(len(arg_list)): + self.args.append(Argument(arg_list[arg_index], inout_state[arg_index])) + self.ret = re.search("^(.*){0}".format(name), self.signature).group(1).strip() + if self.ret == 'double': + self.args.append(Return(self.ret)) + + def args_by_inout(self, inout_filter, prop=None, join=None): + result = [] + for arg in self.args: + if arg.inout_state in inout_filter.split('|'): + if prop is None: + result.append(arg) + else: + result.append(getattr(arg, prop)) + if join is not None: + return join.join(result) + else: + return result diff --git a/cython_npyiter_auto/config.py b/cython_npyiter_auto/config.py new file mode 100644 index 0000000..0bccb8c --- /dev/null +++ b/cython_npyiter_auto/config.py @@ -0,0 +1,26 @@ +#Convention is the following: +# i: input +# o: output +# m: modify +# r: return + +wrappers = {} +#int eraAtco13(double rc, double dc, +# double pr, double pd, double px, double rv, +# double utc1, double utc2, double dut1, +# double elong, double phi, double hm, double xp, double yp, +# double phpa, double tk, double rh, double wl, +# double *aob, double *zob, double *hob, +# double *dob, double *rob, double *eo) +wrappers["eraAtco13"] = "iiiiiiiiiiiiiiiiiioooooo" +#void eraRx(double phi, double r[3][3]) +wrappers["eraRx"] = "im" +#int eraD2dtf(const char *scale, int ndp, double d1, double d2, int *iy, int *im, int *id, int ihmsf[4]) +wrappers["eraD2dtf"] = "iiiioooo" +#void eraAper(double theta, eraASTROM *astrom) +wrappers["eraAper"] = "im" +#void eraAb(double pnat[3], double v[3], double s, double bm1, +# double ppr[3]) +wrappers["eraAb"] = "iiiio" +#void eraIr(double r[3][3]) +wrappers["eraIr"] = "o" \ No newline at end of file diff --git a/cython_npyiter_auto/cython_generator.py b/cython_npyiter_auto/cython_generator.py new file mode 100644 index 0000000..cf1d51e --- /dev/null +++ b/cython_npyiter_auto/cython_generator.py @@ -0,0 +1,54 @@ +import re +from jinja2 import Template, Environment, PackageLoader +from code_analyzer import Function + +import config + +ERFA_SOURCES = "../../erfa/src" + + +#Prepare the jinja2 templating environment + +def prefix(a_list, pre): + return [pre+'{0}'.format(an_element) for an_element in a_list] + +def postfix(a_list, post): + return ['{0}'.format(an_element)+post for an_element in a_list] + +def surround(a_list, pre, post): + return [pre+'{0}'.format(an_element)+post for an_element in a_list] + +env = Environment() +env.filters['prefix'] = prefix +env.filters['postfix'] = postfix +env.filters['surround'] = surround + + +#Extract all the ERFA function names from erfa.h + +with open(ERFA_SOURCES+"/erfa.h", "r") as f: + erfa_h = f.read() +print("Parse ERFA code:") +section_subsection_functions = re.findall('/\* (\w*)/(\w*) \*/\n(.*?)\n\n', erfa_h, flags=re.DOTALL|re.MULTILINE) +funcs = [] +for section, subsection, functions in section_subsection_functions: + print(" {0}.{1}".format(section, subsection)) + func_names = re.findall(' (\w+)\(.*?\);', functions, flags=re.DOTALL) + for name in func_names: + if name in config.wrappers.keys(): + function = Function(name, config.wrappers[name], ERFA_SOURCES) + funcs.append(function) + print(" {0}.{1}.{2} >{3}<...".format(section, subsection, name, function.signature)) +print("Done!") + + +#Render the template + +with open('erfa.pyx.jinja2', 'r') as template_file: + erfa_pyx = env.from_string(template_file.read()).render(funcs=funcs) + + +#Save the result + +with open("erfa.pyx", "w") as f: + f.write(erfa_pyx) diff --git a/cython_npyiter_auto/erfa.pyx.jinja2 b/cython_npyiter_auto/erfa.pyx.jinja2 new file mode 100644 index 0000000..3029490 --- /dev/null +++ b/cython_npyiter_auto/erfa.pyx.jinja2 @@ -0,0 +1,90 @@ +import numpy as np +cimport numpy as np +from cpython.ref cimport PyObject + +np.import_array() + +__all__ = [{{ funcs|map(attribute='pyname')|surround("'","'")|join(", ") }}] + + +cdef extern from "erfam.h": + cdef struct eraASTROM: + pass + +cdef extern from "erfa.h": +{%- for func in funcs %} + {{ func.ret }} {{ func.name }}({{ func.args_by_inout('i|m|o')|map(attribute='ctype_ptr')|join(', ') }}) +{%- endfor %} + +dt_eraASTROM = np.dtype([('pmt','d'), + ('eb','d',(3,)), + ('eh','d',(3,)), + ('em','d'), + ('v','d',(3,)), + ('bm1 ','d'), + ('bpn','d',(3,3)), + ('along','d'), + ('phi','d'), + ('xpl','d'), + ('ypl','d'), + ('sphi','d'), + ('cphi','d'), + ('diurab','d'), + ('eral','d'), + ('refa','d'), + ('refb','d')], align=True) + +ctypedef void NpyIter +ctypedef int (*IterNextFunc)(NpyIter * iter) nogil + +cdef extern from "numpy/arrayobject.h": + IterNextFunc GetIterNext "NpyIter_GetIterNext" (NpyIter *iter, char **) + char** GetDataPtrArray "NpyIter_GetDataPtrArray" (NpyIter* iter) + +ctypedef struct NewNpyArrayIterObject: + PyObject base + NpyIter *iter + +cdef inline NpyIter* GetNpyIter(object iter): + return (iter).iter + + +{% for func in funcs %} +def {{ func.pyname }}({{ func.args_by_inout('i|m')|map(attribute='name')|join(', ') }}): + + #Turn all inputs into arrays + {%- for arg in func.args_by_inout('i|m') %} + {{ arg.name_in }} = np.array({{ arg.name }}, dtype={{ arg.dtype }}, copy=False, subok=True) + {%- endfor %} + + #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed + broadcast = np.broadcast(np.int32(0.0), np.int32(0.0), {{ func.args_by_inout('i|m')|map(attribute='name_in_broadcast')|join(', ') }}) + {%- for arg in func.args_by_inout('m|o|r') %} + {{ arg.name_out }} = np.empty(broadcast.shape+{{ arg.shape }}, dtype={{ arg.dtype }}) + {%- endfor %} + {%- for arg in func.args_by_inout('m') %} + np.copyto({{ arg.name_out }}, {{ arg.name_in }}) + {%- endfor %} + + #Create the iterator, broadcasting on all but the consumed dimensions + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [{{ func.args_by_inout('i')|map(attribute='name_in_broadcast')|join(', ') }}]+[{{ func.args_by_inout('m|o')|map(attribute='name_out_broadcast')|join(', ') }}]] + op_flags = [['readonly']]*{{ func.args_by_inout('i')|count }}+[['readwrite']]*{{ func.args_by_inout('m|o')|count }} + it = np.nditer([{{ func.args_by_inout('i')|map(attribute='name_in')|join(', ') }}]+[{{ func.args_by_inout('m|o')|map(attribute='name_out')|join(', ') }}], op_axes=op_axes, op_flags=op_flags) + + #Iterate + {%- for arg in func.args_by_inout('i|m|o') %} + cdef {{ arg.ctype_ptr }} _{{ arg.name }} + {%- endfor %} + cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) + cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) + cdef int status = 1 + while status: + {%- for arg in func.args_by_inout('i|m|o') %} + _{{ arg.name }} = (<{{ arg.ctype }} *>(dataptrarray[{{ func.args.index(arg) }}])){%- if arg.ctype_ptr[-1] != '*' %}[0]{% endif %} + {%- endfor %} + {{ func.args_by_inout('ret')|map(attribute='name')|surround('_',' = ')|join }}{{ func.name }}({{ func.args_by_inout('i|m|o')|map(attribute='name')|prefix('_')|join(', ') }}) + status = iternext(GetNpyIter(it)) + + return {{ func.args_by_inout('m|o|r')|map(attribute='name_out')|join(', ') }} + +{% endfor %} From f40522c5daaf3a80f63cbd8f45e801b937ae15f3 Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Fri, 7 Nov 2014 11:12:49 +0100 Subject: [PATCH 59/60] Minor cleanups in cython_npyiter To make it match with the automated version and help me splot differences in generated erfa.pyx --- cython_npyiter/erfa.pyx | 297 +++++++++++++++++++++------------------- 1 file changed, 154 insertions(+), 143 deletions(-) diff --git a/cython_npyiter/erfa.pyx b/cython_npyiter/erfa.pyx index b538527..55606ca 100644 --- a/cython_npyiter/erfa.pyx +++ b/cython_npyiter/erfa.pyx @@ -4,26 +4,27 @@ from cpython.ref cimport PyObject np.import_array() -__all__ = ['ab'] +__all__ = ['ab', 'aper', 'atco13', 'd2dtf', 'rx'] +cdef extern from "erfam.h": + cdef struct eraASTROM: + pass + cdef extern from "erfa.h": - void eraRx(double phi, double r[9]) + void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) + void eraAper(double theta, eraASTROM *astrom) int eraAtco13(double rc, double dc, double pr, double pd, double px, double rv, double utc1, double utc2, double dut1, double elong, double phi, double hm, double xp, double yp, - double phpa, double tk, double rh, double wl, + double phpa, double tc, double rh, double wl, double *aob, double *zob, double *hob, double *dob, double *rob, double *eo) - void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) int eraD2dtf(const char *scale, int ndp, double d1, double d2, int *iy, int *im, int *id, int ihmsf[4]) - void eraAper(double theta, eraASTROM *astrom) - -cdef extern from "erfam.h": - cdef struct eraASTROM: - pass + void eraRx(double phi, double r[9]) + void eraIr(double*) dt_eraASTROM = np.dtype([('pmt','d'), ('eb','d',(3,)), @@ -58,82 +59,112 @@ cdef inline NpyIter* GetNpyIter(object iter): return (iter).iter -#void eraRx(double phi, double r[3][3]) -def rx(phi, r): +def ab(pnat, v, s, bm1): #Turn all inputs into arrays - phi_in = np.array(phi, dtype=np.double, copy=False, subok=True) - r_in = np.array(r, dtype=np.double, copy=False, subok=True) + pnat_in = np.array(pnat, dtype=np.double, copy=False, subok=True) + v_in = np.array(v, dtype=np.double, copy=False, subok=True) + s_in = np.array(s, dtype=np.double, copy=False, subok=True) + bm1_in = np.array(bm1, dtype=np.double, copy=False, subok=True) #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed - broadcast = np.broadcast(phi_in, r_in[...,0,0]) - r_out = np.empty(broadcast.shape+(3,3), dtype=np.double) - np.copyto(r_out, r_in) + broadcast = np.broadcast(pnat_in[...,0], v_in[...,0], s_in, bm1_in) + ppr_out = np.empty(broadcast.shape+(3,), dtype=np.double) #Create the iterator, broadcasting on all but the consumed dimensions - op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [phi_in, r_out[...,0,0]]] - op_flags = [['readonly']]*1+[['readwrite']]*1 - it = np.nditer([phi_in, r_out], op_axes=op_axes, op_flags=op_flags) - - cdef double _phi - cdef double* _r + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] + op_flags = [['readonly']]*4+[['readwrite']]*1 + it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes, op_flags=op_flags) #Iterate + cdef double* _pnat + cdef double* _v + cdef double _s + cdef double _bm1 + cdef double* _ppr cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 while status: - _phi = ((dataptrarray[0]))[0] - _r = ((dataptrarray[1])) - eraRx(_phi, _r) + _pnat = ((dataptrarray[0])) + _v = ((dataptrarray[1])) + _s = ((dataptrarray[2]))[0] + _bm1 = ((dataptrarray[3]))[0] + _ppr = ((dataptrarray[4])) + eraAb(_pnat, _v, _s, _bm1, _ppr) status = iternext(GetNpyIter(it)) - return r_out + return ppr_out + + +def aper(theta, astrom): + #Turn all inputs into arrays + theta_in = np.array(theta, dtype=np.double, copy=False, subok=True) + astrom_in = np.array(astrom, dtype=dt_eraASTROM, copy=False, subok=True) + + #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed + broadcast = np.broadcast(theta_in, astrom_in) + astrom_out = np.empty(broadcast.shape, dtype=dt_eraASTROM) + np.copyto(astrom_out, astrom_in) + + #Create the iterator, broadcasting on all but the consumed dimensions + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [theta_in, astrom_out]] + op_flags = [['readonly']]*1+[['readwrite']]*1 + it = np.nditer([theta_in, astrom_out], op_axes=op_axes, op_flags=op_flags) + + #Iterate + cdef double _theta + cdef eraASTROM* _astrom + cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) + cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) + cdef int status = 1 + while status: + _theta = ((dataptrarray[0]))[0] + _astrom = ((dataptrarray[1])) + eraAper(_theta, _astrom) + status = iternext(GetNpyIter(it)) + + return astrom_out + -#int eraAtco13(double rc, double dc, -# double pr, double pd, double px, double rv, -# double utc1, double utc2, double dut1, -# double elong, double phi, double hm, double xp, double yp, -# double phpa, double tk, double rh, double wl, -# double *aob, double *zob, double *hob, -# double *dob, double *rob, double *eo) -def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tk, rh, wl): +def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, phpa, tc, rh, wl): #Turn all inputs into arrays - rc_in = np.array(rc, dtype=np.double, copy=False, subok=True) - dc_in = np.array(dc, dtype=np.double, copy=False, subok=True) - pr_in = np.array(pr, dtype=np.double, copy=False, subok=True) - pd_in = np.array(pd, dtype=np.double, copy=False, subok=True) - px_in = np.array(px, dtype=np.double, copy=False, subok=True) - rv_in = np.array(rv, dtype=np.double, copy=False, subok=True) - utc1_in = np.array(utc1, dtype=np.double, copy=False, subok=True) - utc2_in = np.array(utc2, dtype=np.double, copy=False, subok=True) - dut1_in = np.array(dut1, dtype=np.double, copy=False, subok=True) + rc_in = np.array(rc, dtype=np.double, copy=False, subok=True) + dc_in = np.array(dc, dtype=np.double, copy=False, subok=True) + pr_in = np.array(pr, dtype=np.double, copy=False, subok=True) + pd_in = np.array(pd, dtype=np.double, copy=False, subok=True) + px_in = np.array(px, dtype=np.double, copy=False, subok=True) + rv_in = np.array(rv, dtype=np.double, copy=False, subok=True) + utc1_in = np.array(utc1, dtype=np.double, copy=False, subok=True) + utc2_in = np.array(utc2, dtype=np.double, copy=False, subok=True) + dut1_in = np.array(dut1, dtype=np.double, copy=False, subok=True) elong_in = np.array(elong, dtype=np.double, copy=False, subok=True) - phi_in = np.array(phi, dtype=np.double, copy=False, subok=True) - hm_in = np.array(hm, dtype=np.double, copy=False, subok=True) - xp_in = np.array(xp, dtype=np.double, copy=False, subok=True) - yp_in = np.array(yp, dtype=np.double, copy=False, subok=True) - phpa_in = np.array(phpa, dtype=np.double, copy=False, subok=True) - tk_in = np.array(tk, dtype=np.double, copy=False, subok=True) - rh_in = np.array(rh, dtype=np.double, copy=False, subok=True) - wl_in = np.array(wl, dtype=np.double, copy=False, subok=True) + phi_in = np.array(phi, dtype=np.double, copy=False, subok=True) + hm_in = np.array(hm, dtype=np.double, copy=False, subok=True) + xp_in = np.array(xp, dtype=np.double, copy=False, subok=True) + yp_in = np.array(yp, dtype=np.double, copy=False, subok=True) + phpa_in = np.array(phpa, dtype=np.double, copy=False, subok=True) + tc_in = np.array(tc, dtype=np.double, copy=False, subok=True) + rh_in = np.array(rh, dtype=np.double, copy=False, subok=True) + wl_in = np.array(wl, dtype=np.double, copy=False, subok=True) #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed - broadcast = np.broadcast(rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tk_in, rh_in, wl_in) + broadcast = np.broadcast(rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tc_in, rh_in, wl_in) aob_out = np.empty(broadcast.shape, dtype=np.double) zob_out = np.empty(broadcast.shape, dtype=np.double) hob_out = np.empty(broadcast.shape, dtype=np.double) dob_out = np.empty(broadcast.shape, dtype=np.double) rob_out = np.empty(broadcast.shape, dtype=np.double) - eo_out = np.empty(broadcast.shape, dtype=np.double) + eo_out = np.empty(broadcast.shape, dtype=np.double) #Create the iterator, broadcasting on all but the consumed dimensions - op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tk_in, rh_in, wl_in, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out]] + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tc_in, rh_in, wl_in, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out]] op_flags = [['readonly']]*18+[['readwrite']]*6 - it = np.nditer([rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tk_in, rh_in, wl_in, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out], op_axes=op_axes, op_flags=op_flags) + it = np.nditer([rc_in, dc_in, pr_in, pd_in, px_in, rv_in, utc1_in, utc2_in, dut1_in, elong_in, phi_in, hm_in, xp_in, yp_in, phpa_in, tc_in, rh_in, wl_in, aob_out, zob_out, hob_out, dob_out, rob_out, eo_out], op_axes=op_axes, op_flags=op_flags) + #Iterate cdef double _rc cdef double _dc cdef double _pr @@ -149,52 +180,49 @@ def atco13(rc, dc, pr, pd, px, rv, utc1, utc2, dut1, elong, phi, hm, xp, yp, php cdef double _xp cdef double _yp cdef double _phpa - cdef double _tk + cdef double _tc cdef double _rh cdef double _wl - cdef double *_aob - cdef double *_zob - cdef double *_hob - cdef double *_dob - cdef double *_rob - cdef double *_eo - - #Iterate + cdef double* _aob + cdef double* _zob + cdef double* _hob + cdef double* _dob + cdef double* _rob + cdef double* _eo cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 - while status: - _rc = ((dataptrarray[ 0]))[0] - _dc = ((dataptrarray[ 1]))[0] - _pr = ((dataptrarray[ 2]))[0] - _pd = ((dataptrarray[ 3]))[0] - _px = ((dataptrarray[ 4]))[0] - _rv = ((dataptrarray[ 5]))[0] - _utc1 = ((dataptrarray[ 6]))[0] - _utc2 = ((dataptrarray[ 7]))[0] - _dut1 = ((dataptrarray[ 8]))[0] - _elong = ((dataptrarray[ 9]))[0] - _phi = ((dataptrarray[10]))[0] - _hm = ((dataptrarray[11]))[0] - _xp = ((dataptrarray[12]))[0] - _yp = ((dataptrarray[13]))[0] - _phpa = ((dataptrarray[14]))[0] - _tk = ((dataptrarray[15]))[0] - _rh = ((dataptrarray[16]))[0] - _wl = ((dataptrarray[17]))[0] - _aob = ((dataptrarray[18])) - _zob = ((dataptrarray[19])) - _hob = ((dataptrarray[20])) - _dob = ((dataptrarray[21])) - _rob = ((dataptrarray[22])) - _eo = ((dataptrarray[23])) - ret = eraAtco13(_rc, _dc, _pr, _pd, _px, _rv, _utc1, _utc2, _dut1, _elong, _phi, _hm, _xp, _yp, _phpa, _tk, _rh, _wl, _aob, _zob, _hob, _dob, _rob, _eo) + while status: + _rc = ((dataptrarray[0]))[0] + _dc = ((dataptrarray[1]))[0] + _pr = ((dataptrarray[2]))[0] + _pd = ((dataptrarray[3]))[0] + _px = ((dataptrarray[4]))[0] + _rv = ((dataptrarray[5]))[0] + _utc1 = ((dataptrarray[6]))[0] + _utc2 = ((dataptrarray[7]))[0] + _dut1 = ((dataptrarray[8]))[0] + _elong = ((dataptrarray[9]))[0] + _phi = ((dataptrarray[10]))[0] + _hm = ((dataptrarray[11]))[0] + _xp = ((dataptrarray[12]))[0] + _yp = ((dataptrarray[13]))[0] + _phpa = ((dataptrarray[14]))[0] + _tc = ((dataptrarray[15]))[0] + _rh = ((dataptrarray[16]))[0] + _wl = ((dataptrarray[17]))[0] + _aob = ((dataptrarray[18])) + _zob = ((dataptrarray[19])) + _hob = ((dataptrarray[20])) + _dob = ((dataptrarray[21])) + _rob = ((dataptrarray[22])) + _eo = ((dataptrarray[23])) + ret = eraAtco13(_rc, _dc, _pr, _pd, _px, _rv, _utc1, _utc2, _dut1, _elong, _phi, _hm, _xp, _yp, _phpa, _tc, _rh, _wl, _aob, _zob, _hob, _dob, _rob, _eo) status = iternext(GetNpyIter(it)) return aob_out, zob_out, hob_out, dob_out, rob_out, eo_out -#int eraD2dtf(const char *scale, int ndp, double d1, double d2, int *iy, int *im, int *id, int ihmsf[4]) def d2dtf(scale, ndp, d1, d2): #Turn all inputs into arrays @@ -204,9 +232,9 @@ def d2dtf(scale, ndp, d1, d2): #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed broadcast = np.broadcast(ndp_in, d1_in, d2_in) - iy_out = np.empty(broadcast.shape, dtype=np.int32) - im_out = np.empty(broadcast.shape, dtype=np.int32) - id_out = np.empty(broadcast.shape, dtype=np.int32) + iy_out = np.empty(broadcast.shape, dtype=np.int32) + im_out = np.empty(broadcast.shape, dtype=np.int32) + id_out = np.empty(broadcast.shape, dtype=np.int32) ihmsf_out = np.empty(broadcast.shape+(4,), dtype=np.int32) #Create the iterator, broadcasting on all but the consumed dimensions @@ -214,6 +242,7 @@ def d2dtf(scale, ndp, d1, d2): op_flags = [['readonly']]*3+[['readwrite']]*4 it = np.nditer([ndp_in, d1_in, d2_in, iy_out, im_out, id_out, ihmsf_out], op_axes=op_axes, op_flags=op_flags) + #Iterate cdef char* _scale = scale cdef int _ndp cdef double _d1 @@ -222,93 +251,75 @@ def d2dtf(scale, ndp, d1, d2): cdef int* _im cdef int* _id cdef int* _ihmsf - - #Iterate cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 while status: - _ndp = (< int *>(dataptrarray[0]))[0] - _d1 = ((dataptrarray[1]))[0] - _d2 = ((dataptrarray[2]))[0] - _iy = ( (dataptrarray[3])) - _im = ( (dataptrarray[4])) - _id = ( (dataptrarray[5])) - _ihmsf = ( (dataptrarray[6])) + _ndp = ((dataptrarray[0]))[0] + _d1 = ((dataptrarray[1]))[0] + _d2 = ((dataptrarray[2]))[0] + _iy = ((dataptrarray[3])) + _im = ((dataptrarray[4])) + _id = ((dataptrarray[5])) + _ihmsf = ((dataptrarray[6])) ret = eraD2dtf(_scale, _ndp, _d1, _d2, _iy, _im, _id, _ihmsf) status = iternext(GetNpyIter(it)) return iy_out, im_out, id_out, ihmsf_out -#void eraAper(double theta, eraASTROM *astrom) -def aper(theta, astrom): +def rx(phi, r): #Turn all inputs into arrays - theta_in = np.array(theta, dtype=np.double, copy=False, subok=True) - astrom_in = np.array(astrom, dtype=dt_eraASTROM, copy=False, subok=True) + phi_in = np.array(phi, dtype=np.double, copy=False, subok=True) + r_in = np.array(r, dtype=np.double, copy=False, subok=True) #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed - broadcast = np.broadcast(theta_in, astrom_in) - astrom_out = np.empty(broadcast.shape, dtype=dt_eraASTROM) - np.copyto(astrom_out, astrom_in) + broadcast = np.broadcast(phi_in, r_in[...,0,0]) + r_out = np.empty(broadcast.shape+(3, 3), dtype=np.double) + np.copyto(r_out, r_in) #Create the iterator, broadcasting on all but the consumed dimensions - op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [theta_in, astrom_out]] + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [phi_in, r_out[...,0,0]]] op_flags = [['readonly']]*1+[['readwrite']]*1 - it = np.nditer([theta_in, astrom_out], op_axes=op_axes, op_flags=op_flags) - - cdef double _theta - cdef eraASTROM* _astrom + it = np.nditer([phi_in, r_out], op_axes=op_axes, op_flags=op_flags) #Iterate + cdef double _phi + cdef double* _r cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 while status: - _theta = ( (dataptrarray[0]))[0] - _astrom = ((dataptrarray[1])) - eraAper(_theta, _astrom) + _phi = ((dataptrarray[0]))[0] + _r = ((dataptrarray[1])) + eraRx(_phi, _r) status = iternext(GetNpyIter(it)) - return astrom_out + return r_out -#void eraAb(double pnat[3], double v[3], double s, double bm1, double ppr[3]) -def ab(pnat, v, s, bm1): +def ir(): #Turn all inputs into arrays - pnat_in = np.array(pnat, dtype=np.double, copy=False, subok=True) - v_in = np.array(v, dtype=np.double, copy=False, subok=True) - s_in = np.array(s, dtype=np.double, copy=False, subok=True) - bm1_in = np.array(bm1, dtype=np.double, copy=False, subok=True) #Create the output array, based on the broadcasted shape, adding the generated dimensions if needed - broadcast = np.broadcast(pnat_in[...,0], v_in[...,0], s_in, bm1_in) - ppr_out = np.empty(broadcast.shape+(3,), dtype=np.double) + broadcast = np.broadcast(np.int32(0.0), np.int32(0.0)) + r_out = np.empty(broadcast.shape+(3, 3), dtype=np.double) #Create the iterator, broadcasting on all but the consumed dimensions - op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in [pnat_in[...,0], v_in[...,0], s_in, bm1_in, ppr_out[...,0]]] - op_flags = [['readonly']]*4+[['readwrite']]*1 - it = np.nditer([pnat_in, v_in, s_in, bm1_in, ppr_out], op_axes=op_axes, op_flags=op_flags) - - cdef double* _pnat - cdef double* _v - cdef double _s - cdef double _bm1 - cdef double* _ppr + op_axes = [[-1]*(broadcast.nd-arr.ndim)+range(arr.ndim) for arr in []+[r_out[...,0,0]]] + op_flags = [['readonly']]*0+[['readwrite']]*1 + it = np.nditer([]+[r_out], op_axes=op_axes, op_flags=op_flags) #Iterate + cdef double* _r cdef char** dataptrarray = GetDataPtrArray(GetNpyIter(it)) cdef IterNextFunc iternext = GetIterNext(GetNpyIter(it), NULL) cdef int status = 1 while status: - _pnat = ((dataptrarray[0])) - _v = ((dataptrarray[1])) - _s = ((dataptrarray[2]))[0] - _bm1 = ((dataptrarray[3]))[0] - _ppr = ((dataptrarray[4])) - eraAb(_pnat, _v, _s, _bm1, _ppr) + _r = ((dataptrarray[0])) + eraIr(_r) status = iternext(GetNpyIter(it)) - return ppr_out + return r_out From f80ce83a2bdf77492058ad8372a653d4ea602c9c Mon Sep 17 00:00:00 2001 From: Julien Woillez Date: Fri, 7 Nov 2014 11:13:06 +0100 Subject: [PATCH 60/60] Added cython_npyiter_auto to setup.py --- setup.py | 8 +++++++- 1 file changed, 7 insertions(+), 1 deletion(-) diff --git a/setup.py b/setup.py index 0d84852..8f298bf 100644 --- a/setup.py +++ b/setup.py @@ -22,8 +22,14 @@ include_dirs = [np.get_include(), '/opt/local/include'], library_dirs = ['/opt/local/lib']) +cython_npyiter_auto_ext = Extension("erfa.cython_npyiter_auto", + ["cython_npyiter_auto/erfa.pyx"], + libraries=['erfa'], + include_dirs = [np.get_include(), '/opt/local/include'], + library_dirs = ['/opt/local/lib']) + setup(name = "erfa", cmdclass = { "build_ext": build_ext }, packages = ["erfa"], package_dir = {"erfa":"."}, - ext_modules = [cython_numpy_ext, cython_numpy_auto_ext, cython_npyiter_ext]) + ext_modules = [cython_numpy_ext, cython_numpy_auto_ext, cython_npyiter_ext, cython_npyiter_auto_ext])