GRASS GIS 8 Programmer's Manual 8.4.1(2025)-45ca3179ab
Loading...
Searching...
No Matches
area_ellipse.c
Go to the documentation of this file.
1/*!
2 * \file lib/gis/area_ellipse.c
3 *
4 * \brief GIS Library - Ellipse area routines.
5 *
6 * (C) 2001-2009 by the GRASS Development Team
7 *
8 * This program is free software under the GNU General Public License
9 * (>=v2). Read the file COPYING that comes with GRASS for details.
10 *
11 * \author Original author CERL
12 */
13
14#include <math.h>
15#include <grass/gis.h>
16#include "pi.h"
17
18static struct state {
19 double E;
20 double M;
21} state;
22
23static struct state *st = &state;
24
25/*
26 * a is semi-major axis, e2 is eccentricity squared, s is a scale factor
27 * code will fail if e2==0 (sphere)
28 */
29
30/*!
31 * \brief Begin area calculations for ellipsoid.
32 *
33 * Initializes raster area calculations for an ellipsoid, where <i>a</i>
34 * is the semi-major axis of the ellipse (in meters), <i>e2</i> is the
35 * ellipsoid eccentricity squared, and <i>s</i> is a scale factor to
36 * allow for calculations of part of the zone (<i>s</i>=1.0 is full
37 * zone, <i>s</i>=0.5 is half the zone, and <i>s</i>=360/ew_res is for a
38 * single grid cell).
39 *
40 * <b>Note:</b> <i>e2</i> must be positive. A negative value makes no
41 * sense, and zero implies a sphere.
42 *
43 * \param a semi-major axis
44 * \param e2 ellipsoid eccentricity
45 * \param s scale factor
46 */
47void G_begin_zone_area_on_ellipsoid(double a, double e2, double s)
48{
49 st->E = sqrt(e2);
50 st->M = s * a * a * M_PI * (1 - e2) / st->E;
51}
52
53/*!
54 * \brief Calculate integral for area between two latitudes.
55 *
56 * This routine is part of the integral for the area between two
57 * latitudes.
58 *
59 * \param lat latitude
60 *
61 * \return cell area
62 */
63double G_darea0_on_ellipsoid(double lat)
64{
65 double x;
66
67 x = st->E * sin(Radians(lat));
68
69 return (st->M * (x / (1.0 - x * x) + 0.5 * log((1.0 + x) / (1.0 - x))));
70}
71
72/*!
73 * \brief Calculates area between latitudes.
74 *
75 * This routine shows how to calculate area between two lats, but
76 * isn't efficient for row by row since G_darea0_on_ellipsoid()
77 * will be called twice for the same lat, once as a <i>south</i> then
78 * again as a <i>north</i>.
79 *
80 * Returns the area between latitudes <i>north</i> and <i>south</i>
81 * scaled by the factor <i>s</i> passed to
82 * G_begin_zone_area_on_ellipsoid().
83 *
84 * \param north north coordinate
85 * \param south south coordinate
86 *
87 * \return cell area
88 */
89double G_area_for_zone_on_ellipsoid(double north, double south)
90{
91 return (G_darea0_on_ellipsoid(north) - G_darea0_on_ellipsoid(south));
92}
double G_darea0_on_ellipsoid(double lat)
Calculate integral for area between two latitudes.
double G_area_for_zone_on_ellipsoid(double north, double south)
Calculates area between latitudes.
void G_begin_zone_area_on_ellipsoid(double a, double e2, double s)
Begin area calculations for ellipsoid.
#define Radians(x)
Definition pi.h:6
#define x