00001
00030 #include <math.h>
00031 #include <grass/gis.h>
00032 #include "pi.h"
00033
00034
00035 static int adjust_lat(double *);
00036
00037 #if 0
00038 static int adjust_lon(double *);
00039 #endif
00040
00041 static double TAN_A, TAN1, TAN2, L;
00042 static int parallel;
00043
00044
00057 int G_begin_rhumbline_equation(double lon1, double lat1, double lon2,
00058 double lat2)
00059 {
00060 adjust_lat(&lat1);
00061 adjust_lat(&lat2);
00062
00063 if (lon1 == lon2) {
00064 parallel = 1;
00065 L = lat1;
00066 return 0;
00067 }
00068 if (lat1 == lat2) {
00069 parallel = 1;
00070 L = lat1;
00071 return 1;
00072 }
00073 parallel = 0;
00074 lon1 = Radians(lon1);
00075 lon2 = Radians(lon2);
00076 lat1 = Radians(lat1);
00077 lat2 = Radians(lat2);
00078
00079 TAN1 = tan(M_PI_4 + lat1 / 2.0);
00080 TAN2 = tan(M_PI_4 + lat2 / 2.0);
00081 TAN_A = (lon2 - lon1) / (log(TAN2) - log(TAN1));
00082 L = lon1;
00083
00084 return 1;
00085 }
00086
00087
00097 double G_rhumbline_lat_from_lon(double lon)
00098 {
00099 if (parallel)
00100 return L;
00101
00102 lon = Radians(lon);
00103
00104 return Degrees(2 * atan(exp((lon - L) / TAN_A) * TAN1) - M_PI_2);
00105 }
00106
00107
00108 #if 0
00109 static int adjust_lon(double *lon)
00110 {
00111 while (*lon > 180.0)
00112 *lon -= 360.0;
00113 while (*lon < -180.0)
00114 *lon += 360.0;
00115
00116 return 0;
00117 }
00118 #endif
00119
00120
00121 static int adjust_lat(double *lat)
00122 {
00123 if (*lat > 90.0)
00124 *lat = 90.0;
00125 if (*lat < -90.0)
00126 *lat = -90.0;
00127
00128 return 0;
00129 }