"Fossies" - the Fresh Open Source Software Archive

Member "proj-6.2.1/src/projections/moll.cpp" (6 May 2019, 2723 Bytes) of package /linux/privat/proj-6.2.1.tar.gz:


As a special service "Fossies" has tried to format the requested source page into HTML format using (guessed) C and C++ source code syntax highlighting (style: standard) with prefixed line numbers and code folding option. Alternatively you can here view or download the uninterpreted source code file. For more information about "moll.cpp" see the Fossies "Dox" file reference documentation.

    1 #define PJ_LIB__
    2 
    3 #include <errno.h>
    4 #include <math.h>
    5 
    6 #include "proj.h"
    7 #include "proj_internal.h"
    8 
    9 PROJ_HEAD(moll, "Mollweide") "\n\tPCyl, Sph";
   10 PROJ_HEAD(wag4, "Wagner IV") "\n\tPCyl, Sph";
   11 PROJ_HEAD(wag5, "Wagner V") "\n\tPCyl, Sph";
   12 
   13 #define MAX_ITER    10
   14 #define LOOP_TOL    1e-7
   15 
   16 namespace { // anonymous namespace
   17 struct pj_opaque {
   18     double  C_x, C_y, C_p;
   19 };
   20 } // anonymous namespace
   21 
   22 
   23 static PJ_XY moll_s_forward (PJ_LP lp, PJ *P) {           /* Spheroidal, forward */
   24     PJ_XY xy = {0.0,0.0};
   25     struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque);
   26     double k, V;
   27     int i;
   28 
   29     k = Q->C_p * sin(lp.phi);
   30     for (i = MAX_ITER; i ; --i) {
   31         lp.phi -= V = (lp.phi + sin(lp.phi) - k) /
   32             (1. + cos(lp.phi));
   33         if (fabs(V) < LOOP_TOL)
   34             break;
   35     }
   36     if (!i)
   37         lp.phi = (lp.phi < 0.) ? -M_HALFPI : M_HALFPI;
   38     else
   39         lp.phi *= 0.5;
   40     xy.x = Q->C_x * lp.lam * cos(lp.phi);
   41     xy.y = Q->C_y * sin(lp.phi);
   42     return xy;
   43 }
   44 
   45 
   46 static PJ_LP moll_s_inverse (PJ_XY xy, PJ *P) {           /* Spheroidal, inverse */
   47     PJ_LP lp = {0.0,0.0};
   48     struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque);
   49     lp.phi = aasin(P->ctx, xy.y / Q->C_y);
   50     lp.lam = xy.x / (Q->C_x * cos(lp.phi));
   51         if (fabs(lp.lam) < M_PI) {
   52             lp.phi += lp.phi;
   53             lp.phi = aasin(P->ctx, (lp.phi + sin(lp.phi)) / Q->C_p);
   54         } else {
   55             lp.lam = lp.phi = HUGE_VAL;
   56         }
   57     return lp;
   58 }
   59 
   60 
   61 static PJ * setup(PJ *P, double p) {
   62     struct pj_opaque *Q = static_cast<struct pj_opaque*>(P->opaque);
   63     double r, sp, p2 = p + p;
   64 
   65     P->es = 0;
   66     sp = sin(p);
   67     r = sqrt(M_TWOPI * sp / (p2 + sin(p2)));
   68 
   69     Q->C_x = 2. * r / M_PI;
   70     Q->C_y = r / sp;
   71     Q->C_p = p2 + sin(p2);
   72 
   73     P->inv = moll_s_inverse;
   74     P->fwd = moll_s_forward;
   75     return P;
   76 }
   77 
   78 
   79 PJ *PROJECTION(moll) {
   80     struct pj_opaque *Q = static_cast<struct pj_opaque*>(pj_calloc (1, sizeof (struct pj_opaque)));
   81     if (nullptr==Q)
   82         return pj_default_destructor (P, ENOMEM);
   83     P->opaque = Q;
   84 
   85     return setup(P, M_HALFPI);
   86 }
   87 
   88 
   89 PJ *PROJECTION(wag4) {
   90     struct pj_opaque *Q = static_cast<struct pj_opaque*>(pj_calloc (1, sizeof (struct pj_opaque)));
   91     if (nullptr==Q)
   92         return pj_default_destructor (P, ENOMEM);
   93     P->opaque = Q;
   94 
   95     return setup(P, M_PI/3.);
   96 }
   97 
   98 PJ *PROJECTION(wag5) {
   99     struct pj_opaque *Q = static_cast<struct pj_opaque*>(pj_calloc (1, sizeof (struct pj_opaque)));
  100     if (nullptr==Q)
  101         return pj_default_destructor (P, ENOMEM);
  102     P->opaque = Q;
  103 
  104     P->es = 0;
  105     Q->C_x = 0.90977;
  106     Q->C_y = 1.65014;
  107     Q->C_p = 3.00896;
  108 
  109     P->inv = moll_s_inverse;
  110     P->fwd = moll_s_forward;
  111 
  112     return P;
  113 }