diff options
| author | Even Rouault <even.rouault@spatialys.com> | 2018-12-18 20:24:11 +0100 |
|---|---|---|
| committer | Even Rouault <even.rouault@spatialys.com> | 2018-12-26 10:08:53 +0100 |
| commit | 610957f7035242f15743c399ffd429b92bc36206 (patch) | |
| tree | 73f0d51147e2f4860c4bfc875f7a4bf9359386d4 /src/biveval.cpp | |
| parent | 355d681ed88019e97742344bd642c2fd97e700a1 (diff) | |
| download | PROJ-610957f7035242f15743c399ffd429b92bc36206.tar.gz PROJ-610957f7035242f15743c399ffd429b92bc36206.zip | |
cpp conversion: minimal steps to fix compilation errors, not warnings
Diffstat (limited to 'src/biveval.cpp')
| -rw-r--r-- | src/biveval.cpp | 87 |
1 files changed, 87 insertions, 0 deletions
diff --git a/src/biveval.cpp b/src/biveval.cpp new file mode 100644 index 00000000..05e83e21 --- /dev/null +++ b/src/biveval.cpp @@ -0,0 +1,87 @@ +/* procedures for evaluating Tseries */ +# include "projects.h" +# define NEAR_ONE 1.00001 +static double ceval(struct PW_COEF *C, int n, projUV w, projUV w2) { + double d=0, dd=0, vd, vdd, tmp, *c; + int j; + + for (C += n ; n-- ; --C ) { + if ( (j = C->m) != 0) { + vd = vdd = 0.; + for (c = C->c + --j; j ; --j ) { + vd = w2.v * (tmp = vd) - vdd + *c--; + vdd = tmp; + } + d = w2.u * (tmp = d) - dd + w.v * vd - vdd + 0.5 * *c; + } else + d = w2.u * (tmp = d) - dd; + dd = tmp; + } + if ( (j = C->m) != 0 ) { + vd = vdd = 0.; + for (c = C->c + --j; j ; --j ) { + vd = w2.v * (tmp = vd) - vdd + *c--; + vdd = tmp; + } + return (w.u * d - dd + 0.5 * ( w.v * vd - vdd + 0.5 * *c )); + } else + return (w.u * d - dd); +} + +projUV /* bivariate Chebyshev polynomial entry point */ +bcheval(projUV in, Tseries *T) { + projUV w2, w; + projUV out; + /* scale to +-1 */ + w.u = ( in.u + in.u - T->a.u ) * T->b.u; + w.v = ( in.v + in.v - T->a.v ) * T->b.v; + if (fabs(w.u) > NEAR_ONE || fabs(w.v) > NEAR_ONE) { + out.u = out.v = HUGE_VAL; + pj_errno = PJD_ERR_TCHEBY_VAL_OUT_OF_RANGE; + } else { /* double evaluation */ + w2.u = w.u + w.u; + w2.v = w.v + w.v; + out.u = ceval(T->cu, T->mu, w, w2); + out.v = ceval(T->cv, T->mv, w, w2); + } + return out; +} + +projUV /* bivariate power polynomial entry point */ +bpseval(projUV in, Tseries *T) { + projUV out; + double *c, row; + int i, m; + + out.u = out.v = 0.; + for (i = T->mu; i >= 0; --i) { + row = 0.; + if ((m = T->cu[i].m) != 0) { + c = T->cu[i].c + m; + while (m--) + row = *--c + in.v * row; + } + out.u = row + in.u * out.u; + } + for (i = T->mv; i >= 0; --i) { + row = 0.; + if ((m = T->cv[i].m) != 0) { + c = T->cv[i].c + m; + while (m--) + row = *--c + in.v * row; + } + out.v = row + in.u * out.v; + } + return out; +} + +projUV /* general entry point selecting evaluation mode */ +biveval(projUV in, Tseries *T) { + + if (T->power) { + return bpseval(in, T); + } else { + return bcheval(in, T); + } +} + |
