OB.DAAC Logo
NASA Logo
Ocean Color Science Software

ocssw V2022
proj_molwfor.c
Go to the documentation of this file.
1 /*******************************************************************************
2 NAME MOLLWEIDE
3 
4 PURPOSE: Transforms input longitude and latitude to Easting and
5  Northing for the MOllweide projection. The
6  longitude and latitude must be in radians. The Easting
7  and Northing values will be returned in meters.
8 
9 PROGRAMMER DATE
10 ---------- ----
11 D. Steinwand, EROS May, 1991; Updated Sept, 1992; Updated Feb, 1993
12 S. Nelson, EDC Jun, 2993; Made corrections in precision and
13  number of iterations.
14 
15 ALGORITHM REFERENCES
16 
17 1. Snyder, John P. and Voxland, Philip M., "An Album of Map Projections",
18  U.S. Geological Survey Professional Paper 1453 , United State Government
19  Printing Office, Washington D.C., 1989.
20 
21 2. Snyder, John P., "Map Projections--A Working Manual", U.S. Geological
22  Survey Professional Paper 1395 (Supersedes USGS Bulletin 1532), United
23  State Government Printing Office, Washington D.C., 1987.
24  *******************************************************************************/
25 #include "proj_cproj.h"
26 #include <stdint.h>
27 
28 /* Variables common to all subroutines in this code file
29  -----------------------------------------------------*/
30 static double lon_center; /* Center longitude (projection center) */
31 static double R; /* Radius of the earth (sphere) */
32 static double false_easting; /* x offset in meters */
33 static double false_northing; /* y offset in meters */
34 
35 /* Initialize the Mollweide projection
36  ------------------------------------*/
37 int molwforint(double r, double center_long, double false_east, double false_north)
38 /* r (I) Radius of the earth (sphere) */
39 /* center_long (I) Center longitude */
40 /* false_east x offset in meters */
41 /* false_north y offset in meters */ {
42  /* Place parameters in static storage for common use
43  -------------------------------------------------*/
44  false_easting = false_east;
45  false_northing = false_north;
46  R = r;
47  lon_center = center_long;
48  return (OK);
49 }
50 
51 /* Mollweide forward equations--mapping lat,long to x,y
52  ----------------------------------------------------*/
53 int molwfor(double lon, double lat, double *x, double *y)
54 /* lon (I) Longitude */
55 /* lat (I) Latitude */
56 /* x (O) X projection coordinate */
57 /* y (O) Y projection coordinate */ {
58  double adjust_lon(double lon); /* Function to adjust longitude to -180 - 180 */
59  double delta_lon; /* Delta longitude (Given longitude - center */
60  double theta;
61  double delta_theta;
62  double con;
63  int32_t i;
64 
65  /* Forward equations
66  -----------------*/
67  delta_lon = adjust_lon(lon - lon_center);
68  theta = lat;
69  con = PI * sin(lat);
70 
71  /* Iterate using the Newton-Raphson method to find theta
72  -----------------------------------------------------*/
73  for (i = 0;; i++) {
74  delta_theta = -(theta + sin(theta) - con) / (1.0 + cos(theta));
75  theta += delta_theta;
76  if (fabs(delta_theta) < EPSLN) break;
77  if (i >= 50) {
78  /*
79  p_error("Iteration failed to converge","Mollweide-forward");
80  */
81  return (241);
82  }
83  }
84  theta /= 2.0;
85 
86  /* If the latitude is 90 deg, force the x coordinate to be "0 + false easting"
87  this is done here because of precision problems with "cos(theta)"
88  --------------------------------------------------------------------------*/
89  if (PI / 2 - fabs(lat) < EPSLN)
90  delta_lon = 0;
91  *x = 0.900316316158 * R * delta_lon * cos(theta) + false_easting;
92  *y = 1.4142135623731 * R * sin(theta) + false_northing;
93  return (OK);
94 }
int r
Definition: decode_rs.h:73
int molwforint(double r, double center_long, double false_east, double false_north)
Definition: proj_molwfor.c:37
float * lat
double adjust_lon(double x)
Definition: proj_cproj.c:349
#define PI
Definition: l3_get_org.c:6
#define OK
Definition: ancil.h:30
#define fabs(a)
Definition: misc.h:93
int molwfor(double lon, double lat, double *x, double *y)
Definition: proj_molwfor.c:53
float * lon
#define R
Definition: make_L3_v1.1.c:96
int i
Definition: decode_rs.h:71
#define EPSLN
Definition: proj_define.h:86