86
86
<a name="l00061"></a>00061 <span class="comment"> * 24-May-99 Original Code</span>
87
87
<a name="l00062"></a>00062 <span class="comment"> * 09-Jan-06 Added new geoid height interpolation methods</span>
88
88
<a name="l00063"></a>00063 <span class="comment"> * 03-14-07 Original C++ Code</span>
89
<a name="l00064"></a>00064 <span class="comment"> *</span>
90
<a name="l00065"></a>00065 <span class="comment"> */</span>
91
<a name="l00066"></a>00066
92
<a name="l00067"></a>00067
93
<a name="l00068"></a>00068 <span class="comment">/***************************************************************************/</span>
94
<a name="l00069"></a>00069 <span class="comment">/*</span>
95
<a name="l00070"></a>00070 <span class="comment"> * INCLUDES</span>
96
<a name="l00071"></a>00071 <span class="comment"> */</span>
97
<a name="l00072"></a>00072
98
<a name="l00073"></a>00073 <span class="preprocessor">#include <string.h></span>
99
<a name="l00074"></a>00074 <span class="preprocessor">#include <stdlib.h></span>
100
<a name="l00075"></a>00075 <span class="preprocessor">#include <stdio.h></span>
101
<a name="l00076"></a>00076 <span class="preprocessor">#include "<a class="code" href="_geoid_library_8h.html">GeoidLibrary.h</a>"</span>
102
<a name="l00077"></a>00077 <span class="preprocessor">#include "<a class="code" href="threads_8h.html">threads.h</a>"</span>
103
<a name="l00078"></a>00078 <span class="preprocessor">#include "<a class="code" href="_coordinate_conversion_exception_8h.html">CoordinateConversionException.h</a>"</span>
104
<a name="l00079"></a>00079 <span class="preprocessor">#include "<a class="code" href="_error_messages_8h.html">ErrorMessages.h</a>"</span>
105
<a name="l00080"></a>00080
106
<a name="l00081"></a>00081 <span class="comment">/*</span>
107
<a name="l00082"></a>00082 <span class="comment"> * string.h - standard C string handling library</span>
108
<a name="l00083"></a>00083 <span class="comment"> * stdio.h - standard C input/output library</span>
109
<a name="l00084"></a>00084 <span class="comment"> * stdlib.h - standard C general utilities library</span>
110
<a name="l00085"></a>00085 <span class="comment"> * GeoidLibrary.h - prototype error checking and error codes</span>
111
<a name="l00086"></a>00086 <span class="comment"> * threads.h - used for thread safety</span>
112
<a name="l00087"></a>00087 <span class="comment"> * CoordinateConversionException.h - Exception handler</span>
113
<a name="l00088"></a>00088 <span class="comment"> * ErrorMessages.h - Contains exception messages</span>
114
<a name="l00089"></a>00089 <span class="comment"> */</span>
89
<a name="l00064"></a>00064 <span class="comment"> * 05-12-10 S. Gillis, BAEts26542, MSP TS MSL-HAE conversion </span>
90
<a name="l00065"></a>00065 <span class="comment"> * should use CCS </span>
91
<a name="l00066"></a>00066 <span class="comment"> * 06-11-10 S. Gillis, BAEts26724, Fixed memory error problem</span>
92
<a name="l00067"></a>00067 <span class="comment"> * when MSPCCS_DATA is not set</span>
93
<a name="l00068"></a>00068 <span class="comment"> * 06-16-10 S. Gillis, BAEts27144, Fixed memory error problem</span>
94
<a name="l00069"></a>00069 <span class="comment"> * when MSPCCS_DATA is not set</span>
95
<a name="l00070"></a>00070 <span class="comment"> * 07-07-10 K.Lam, BAEts27269, Replace C functions in threads.h</span>
96
<a name="l00071"></a>00071 <span class="comment"> * with C++ methods in classes CCSThreadMutex</span>
97
<a name="l00072"></a>00072 <span class="comment"> */</span>
98
<a name="l00073"></a>00073
99
<a name="l00074"></a>00074
100
<a name="l00075"></a>00075 <span class="comment">/***************************************************************************/</span>
101
<a name="l00076"></a>00076 <span class="comment">/*</span>
102
<a name="l00077"></a>00077 <span class="comment"> * INCLUDES</span>
103
<a name="l00078"></a>00078 <span class="comment"> */</span>
104
<a name="l00079"></a>00079
105
<a name="l00080"></a>00080 <span class="preprocessor">#include <string.h></span>
106
<a name="l00081"></a>00081 <span class="preprocessor">#include <stdlib.h></span>
107
<a name="l00082"></a>00082 <span class="preprocessor">#include <stdio.h></span>
108
<a name="l00083"></a>00083 <span class="preprocessor">#include "<a class="code" href="_geoid_library_8h.html">GeoidLibrary.h</a>"</span>
109
<a name="l00084"></a>00084 <span class="preprocessor">#include "<a class="code" href="_coordinate_conversion_exception_8h.html">CoordinateConversionException.h</a>"</span>
110
<a name="l00085"></a>00085 <span class="preprocessor">#include "<a class="code" href="_error_messages_8h.html">ErrorMessages.h</a>"</span>
111
<a name="l00086"></a>00086 <span class="preprocessor">#include "<a class="code" href="_c_c_s_thread_mutex_8h.html">CCSThreadMutex.h</a>"</span>
112
<a name="l00087"></a>00087 <span class="preprocessor">#include "<a class="code" href="_c_c_s_thread_lock_8h.html">CCSThreadLock.h</a>"</span>
113
<a name="l00088"></a>00088
114
<a name="l00089"></a>00089 <span class="preprocessor">#include <vector></span>
115
115
<a name="l00090"></a>00090
116
<a name="l00091"></a>00091
117
<a name="l00092"></a>00092 <span class="keyword">using namespace </span>MSP::CCS;
118
<a name="l00093"></a>00093
119
<a name="l00094"></a>00094
120
<a name="l00095"></a>00095 <span class="comment">/***************************************************************************/</span>
121
<a name="l00096"></a>00096 <span class="comment">/* DEFINES </span>
122
<a name="l00097"></a>00097 <span class="comment"> *</span>
123
<a name="l00098"></a>00098 <span class="comment"> */</span>
124
<a name="l00099"></a>00099
125
<a name="l00100"></a><a class="code" href="_geoid_library_8cpp.html#a952eac791b596a61bba0a133a3bb439f">00100</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> = 3.14159265358979323e0; <span class="comment">/* PI */</span>
126
<a name="l00101"></a><a class="code" href="_geoid_library_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">00101</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> = <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> / 2.0;
127
<a name="l00102"></a><a class="code" href="_geoid_library_8cpp.html#a820597b124334bf7fb85214114f4876d">00102</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_albers_equal_area_conic_8cpp.html#a820597b124334bf7fb85214114f4876d">TWO_PI</a> = 2.0 * <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a>;
128
<a name="l00103"></a><a class="code" href="_geoid_library_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">00103</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a> = (180.0 / <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a>);
129
<a name="l00104"></a><a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">00104</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a> = 1441; <span class="comment">/* 360 degrees of longitude at 15 minute spacing */</span>
130
<a name="l00105"></a><a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">00105</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a> = 721; <span class="comment">/* 180 degrees of latitude at 15 minute spacing */</span>
131
<a name="l00106"></a><a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">00106</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a> = 37; <span class="comment">/* 360 degrees of longitude at 10 degree spacing */</span>
132
<a name="l00107"></a><a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">00107</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a> = 19; <span class="comment">/* 180 degrees of latitude at 10 degree spacing */</span>
133
<a name="l00108"></a><a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">00108</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">EGM96_HEADER_ITEMS</a> = 6; <span class="comment">/* min, max lat, min, max long, lat, long spacing*/</span>
134
<a name="l00109"></a><a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">00109</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a> = .25; <span class="comment">/* 4 grid cells per degree at 15 minute spacing */</span>
135
<a name="l00110"></a><a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">00110</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a> = 10; <span class="comment">/* 1 / 10.0 grid cells per degree at 10 degree spacing */</span>
136
<a name="l00111"></a><a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">00111</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a> = .5; <span class="comment">/* 2 grid cells per degree at 30 minute spacing */</span>
137
<a name="l00112"></a><a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">00112</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a> = 1; <span class="comment">/* 1 grid cell per degree at 1 degree spacing */</span>
138
<a name="l00113"></a><a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">00113</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a> = 2; <span class="comment">/* 1 / 2 grid cells per degree at 2 degree spacing */</span>
139
<a name="l00114"></a><a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">00114</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a> = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a> * <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
140
<a name="l00115"></a><a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">00115</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a> = <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a> * <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>;
141
<a name="l00116"></a><a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">00116</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a> = 53;
142
<a name="l00117"></a>00117
143
<a name="l00118"></a>00118
144
<a name="l00119"></a>00119 <span class="comment">/* defines the egm96 variable grid */</span>
145
<a name="l00120"></a><a class="code" href="struct_e_g_m96___variable___grid.html">00120</a> <span class="keyword">struct </span><a class="code" href="struct_e_g_m96___variable___grid.html">EGM96_Variable_Grid</a>
146
<a name="l00121"></a>00121 {
147
<a name="l00122"></a><a class="code" href="struct_e_g_m96___variable___grid.html#adc0cb102acb3063f76beba066fb0669c">00122</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#adc0cb102acb3063f76beba066fb0669c">min_lat</a>; <span class="comment">/* Min cell latitude */</span>
148
<a name="l00123"></a><a class="code" href="struct_e_g_m96___variable___grid.html#afe14a01a2dcaf81a06e7c9aeac96dcf9">00123</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#afe14a01a2dcaf81a06e7c9aeac96dcf9">max_lat</a>; <span class="comment">/* Max cell latitude */</span>
149
<a name="l00124"></a><a class="code" href="struct_e_g_m96___variable___grid.html#a62917673950dfbb44e78d36e76d97795">00124</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#a62917673950dfbb44e78d36e76d97795">min_lon</a>; <span class="comment">/* Min cell longitude */</span>
150
<a name="l00125"></a><a class="code" href="struct_e_g_m96___variable___grid.html#a59e9e4f3c67a6701f968cc7aa0a3f235">00125</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#a59e9e4f3c67a6701f968cc7aa0a3f235">max_lon</a>; <span class="comment">/* Max cell latitude */</span>
151
<a name="l00126"></a>00126 };
152
<a name="l00127"></a>00127
153
<a name="l00128"></a>00128 <span class="comment">/* Table of EGM96 variable grid 30 minute inset areas */</span>
154
<a name="l00129"></a><a class="code" href="_geoid_library_8cpp.html#af4dbe9967590d5fd26163e5f66e19d3d">00129</a> <span class="keyword">const</span> <a class="code" href="struct_e_g_m96___variable___grid.html">EGM96_Variable_Grid</a> <a class="code" href="_geoid_library_8cpp.html#af4dbe9967590d5fd26163e5f66e19d3d">EGM96_Variable_Grid_Table</a>[<a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a>] =
155
<a name="l00130"></a>00130 {{74.5, 75.5, 273.5, 280.0},
156
<a name="l00131"></a>00131 {66.5, 67.5, 293.5, 295.0},
157
<a name="l00132"></a>00132 {62.5, 64.0, 133.0, 136.5},
158
<a name="l00133"></a>00133 {60.5, 61.5, 208.5, 210.0},
159
<a name="l00134"></a>00134 {60.5, 61.0, 219.0, 220.5},
160
<a name="l00135"></a>00135 {51.0, 53.0, 172.0, 174.5},
161
<a name="l00136"></a>00136 {52.0, 53.0, 192.5, 194.0},
162
<a name="l00137"></a>00137 {51.0, 52.0, 188.5, 191.0},
163
<a name="l00138"></a>00138 {50.0, 52.0, 178.0, 182.5},
164
<a name="l00139"></a>00139 {43.0, 46.0, 148.0, 153.5},
165
<a name="l00140"></a>00140 {43.0, 45.0, 84.0, 89.5},
166
<a name="l00141"></a>00141 {40.0, 41.0, 70.5, 72.0},
167
<a name="l00142"></a>00142 {36.5, 37.0, 78.5, 79.0},
168
<a name="l00143"></a>00143 {36.0, 37.0, 348.0, 349.5},
169
<a name="l00144"></a>00144 {35.0, 36.0, 171.0, 172.5},
170
<a name="l00145"></a>00145 {34.0, 35.0, 140.5, 142.0},
171
<a name="l00146"></a>00146 {29.5, 31.0, 78.5, 81.0},
172
<a name="l00147"></a>00147 {28.5, 30.0, 81.5, 83.0},
173
<a name="l00148"></a>00148 {26.5, 30.0, 142.0, 143.5},
174
<a name="l00149"></a>00149 {26.0, 29.0, 91.5, 96.0},
175
<a name="l00150"></a>00150 {27.5, 29.0, 84.0, 86.5},
176
<a name="l00151"></a>00151 {28.0, 29.0, 342.5, 344.0},
177
<a name="l00152"></a>00152 {26.5, 28.0, 88.5, 90.0},
178
<a name="l00153"></a>00153 {25.0, 26.0, 189.0, 190.5},
179
<a name="l00154"></a>00154 {23.0, 24.0, 195.0, 196.5},
180
<a name="l00155"></a>00155 {21.0, 21.5, 204.0, 204.5},
181
<a name="l00156"></a>00156 {20.0, 21.0, 283.5, 288.0},
182
<a name="l00157"></a>00157 {18.5, 20.5, 204.0, 205.5},
183
<a name="l00158"></a>00158 {18.0, 20.0, 291.0, 296.5},
184
<a name="l00159"></a>00159 {17.0, 18.0, 298.0, 299.5},
185
<a name="l00160"></a>00160 {15.0, 16.0, 122.0, 123.5},
186
<a name="l00161"></a>00161 {12.0, 14.0, 144.5, 147.0},
187
<a name="l00162"></a>00162 {11.0, 12.0, 141.5, 144.0},
188
<a name="l00163"></a>00163 {9.5, 11.5, 125.0, 127.5},
189
<a name="l00164"></a>00164 {10.0, 11.0, 286.0, 287.5},
190
<a name="l00165"></a>00165 {6.0, 9.5, 287.0, 289.5},
191
<a name="l00166"></a>00166 {5.0, 7.0, 124.0, 128.5},
192
<a name="l00167"></a>00167 {-1.0, 1.0, 125.0, 128.5},
193
<a name="l00168"></a>00168 {-3.0, -1.5, 281.0, 282.5},
194
<a name="l00169"></a>00169 {-7.0, -5.0, 150.5, 155.0},
195
<a name="l00170"></a>00170 {-8.0, -7.0, 107.0, 108.5},
196
<a name="l00171"></a>00171 {-9.0, -7.0, 147.0, 149.5},
197
<a name="l00172"></a>00172 {-11.0, -10.0, 161.5, 163.0},
198
<a name="l00173"></a>00173 {-14.5, -13.5, 166.0, 167.5},
199
<a name="l00174"></a>00174 {-18.5, -17.0, 186.5, 188.0},
200
<a name="l00175"></a>00175 {-20.5, -20.0, 168.0, 169.5},
201
<a name="l00176"></a>00176 {-23.0, -20.0, 184.5, 187.0},
202
<a name="l00177"></a>00177 {-27.0, -24.0, 288.0, 290.5},
203
<a name="l00178"></a>00178 {-53.0, -52.0, 312.0, 313.5},
204
<a name="l00179"></a>00179 {-56.0, -55.0, 333.0, 334.5},
205
<a name="l00180"></a>00180 {-61.5, -60.0, 312.5, 317.0},
206
<a name="l00181"></a>00181 {-61.5, -60.5, 300.5, 303.0},
207
<a name="l00182"></a>00182 {-73.0, -72.0, 24.5, 26.0}};
208
<a name="l00183"></a>00183
209
<a name="l00184"></a>00184
210
<a name="l00185"></a>00185 <span class="comment">/************************************************************************/</span>
211
<a name="l00186"></a>00186 <span class="comment">/* LOCAL FUNCTIONS </span>
212
<a name="l00187"></a>00187 <span class="comment"> *</span>
213
<a name="l00188"></a>00188 <span class="comment"> */</span>
214
<a name="l00189"></a>00189
215
<a name="l00190"></a><a class="code" href="_geoid_library_8cpp.html#ac0f3a476600c87d7296394be3c5b047d">00190</a> <span class="keywordtype">float</span> <a class="code" href="_geoid_library_8cpp.html#ac0f3a476600c87d7296394be3c5b047d">readFloat</a>( <span class="keywordtype">int</span> *numRead, FILE* file )
216
<a name="l00191"></a>00191 <span class="comment">/*</span>
217
<a name="l00192"></a>00192 <span class="comment"> * The private function readFloat returns the geoid height</span>
218
<a name="l00193"></a>00193 <span class="comment"> * read from the geoid file. 4 bytes are read from the file and,</span>
219
<a name="l00194"></a>00194 <span class="comment"> * if necessary, the bytes are swapped.</span>
220
<a name="l00195"></a>00195 <span class="comment"> *</span>
221
<a name="l00196"></a>00196 <span class="comment"> * numRead : Number of heights read from file (output)</span>
222
<a name="l00197"></a>00197 <span class="comment"> *</span>
223
<a name="l00198"></a>00198 <span class="comment"> */</span>
224
<a name="l00199"></a>00199 {
225
<a name="l00200"></a>00200 <span class="keywordtype">float</span> result;
226
<a name="l00201"></a>00201 <span class="keywordtype">char</span>* swap = ( <span class="keywordtype">char</span>* )&result;
227
<a name="l00202"></a>00202 <span class="keywordtype">char</span> temp;
228
<a name="l00203"></a>00203
229
<a name="l00204"></a>00204 *numRead = fread( swap, 4, 1, file );
230
<a name="l00205"></a>00205 <span class="preprocessor">#ifdef LITTLE_ENDIAN</span>
231
<a name="l00206"></a>00206 <span class="preprocessor"></span> temp = swap[0];
232
<a name="l00207"></a>00207 swap[0] = swap[3];
233
<a name="l00208"></a>00208 swap[3] = temp;
234
<a name="l00209"></a>00209 temp = swap[1];
235
<a name="l00210"></a>00210 swap[1] = swap[2];
236
<a name="l00211"></a>00211 swap[2] = temp;
237
<a name="l00212"></a>00212 <span class="preprocessor">#endif</span>
238
<a name="l00213"></a>00213 <span class="preprocessor"></span>
239
<a name="l00214"></a>00214 <span class="keywordflow">return</span> result;
240
<a name="l00215"></a>00215 }
241
<a name="l00216"></a>00216
242
<a name="l00217"></a>00217
243
<a name="l00218"></a>00218 <span class="comment">/************************************************************************/</span>
244
<a name="l00219"></a>00219 <span class="comment">/* FUNCTIONS </span>
245
<a name="l00220"></a>00220 <span class="comment"> *</span>
246
<a name="l00221"></a>00221 <span class="comment"> */</span>
116
<a name="l00091"></a>00091 <span class="comment">/*</span>
117
<a name="l00092"></a>00092 <span class="comment"> * string.h - standard C string handling library</span>
118
<a name="l00093"></a>00093 <span class="comment"> * stdio.h - standard C input/output library</span>
119
<a name="l00094"></a>00094 <span class="comment"> * stdlib.h - standard C general utilities library</span>
120
<a name="l00095"></a>00095 <span class="comment"> * GeoidLibrary.h - prototype error checking and error codes</span>
121
<a name="l00096"></a>00096 <span class="comment"> * threads.h - used for thread safety</span>
122
<a name="l00097"></a>00097 <span class="comment"> * CoordinateConversionException.h - Exception handler</span>
123
<a name="l00098"></a>00098 <span class="comment"> * ErrorMessages.h - Contains exception messages</span>
124
<a name="l00099"></a>00099 <span class="comment"> */</span>
125
<a name="l00100"></a>00100
126
<a name="l00101"></a>00101
127
<a name="l00102"></a>00102 <span class="keyword">using namespace </span>MSP::CCS;
128
<a name="l00103"></a>00103 <span class="keyword">using</span> <a class="code" href="class_m_s_p_1_1_c_c_s_thread_mutex.html">MSP::CCSThreadMutex</a>;
129
<a name="l00104"></a>00104 <span class="keyword">using</span> <a class="code" href="class_m_s_p_1_1_c_c_s_thread_lock.html">MSP::CCSThreadLock</a>;
130
<a name="l00105"></a>00105
131
<a name="l00106"></a>00106
132
<a name="l00107"></a>00107 <span class="comment">/***************************************************************************/</span>
133
<a name="l00108"></a>00108 <span class="comment">/* DEFINES </span>
134
<a name="l00109"></a>00109 <span class="comment"> *</span>
135
<a name="l00110"></a>00110 <span class="comment"> */</span>
136
<a name="l00111"></a>00111
137
<a name="l00112"></a><a class="code" href="_geoid_library_8cpp.html#a952eac791b596a61bba0a133a3bb439f">00112</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> = 3.14159265358979323e0; <span class="comment">/* PI */</span>
138
<a name="l00113"></a><a class="code" href="_geoid_library_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">00113</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> = <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> / 2.0;
139
<a name="l00114"></a><a class="code" href="_geoid_library_8cpp.html#a820597b124334bf7fb85214114f4876d">00114</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_albers_equal_area_conic_8cpp.html#a820597b124334bf7fb85214114f4876d">TWO_PI</a> = 2.0 * <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a>;
140
<a name="l00115"></a><a class="code" href="_geoid_library_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">00115</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a> = (180.0 / <a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a>);
141
<a name="l00116"></a><a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">00116</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a> = 1441; <span class="comment">/* 360 degrees of longitude at 15 minute spacing */</span>
142
<a name="l00117"></a><a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">00117</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a> = 721; <span class="comment">/* 180 degrees of latitude at 15 minute spacing */</span>
143
<a name="l00118"></a><a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">00118</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a> = 37; <span class="comment">/* 360 degrees of longitude at 10 degree spacing */</span>
144
<a name="l00119"></a><a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">00119</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a> = 19; <span class="comment">/* 180 degrees of latitude at 10 degree spacing */</span>
145
<a name="l00120"></a><a class="code" href="_geoid_library_8cpp.html#a580782fa058e200088161c1680921220">00120</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a580782fa058e200088161c1680921220">EGM84_30_MIN_COLS</a> = 721; <span class="comment">/* 360 degrees of longitude at 30 minute spacing */</span>
146
<a name="l00121"></a><a class="code" href="_geoid_library_8cpp.html#a1ea8fba0f5bc5035bc21da731b6d8799">00121</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a1ea8fba0f5bc5035bc21da731b6d8799">EGM84_30_MIN_ROWS</a> = 361; <span class="comment">/* 180 degrees of latitude at 30 minute spacing */</span>
147
<a name="l00122"></a><a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">00122</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">EGM96_HEADER_ITEMS</a> = 6; <span class="comment">/* min, max lat, min, max long, lat, long spacing*/</span>
148
<a name="l00123"></a><a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">00123</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a> = .25; <span class="comment">/* 4 grid cells per degree at 15 minute spacing */</span>
149
<a name="l00124"></a><a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">00124</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a> = 10; <span class="comment">/* 1 / 10.0 grid cells per degree at 10 degree spacing */</span>
150
<a name="l00125"></a><a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">00125</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a> = .5; <span class="comment">/* 2 grid cells per degree at 30 minute spacing */</span>
151
<a name="l00126"></a><a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">00126</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a> = 1; <span class="comment">/* 1 grid cell per degree at 1 degree spacing */</span>
152
<a name="l00127"></a><a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">00127</a> <span class="keyword">const</span> <span class="keywordtype">double</span> <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a> = 2; <span class="comment">/* 1 / 2 grid cells per degree at 2 degree spacing */</span>
153
<a name="l00128"></a><a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">00128</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a> = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a> * <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
154
<a name="l00129"></a><a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">00129</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a> = <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a> * <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>;
155
<a name="l00130"></a><a class="code" href="_geoid_library_8cpp.html#aa3581e7c236f1d903b4b3ccf242d1088">00130</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#aa3581e7c236f1d903b4b3ccf242d1088">EGM84_30_MIN_ELEVATIONS</a> = <a class="code" href="_geoid_library_8cpp.html#a580782fa058e200088161c1680921220">EGM84_30_MIN_COLS</a> * <a class="code" href="_geoid_library_8cpp.html#a1ea8fba0f5bc5035bc21da731b6d8799">EGM84_30_MIN_ROWS</a>;
156
<a name="l00131"></a><a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">00131</a> <span class="keyword">const</span> <span class="keywordtype">int</span> <a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a> = 53;
157
<a name="l00132"></a>00132
158
<a name="l00133"></a>00133
159
<a name="l00134"></a>00134 <span class="comment">/* defines the egm96 variable grid */</span>
160
<a name="l00135"></a><a class="code" href="struct_e_g_m96___variable___grid.html">00135</a> <span class="keyword">struct </span><a class="code" href="struct_e_g_m96___variable___grid.html">EGM96_Variable_Grid</a>
161
<a name="l00136"></a>00136 {
162
<a name="l00137"></a><a class="code" href="struct_e_g_m96___variable___grid.html#adc0cb102acb3063f76beba066fb0669c">00137</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#adc0cb102acb3063f76beba066fb0669c">min_lat</a>; <span class="comment">/* Min cell latitude */</span>
163
<a name="l00138"></a><a class="code" href="struct_e_g_m96___variable___grid.html#afe14a01a2dcaf81a06e7c9aeac96dcf9">00138</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#afe14a01a2dcaf81a06e7c9aeac96dcf9">max_lat</a>; <span class="comment">/* Max cell latitude */</span>
164
<a name="l00139"></a><a class="code" href="struct_e_g_m96___variable___grid.html#a62917673950dfbb44e78d36e76d97795">00139</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#a62917673950dfbb44e78d36e76d97795">min_lon</a>; <span class="comment">/* Min cell longitude */</span>
165
<a name="l00140"></a><a class="code" href="struct_e_g_m96___variable___grid.html#a59e9e4f3c67a6701f968cc7aa0a3f235">00140</a> <span class="keywordtype">double</span> <a class="code" href="struct_e_g_m96___variable___grid.html#a59e9e4f3c67a6701f968cc7aa0a3f235">max_lon</a>; <span class="comment">/* Max cell latitude */</span>
166
<a name="l00141"></a>00141 };
167
<a name="l00142"></a>00142
168
<a name="l00143"></a>00143 <span class="comment">/* Table of EGM96 variable grid 30 minute inset areas */</span>
169
<a name="l00144"></a><a class="code" href="_geoid_library_8cpp.html#af4dbe9967590d5fd26163e5f66e19d3d">00144</a> <span class="keyword">const</span> <a class="code" href="struct_e_g_m96___variable___grid.html">EGM96_Variable_Grid</a> <a class="code" href="_geoid_library_8cpp.html#af4dbe9967590d5fd26163e5f66e19d3d">EGM96_Variable_Grid_Table</a>[<a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a>] =
170
<a name="l00145"></a>00145 {{74.5, 75.5, 273.5, 280.0},
171
<a name="l00146"></a>00146 {66.5, 67.5, 293.5, 295.0},
172
<a name="l00147"></a>00147 {62.5, 64.0, 133.0, 136.5},
173
<a name="l00148"></a>00148 {60.5, 61.5, 208.5, 210.0},
174
<a name="l00149"></a>00149 {60.5, 61.0, 219.0, 220.5},
175
<a name="l00150"></a>00150 {51.0, 53.0, 172.0, 174.5},
176
<a name="l00151"></a>00151 {52.0, 53.0, 192.5, 194.0},
177
<a name="l00152"></a>00152 {51.0, 52.0, 188.5, 191.0},
178
<a name="l00153"></a>00153 {50.0, 52.0, 178.0, 182.5},
179
<a name="l00154"></a>00154 {43.0, 46.0, 148.0, 153.5},
180
<a name="l00155"></a>00155 {43.0, 45.0, 84.0, 89.5},
181
<a name="l00156"></a>00156 {40.0, 41.0, 70.5, 72.0},
182
<a name="l00157"></a>00157 {36.5, 37.0, 78.5, 79.0},
183
<a name="l00158"></a>00158 {36.0, 37.0, 348.0, 349.5},
184
<a name="l00159"></a>00159 {35.0, 36.0, 171.0, 172.5},
185
<a name="l00160"></a>00160 {34.0, 35.0, 140.5, 142.0},
186
<a name="l00161"></a>00161 {29.5, 31.0, 78.5, 81.0},
187
<a name="l00162"></a>00162 {28.5, 30.0, 81.5, 83.0},
188
<a name="l00163"></a>00163 {26.5, 30.0, 142.0, 143.5},
189
<a name="l00164"></a>00164 {26.0, 29.0, 91.5, 96.0},
190
<a name="l00165"></a>00165 {27.5, 29.0, 84.0, 86.5},
191
<a name="l00166"></a>00166 {28.0, 29.0, 342.5, 344.0},
192
<a name="l00167"></a>00167 {26.5, 28.0, 88.5, 90.0},
193
<a name="l00168"></a>00168 {25.0, 26.0, 189.0, 190.5},
194
<a name="l00169"></a>00169 {23.0, 24.0, 195.0, 196.5},
195
<a name="l00170"></a>00170 {21.0, 21.5, 204.0, 204.5},
196
<a name="l00171"></a>00171 {20.0, 21.0, 283.5, 288.0},
197
<a name="l00172"></a>00172 {18.5, 20.5, 204.0, 205.5},
198
<a name="l00173"></a>00173 {18.0, 20.0, 291.0, 296.5},
199
<a name="l00174"></a>00174 {17.0, 18.0, 298.0, 299.5},
200
<a name="l00175"></a>00175 {15.0, 16.0, 122.0, 123.5},
201
<a name="l00176"></a>00176 {12.0, 14.0, 144.5, 147.0},
202
<a name="l00177"></a>00177 {11.0, 12.0, 141.5, 144.0},
203
<a name="l00178"></a>00178 {9.5, 11.5, 125.0, 127.5},
204
<a name="l00179"></a>00179 {10.0, 11.0, 286.0, 287.5},
205
<a name="l00180"></a>00180 {6.0, 9.5, 287.0, 289.5},
206
<a name="l00181"></a>00181 {5.0, 7.0, 124.0, 128.5},
207
<a name="l00182"></a>00182 {-1.0, 1.0, 125.0, 128.5},
208
<a name="l00183"></a>00183 {-3.0, -1.5, 281.0, 282.5},
209
<a name="l00184"></a>00184 {-7.0, -5.0, 150.5, 155.0},
210
<a name="l00185"></a>00185 {-8.0, -7.0, 107.0, 108.5},
211
<a name="l00186"></a>00186 {-9.0, -7.0, 147.0, 149.5},
212
<a name="l00187"></a>00187 {-11.0, -10.0, 161.5, 163.0},
213
<a name="l00188"></a>00188 {-14.5, -13.5, 166.0, 167.5},
214
<a name="l00189"></a>00189 {-18.5, -17.0, 186.5, 188.0},
215
<a name="l00190"></a>00190 {-20.5, -20.0, 168.0, 169.5},
216
<a name="l00191"></a>00191 {-23.0, -20.0, 184.5, 187.0},
217
<a name="l00192"></a>00192 {-27.0, -24.0, 288.0, 290.5},
218
<a name="l00193"></a>00193 {-53.0, -52.0, 312.0, 313.5},
219
<a name="l00194"></a>00194 {-56.0, -55.0, 333.0, 334.5},
220
<a name="l00195"></a>00195 {-61.5, -60.0, 312.5, 317.0},
221
<a name="l00196"></a>00196 {-61.5, -60.5, 300.5, 303.0},
222
<a name="l00197"></a>00197 {-73.0, -72.0, 24.5, 26.0}};
223
<a name="l00198"></a>00198
224
<a name="l00199"></a>00199
225
<a name="l00200"></a>00200 <span class="comment">/************************************************************************/</span>
226
<a name="l00201"></a>00201 <span class="comment">/* LOCAL FUNCTIONS </span>
227
<a name="l00202"></a>00202 <span class="comment"> *</span>
228
<a name="l00203"></a>00203 <span class="comment"> */</span>
229
<a name="l00204"></a>00204
230
<a name="l00205"></a><a class="code" href="_geoid_library_8cpp.html#a7b114915125b927e278ce213af694a31">00205</a> <span class="keywordtype">void</span> <a class="code" href="_geoid_library_8cpp.html#a7b114915125b927e278ce213af694a31">swapBytes</a>(
231
<a name="l00206"></a>00206 <span class="keywordtype">void</span> *buffer,
232
<a name="l00207"></a>00207 <span class="keywordtype">size_t</span> size,
233
<a name="l00208"></a>00208 <span class="keywordtype">size_t</span> count)
234
<a name="l00209"></a>00209 {
235
<a name="l00210"></a>00210 <span class="keywordtype">char</span> *b = (<span class="keywordtype">char</span> *) buffer;
236
<a name="l00211"></a>00211 <span class="keywordtype">char</span> temp;
237
<a name="l00212"></a>00212 <span class="keywordflow">for</span> (<span class="keywordtype">size_t</span> c = 0; c < count; c ++)
238
<a name="l00213"></a>00213 {
239
<a name="l00214"></a>00214 <span class="keywordflow">for</span> (<span class="keywordtype">size_t</span> s = 0; s < size / 2; s ++)
240
<a name="l00215"></a>00215 {
241
<a name="l00216"></a>00216 temp = b[c*size + s];
242
<a name="l00217"></a>00217 b[c*size + s] = b[c*size + size - s - 1];
243
<a name="l00218"></a>00218 b[c*size + size - s - 1] = temp;
244
<a name="l00219"></a>00219 }
245
<a name="l00220"></a>00220 }
246
<a name="l00221"></a>00221 }
247
247
<a name="l00222"></a>00222
248
<a name="l00223"></a>00223 <span class="comment">/* This class is a safeguard to make sure the singleton gets deleted</span>
249
<a name="l00224"></a>00224 <span class="comment"> * when the application exits</span>
250
<a name="l00225"></a>00225 <span class="comment"> */</span>
251
<a name="l00226"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html">00226</a> <span class="keyword">class </span><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html">MSP::CCS::GeoidLibraryCleaner</a>
252
<a name="l00227"></a>00227 {
253
<a name="l00228"></a>00228 <span class="keyword">public</span>:
254
<a name="l00229"></a>00229
255
<a name="l00230"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html#a81b8f378199b3b6ad0b4a0dddb30fdca">00230</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html#a81b8f378199b3b6ad0b4a0dddb30fdca">~GeoidLibraryCleaner</a>()
256
<a name="l00231"></a>00231 {
257
<a name="l00232"></a>00232 GeoidLibrary::deleteInstance();
258
<a name="l00233"></a>00233 }
259
<a name="l00234"></a>00234
260
<a name="l00235"></a>00235 } <a class="code" href="_geoid_library_8cpp.html#ad6b9344cfa44b6798b0959c9b97929c6">geoidLibraryCleanerInstance</a>;
261
<a name="l00236"></a>00236
262
<a name="l00237"></a>00237
263
<a name="l00238"></a>00238 <span class="comment">// Make this class a singleton, so the data files are only initialized once</span>
264
<a name="l00239"></a>00239 <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>* GeoidLibrary::instance = 0;
265
<a name="l00240"></a>00240 <span class="keywordtype">int</span> GeoidLibrary::instanceCount = 0;
266
<a name="l00241"></a>00241
267
<a name="l00242"></a>00242
268
<a name="l00243"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a0c0c5e06957b3159e0cdadb3f929edd8">00243</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>* GeoidLibrary::getInstance()
269
<a name="l00244"></a>00244 {
270
<a name="l00245"></a>00245 <span class="keywordflow">if</span>( instance == 0 )
271
<a name="l00246"></a>00246 instance = <span class="keyword">new</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>;
272
<a name="l00247"></a>00247
273
<a name="l00248"></a>00248 instanceCount++;
274
<a name="l00249"></a>00249
275
<a name="l00250"></a>00250 <span class="keywordflow">return</span> instance;
276
<a name="l00251"></a>00251 }
277
<a name="l00252"></a>00252
278
<a name="l00253"></a>00253
279
<a name="l00254"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a1bc750a31ff20908209cb97067102e7e">00254</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a1bc750a31ff20908209cb97067102e7e">GeoidLibrary::removeInstance</a>()
248
<a name="l00223"></a>00223
249
<a name="l00224"></a><a class="code" href="_geoid_library_8cpp.html#a86b4187eeab007da13367ff3ac238250">00224</a> <span class="keywordtype">size_t</span> <a class="code" href="_geoid_library_8cpp.html#a86b4187eeab007da13367ff3ac238250">readBinary</a>(
250
<a name="l00225"></a>00225 <span class="keywordtype">void</span> *buffer,
251
<a name="l00226"></a>00226 <span class="keywordtype">size_t</span> size,
252
<a name="l00227"></a>00227 <span class="keywordtype">size_t</span> count,
253
<a name="l00228"></a>00228 FILE *stream )
254
<a name="l00229"></a>00229 {
255
<a name="l00230"></a>00230 <span class="keywordtype">size_t</span> actual_count = fread(buffer, size, count, stream);
256
<a name="l00231"></a>00231
257
<a name="l00232"></a>00232 <span class="keywordflow">if</span>(ferror(stream) || actual_count < count )
258
<a name="l00233"></a>00233 {
259
<a name="l00234"></a>00234 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
260
<a name="l00235"></a>00235 strcpy( message, <span class="stringliteral">"Error reading binary file."</span> );
261
<a name="l00236"></a>00236 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
262
<a name="l00237"></a>00237 }
263
<a name="l00238"></a>00238
264
<a name="l00239"></a>00239 <span class="preprocessor">#ifdef LITTLE_ENDIAN</span>
265
<a name="l00240"></a>00240 <span class="preprocessor"></span> <a class="code" href="_geoid_library_8cpp.html#a7b114915125b927e278ce213af694a31">swapBytes</a>(buffer, size, count);
266
<a name="l00241"></a>00241 <span class="preprocessor">#endif</span>
267
<a name="l00242"></a>00242 <span class="preprocessor"></span>
268
<a name="l00243"></a>00243 <span class="keywordflow">return</span> actual_count;
269
<a name="l00244"></a>00244 }
270
<a name="l00245"></a>00245
271
<a name="l00246"></a>00246 <span class="comment">/************************************************************************/</span>
272
<a name="l00247"></a>00247 <span class="comment">/* FUNCTIONS </span>
273
<a name="l00248"></a>00248 <span class="comment"> *</span>
274
<a name="l00249"></a>00249 <span class="comment"> */</span>
275
<a name="l00250"></a>00250
276
<a name="l00251"></a>00251 <span class="comment">/* This class is a safeguard to make sure the singleton gets deleted</span>
277
<a name="l00252"></a>00252 <span class="comment"> * when the application exits</span>
278
<a name="l00253"></a>00253 <span class="comment"> */</span>
279
<a name="l00254"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html">00254</a> <span class="keyword">class </span><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html">MSP::CCS::GeoidLibraryCleaner</a>
280
280
<a name="l00255"></a>00255 {
281
<a name="l00256"></a>00256 <span class="comment">/*</span>
282
<a name="l00257"></a>00257 <span class="comment"> * The function removeInstance removes this GeoidLibrary instance from the</span>
283
<a name="l00258"></a>00258 <span class="comment"> * total number of instances. </span>
284
<a name="l00259"></a>00259 <span class="comment"> */</span>
285
<a name="l00260"></a>00260 <span class="keywordflow">if</span>( --instanceCount < 1 )
286
<a name="l00261"></a>00261 {
287
<a name="l00262"></a>00262 deleteInstance();
288
<a name="l00263"></a>00263 }
289
<a name="l00264"></a>00264 }
281
<a name="l00256"></a>00256 <span class="keyword">public</span>:
282
<a name="l00257"></a>00257
283
<a name="l00258"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html#a81b8f378199b3b6ad0b4a0dddb30fdca">00258</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library_cleaner.html#a81b8f378199b3b6ad0b4a0dddb30fdca">~GeoidLibraryCleaner</a>()
284
<a name="l00259"></a>00259 {
285
<a name="l00260"></a>00260 <a class="code" href="class_m_s_p_1_1_c_c_s_thread_lock.html">CCSThreadLock</a> lock(&GeoidLibrary::mutex);
286
<a name="l00261"></a>00261 GeoidLibrary::deleteInstance();
287
<a name="l00262"></a>00262 }
288
<a name="l00263"></a>00263
289
<a name="l00264"></a>00264 } <a class="code" href="_geoid_library_8cpp.html#ad6b9344cfa44b6798b0959c9b97929c6">geoidLibraryCleanerInstance</a>;
290
290
<a name="l00265"></a>00265
291
291
<a name="l00266"></a>00266
292
<a name="l00267"></a>00267 <span class="keywordtype">void</span> GeoidLibrary::deleteInstance()
293
<a name="l00268"></a>00268 {
294
<a name="l00269"></a>00269 <span class="comment">/*</span>
295
<a name="l00270"></a>00270 <span class="comment"> * Delete the singleton.</span>
296
<a name="l00271"></a>00271 <span class="comment"> */</span>
292
<a name="l00267"></a>00267 <span class="comment">// Make this class a singleton, so the data files are only initialized once</span>
293
<a name="l00268"></a>00268 <a class="code" href="class_m_s_p_1_1_c_c_s_thread_mutex.html">CCSThreadMutex</a> GeoidLibrary::mutex;
294
<a name="l00269"></a>00269 <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>* GeoidLibrary::instance = 0;
295
<a name="l00270"></a>00270 <span class="keywordtype">int</span> GeoidLibrary::instanceCount = 0;
296
<a name="l00271"></a>00271
297
297
<a name="l00272"></a>00272
298
<a name="l00273"></a>00273 <span class="keywordflow">if</span>( instance != 0 )
299
<a name="l00274"></a>00274 {
300
<a name="l00275"></a>00275 <span class="keyword">delete</span> instance;
301
<a name="l00276"></a>00276 instance = 0;
302
<a name="l00277"></a>00277 }
303
<a name="l00278"></a>00278 }
304
<a name="l00279"></a>00279
298
<a name="l00273"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a0c0c5e06957b3159e0cdadb3f929edd8">00273</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>* GeoidLibrary::getInstance()
299
<a name="l00274"></a>00274 {
300
<a name="l00275"></a>00275 <a class="code" href="class_m_s_p_1_1_c_c_s_thread_lock.html">CCSThreadLock</a> lock(&mutex);
301
<a name="l00276"></a>00276 <span class="keywordflow">if</span>( instance == 0 )
302
<a name="l00277"></a>00277 instance = <span class="keyword">new</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>;
303
<a name="l00278"></a>00278
304
<a name="l00279"></a>00279 instanceCount++;
305
305
<a name="l00280"></a>00280
306
<a name="l00281"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a04ed88a0fbf5712724bf167e53461640">00281</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a04ed88a0fbf5712724bf167e53461640">GeoidLibrary::GeoidLibrary</a>()
307
<a name="l00282"></a>00282 {
308
<a name="l00283"></a>00283 <span class="comment">/*</span>
309
<a name="l00284"></a>00284 <span class="comment"> * The constructor creates an empty list which is used to store the geoid separation data</span>
310
<a name="l00285"></a>00285 <span class="comment"> * contained in the data files egm84.grd and egm96.grd</span>
311
<a name="l00286"></a>00286 <span class="comment"> */</span>
312
<a name="l00287"></a>00287
313
<a name="l00288"></a>00288 egm96GeoidList.reserve( <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a> );
314
<a name="l00289"></a>00289 egm84GeoidList.reserve( <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a> );
315
<a name="l00290"></a>00290
316
<a name="l00291"></a>00291 loadGeoids();
317
<a name="l00292"></a>00292 }
318
<a name="l00293"></a>00293
319
<a name="l00294"></a>00294
320
<a name="l00295"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ac1b892be471f2fa4b254d5eb6d0055af">00295</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a04ed88a0fbf5712724bf167e53461640">GeoidLibrary::GeoidLibrary</a>( <span class="keyword">const</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a> &gl )
321
<a name="l00296"></a>00296 {
322
<a name="l00297"></a>00297 egm96GeoidList = gl.egm96GeoidList;
323
<a name="l00298"></a>00298 egm84GeoidList = gl.egm84GeoidList;
324
<a name="l00299"></a>00299 }
325
<a name="l00300"></a>00300
326
<a name="l00301"></a>00301
327
<a name="l00302"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ae9391467ad025e2f8a7e35bd39d9e41d">00302</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ae9391467ad025e2f8a7e35bd39d9e41d">GeoidLibrary::~GeoidLibrary</a>()
328
<a name="l00303"></a>00303 {
329
<a name="l00304"></a>00304 egm96GeoidList.clear();
330
<a name="l00305"></a>00305 egm84GeoidList.clear();
331
<a name="l00306"></a>00306 }
332
<a name="l00307"></a>00307
333
<a name="l00308"></a>00308
334
<a name="l00309"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a49f4a11d8d5e91e55bb3221f92f0901d">00309</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>& <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a49f4a11d8d5e91e55bb3221f92f0901d">GeoidLibrary::operator=</a>( <span class="keyword">const</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a> &gl )
335
<a name="l00310"></a>00310 {
336
<a name="l00311"></a>00311 <span class="keywordflow">if</span> ( &gl == <span class="keyword">this</span> )
337
<a name="l00312"></a>00312 <span class="keywordflow">return</span> *<span class="keyword">this</span>;
338
<a name="l00313"></a>00313
339
<a name="l00314"></a>00314 egm96GeoidList = gl.egm96GeoidList;
340
<a name="l00315"></a>00315 egm84GeoidList = gl.egm84GeoidList;
341
<a name="l00316"></a>00316
342
<a name="l00317"></a>00317 <span class="keywordflow">return</span> *<span class="keyword">this</span>;
343
<a name="l00318"></a>00318 }
344
<a name="l00319"></a>00319
345
<a name="l00320"></a>00320
346
<a name="l00321"></a>00321 <span class="keywordtype">void</span> GeoidLibrary::loadGeoids()
347
<a name="l00322"></a>00322 {
348
<a name="l00323"></a>00323 <span class="comment">/*</span>
349
<a name="l00324"></a>00324 <span class="comment"> * The function loadGeoids reads geoid separation data from a file in</span>
350
<a name="l00325"></a>00325 <span class="comment"> * the current directory and builds the geoid separation table from it.</span>
351
<a name="l00326"></a>00326 <span class="comment"> * If the separation file can not be found or accessed, an error code of</span>
352
<a name="l00327"></a>00327 <span class="comment"> * GEOID_FILE_OPEN_ERROR is returned, If the separation file is incomplete</span>
353
<a name="l00328"></a>00328 <span class="comment"> * or improperly formatted, an error code of GEOID_INITIALIZE_ERROR is returned,</span>
354
<a name="l00329"></a>00329 <span class="comment"> * otherwise GEOID_NO_ERROR is returned. This function must be called before </span>
355
<a name="l00330"></a>00330 <span class="comment"> * any of the other functions in this component.</span>
356
<a name="l00331"></a>00331 <span class="comment"> */</span>
306
<a name="l00281"></a>00281 <span class="keywordflow">return</span> instance;
307
<a name="l00282"></a>00282 }
308
<a name="l00283"></a>00283
309
<a name="l00284"></a>00284
310
<a name="l00285"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a1bc750a31ff20908209cb97067102e7e">00285</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a1bc750a31ff20908209cb97067102e7e">GeoidLibrary::removeInstance</a>()
311
<a name="l00286"></a>00286 {
312
<a name="l00287"></a>00287 <span class="comment">/*</span>
313
<a name="l00288"></a>00288 <span class="comment"> * The function removeInstance removes this GeoidLibrary instance from the</span>
314
<a name="l00289"></a>00289 <span class="comment"> * total number of instances. </span>
315
<a name="l00290"></a>00290 <span class="comment"> */</span>
316
<a name="l00291"></a>00291 <a class="code" href="class_m_s_p_1_1_c_c_s_thread_lock.html">CCSThreadLock</a> lock(&mutex);
317
<a name="l00292"></a>00292 <span class="keywordflow">if</span> ( --instanceCount < 1 )
318
<a name="l00293"></a>00293 {
319
<a name="l00294"></a>00294 deleteInstance();
320
<a name="l00295"></a>00295 }
321
<a name="l00296"></a>00296 }
322
<a name="l00297"></a>00297
323
<a name="l00298"></a>00298
324
<a name="l00299"></a>00299 <span class="keywordtype">void</span> GeoidLibrary::deleteInstance()
325
<a name="l00300"></a>00300 {
326
<a name="l00301"></a>00301 <span class="comment">/*</span>
327
<a name="l00302"></a>00302 <span class="comment"> * Delete the singleton.</span>
328
<a name="l00303"></a>00303 <span class="comment"> */</span>
329
<a name="l00304"></a>00304
330
<a name="l00305"></a>00305 <span class="keywordflow">if</span>( instance != 0 )
331
<a name="l00306"></a>00306 {
332
<a name="l00307"></a>00307 <span class="keyword">delete</span> instance;
333
<a name="l00308"></a>00308 instance = 0;
334
<a name="l00309"></a>00309 }
335
<a name="l00310"></a>00310 }
336
<a name="l00311"></a>00311
337
<a name="l00312"></a>00312
338
<a name="l00313"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a04ed88a0fbf5712724bf167e53461640">00313</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a04ed88a0fbf5712724bf167e53461640">GeoidLibrary::GeoidLibrary</a>()
339
<a name="l00314"></a>00314 {
340
<a name="l00315"></a>00315 loadGeoids();
341
<a name="l00316"></a>00316 }
342
<a name="l00317"></a>00317
343
<a name="l00318"></a>00318
344
<a name="l00319"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ac1b892be471f2fa4b254d5eb6d0055af">00319</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a04ed88a0fbf5712724bf167e53461640">GeoidLibrary::GeoidLibrary</a>( <span class="keyword">const</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a> &gl )
345
<a name="l00320"></a>00320 {
346
<a name="l00321"></a>00321 *<span class="keyword">this</span> = gl;
347
<a name="l00322"></a>00322 }
348
<a name="l00323"></a>00323
349
<a name="l00324"></a>00324
350
<a name="l00325"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ae9391467ad025e2f8a7e35bd39d9e41d">00325</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ae9391467ad025e2f8a7e35bd39d9e41d">GeoidLibrary::~GeoidLibrary</a>()
351
<a name="l00326"></a>00326 {
352
<a name="l00327"></a>00327 <span class="keyword">delete</span> [] egm96GeoidList;
353
<a name="l00328"></a>00328 <span class="keyword">delete</span> [] egm84GeoidList;
354
<a name="l00329"></a>00329 <span class="keyword">delete</span> [] egm84ThirtyMinGeoidList;
355
<a name="l00330"></a>00330 }
356
<a name="l00331"></a>00331
357
357
<a name="l00332"></a>00332
358
<a name="l00333"></a>00333 initializeEGM96Geoid();
359
<a name="l00334"></a>00334
360
<a name="l00335"></a>00335 initializeEGM84Geoid();
361
<a name="l00336"></a>00336 }
358
<a name="l00333"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a49f4a11d8d5e91e55bb3221f92f0901d">00333</a> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a>& <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a49f4a11d8d5e91e55bb3221f92f0901d">GeoidLibrary::operator=</a>( <span class="keyword">const</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html">GeoidLibrary</a> &gl )
359
<a name="l00334"></a>00334 {
360
<a name="l00335"></a>00335 <span class="keywordflow">if</span> ( &gl == <span class="keyword">this</span> )
361
<a name="l00336"></a>00336 <span class="keywordflow">return</span> *<span class="keyword">this</span>;
362
362
<a name="l00337"></a>00337
363
<a name="l00338"></a>00338
364
<a name="l00339"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#abc13e044252f3e6a033fd1596fff50a8">00339</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#abc13e044252f3e6a033fd1596fff50a8">GeoidLibrary::convertEllipsoidToEGM96FifteenMinBilinearGeoidHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> ellipsoidHeight, <span class="keywordtype">double</span> *geoidHeight )
365
<a name="l00340"></a>00340 {
366
<a name="l00341"></a>00341 <span class="comment">/*</span>
367
<a name="l00342"></a>00342 <span class="comment"> * The function convertEllipsoidToEGM96FifteenMinBilinearGeoidHeight converts the specified WGS84</span>
368
<a name="l00343"></a>00343 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
369
<a name="l00344"></a>00344 <span class="comment"> * geoid height, using the EGM96 gravity model and the bilinear interpolation method.</span>
370
<a name="l00345"></a>00345 <span class="comment"> *</span>
371
<a name="l00346"></a>00346 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
372
<a name="l00347"></a>00347 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
373
<a name="l00348"></a>00348 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
374
<a name="l00349"></a>00349 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
375
<a name="l00350"></a>00350 <span class="comment"> *</span>
376
<a name="l00351"></a>00351 <span class="comment"> */</span>
377
<a name="l00352"></a>00352
378
<a name="l00353"></a>00353 <span class="keywordtype">double</span> delta_height;
379
<a name="l00354"></a>00354
380
<a name="l00355"></a>00355 bilinearInterpolate( longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>, <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>, egm96GeoidList, &delta_height );
381
<a name="l00356"></a>00356 *geoidHeight = ellipsoidHeight - delta_height;
382
<a name="l00357"></a>00357 }
383
<a name="l00358"></a>00358
363
<a name="l00338"></a>00338 egm96GeoidList = <span class="keyword">new</span> <span class="keywordtype">float</span>[<a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>];
364
<a name="l00339"></a>00339 egm84GeoidList = <span class="keyword">new</span> <span class="keywordtype">float</span>[<a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>];
365
<a name="l00340"></a>00340 egm84ThirtyMinGeoidList = <span class="keyword">new</span> <span class="keywordtype">double</span>[<a class="code" href="_geoid_library_8cpp.html#aa3581e7c236f1d903b4b3ccf242d1088">EGM84_30_MIN_ELEVATIONS</a>];
366
<a name="l00341"></a>00341
367
<a name="l00342"></a>00342 <span class="keywordflow">for</span>( <span class="keywordtype">int</span> i = 0; i < <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>; i++ )
368
<a name="l00343"></a>00343 {
369
<a name="l00344"></a>00344 egm96GeoidList[i] = gl.egm96GeoidList[i];
370
<a name="l00345"></a>00345 }
371
<a name="l00346"></a>00346
372
<a name="l00347"></a>00347 <span class="keywordflow">for</span>( <span class="keywordtype">int</span> j = 0; j < <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>; j++ )
373
<a name="l00348"></a>00348 {
374
<a name="l00349"></a>00349 egm84GeoidList[j] = gl.egm84GeoidList[j];
375
<a name="l00350"></a>00350 }
376
<a name="l00351"></a>00351
377
<a name="l00352"></a>00352 <span class="keywordflow">for</span>( <span class="keywordtype">int</span> k = 0; k < <a class="code" href="_geoid_library_8cpp.html#aa3581e7c236f1d903b4b3ccf242d1088">EGM84_30_MIN_ELEVATIONS</a>; k++ )
378
<a name="l00353"></a>00353 {
379
<a name="l00354"></a>00354 egm84ThirtyMinGeoidList[k] = gl.egm84ThirtyMinGeoidList[k];
380
<a name="l00355"></a>00355 }
381
<a name="l00356"></a>00356
382
<a name="l00357"></a>00357 <span class="keywordflow">return</span> *<span class="keyword">this</span>;
383
<a name="l00358"></a>00358 }
384
384
<a name="l00359"></a>00359
385
<a name="l00360"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ada3967bed409cf6c185849da0fba7464">00360</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ada3967bed409cf6c185849da0fba7464">GeoidLibrary::convertEllipsoidToEGM96VariableNaturalSplineHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> ellipsoidHeight, <span class="keywordtype">double</span> *geoidHeight )
386
<a name="l00361"></a>00361 {
387
<a name="l00362"></a>00362 <span class="comment">/*</span>
388
<a name="l00363"></a>00363 <span class="comment"> * The function convertEllipsoidToEGM96VariableNaturalSplineHeight converts the specified WGS84</span>
389
<a name="l00364"></a>00364 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
390
<a name="l00365"></a>00365 <span class="comment"> * geoid height, using the EGM96 gravity model and the natural spline interpolation method..</span>
391
<a name="l00366"></a>00366 <span class="comment"> *</span>
392
<a name="l00367"></a>00367 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
393
<a name="l00368"></a>00368 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
394
<a name="l00369"></a>00369 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
395
<a name="l00370"></a>00370 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
396
<a name="l00371"></a>00371 <span class="comment"> *</span>
397
<a name="l00372"></a>00372 <span class="comment"> */</span>
398
<a name="l00373"></a>00373
399
<a name="l00374"></a>00374 <span class="keywordtype">int</span> i = 0;
400
<a name="l00375"></a>00375 <span class="keywordtype">int</span> num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
401
<a name="l00376"></a>00376 <span class="keywordtype">int</span> num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
402
<a name="l00377"></a>00377 <span class="keywordtype">double</span> latitude_degrees = latitude * <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a>;
403
<a name="l00378"></a>00378 <span class="keywordtype">double</span> longitude_degrees = longitude * _180_OVER_PI;
404
<a name="l00379"></a>00379 <span class="keywordtype">double</span> scale_factor = <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>;
405
<a name="l00380"></a>00380 <span class="keywordtype">double</span> delta_height;
406
<a name="l00381"></a>00381 <span class="keywordtype">long</span> found = 0;
407
<a name="l00382"></a>00382
408
<a name="l00383"></a>00383 <span class="keywordflow">if</span>( longitude_degrees < 0.0 )
409
<a name="l00384"></a>00384 longitude_degrees += 360.0;
410
<a name="l00385"></a>00385
411
<a name="l00386"></a>00386 <span class="keywordflow">while</span>( !found && i < <a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a> )
412
<a name="l00387"></a>00387 {
413
<a name="l00388"></a>00388 <span class="keywordflow">if</span>( ( latitude_degrees >= EGM96_Variable_Grid_Table[i].min_lat ) && ( longitude_degrees >= EGM96_Variable_Grid_Table[i].min_lon ) &&
414
<a name="l00389"></a>00389 ( latitude_degrees < EGM96_Variable_Grid_Table[i].max_lat ) && ( longitude_degrees < EGM96_Variable_Grid_Table[i].max_lon ) )
415
<a name="l00390"></a>00390 {
416
<a name="l00391"></a>00391 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a>; <span class="comment">// use 30 minute by 30 minute grid</span>
417
<a name="l00392"></a>00392 num_cols = 721;
418
<a name="l00393"></a>00393 num_rows = 361;
419
<a name="l00394"></a>00394 found = 1;
420
<a name="l00395"></a>00395 }
385
<a name="l00360"></a>00360
386
<a name="l00361"></a>00361 <span class="keywordtype">void</span> GeoidLibrary::loadGeoids()
387
<a name="l00362"></a>00362 {
388
<a name="l00363"></a>00363 <span class="comment">/*</span>
389
<a name="l00364"></a>00364 <span class="comment"> * The function loadGeoids reads geoid separation data from a file in</span>
390
<a name="l00365"></a>00365 <span class="comment"> * the current directory and builds the geoid separation table from it.</span>
391
<a name="l00366"></a>00366 <span class="comment"> * If the separation file can not be found or accessed, an error code of</span>
392
<a name="l00367"></a>00367 <span class="comment"> * GEOID_FILE_OPEN_ERROR is returned, If the separation file is incomplete</span>
393
<a name="l00368"></a>00368 <span class="comment"> * or improperly formatted, an error code of GEOID_INITIALIZE_ERROR is returned,</span>
394
<a name="l00369"></a>00369 <span class="comment"> * otherwise GEOID_NO_ERROR is returned. This function must be called before </span>
395
<a name="l00370"></a>00370 <span class="comment"> * any of the other functions in this component.</span>
396
<a name="l00371"></a>00371 <span class="comment"> */</span>
397
<a name="l00372"></a>00372
398
<a name="l00373"></a>00373 initializeEGM96Geoid();
399
<a name="l00374"></a>00374 initializeEGM84Geoid();
400
<a name="l00375"></a>00375 initializeEGM84ThirtyMinGeoid();
401
<a name="l00376"></a>00376 }
402
<a name="l00377"></a>00377
403
<a name="l00378"></a>00378
404
<a name="l00379"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#abc13e044252f3e6a033fd1596fff50a8">00379</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#abc13e044252f3e6a033fd1596fff50a8">GeoidLibrary::convertEllipsoidToEGM96FifteenMinBilinearGeoidHeight</a>(
405
<a name="l00380"></a>00380 <span class="keywordtype">double</span> longitude,
406
<a name="l00381"></a>00381 <span class="keywordtype">double</span> latitude,
407
<a name="l00382"></a>00382 <span class="keywordtype">double</span> ellipsoidHeight,
408
<a name="l00383"></a>00383 <span class="keywordtype">double</span> *geoidHeight )
409
<a name="l00384"></a>00384 {
410
<a name="l00385"></a>00385 <span class="comment">/*</span>
411
<a name="l00386"></a>00386 <span class="comment"> * The function convertEllipsoidToEGM96FifteenMinBilinearGeoidHeight converts the specified WGS84</span>
412
<a name="l00387"></a>00387 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
413
<a name="l00388"></a>00388 <span class="comment"> * geoid height, using the EGM96 gravity model and the bilinear interpolation method.</span>
414
<a name="l00389"></a>00389 <span class="comment"> *</span>
415
<a name="l00390"></a>00390 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
416
<a name="l00391"></a>00391 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
417
<a name="l00392"></a>00392 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
418
<a name="l00393"></a>00393 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
419
<a name="l00394"></a>00394 <span class="comment"> *</span>
420
<a name="l00395"></a>00395 <span class="comment"> */</span>
421
421
<a name="l00396"></a>00396
422
<a name="l00397"></a>00397 i++;
423
<a name="l00398"></a>00398 }
424
<a name="l00399"></a>00399
425
<a name="l00400"></a>00400 <span class="keywordflow">if</span>( !found )
426
<a name="l00401"></a>00401 {
427
<a name="l00402"></a>00402 <span class="keywordflow">if</span>( latitude_degrees >= -60.0 && latitude_degrees < 60.0 )
428
<a name="l00403"></a>00403 {
429
<a name="l00404"></a>00404 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a>; <span class="comment">// use 1 degree by 1 degree grid</span>
430
<a name="l00405"></a>00405 num_cols = 361;
431
<a name="l00406"></a>00406 num_rows = 181;
432
<a name="l00407"></a>00407 }
433
<a name="l00408"></a>00408 <span class="keywordflow">else</span>
434
<a name="l00409"></a>00409 {
435
<a name="l00410"></a>00410 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a>; <span class="comment">// use 2 degree by 2 degree grid</span>
436
<a name="l00411"></a>00411 num_cols = 181;
437
<a name="l00412"></a>00412 num_rows = 91;
438
<a name="l00413"></a>00413 }
439
<a name="l00414"></a>00414 }
440
<a name="l00415"></a>00415
441
<a name="l00416"></a>00416 naturalSplineInterpolate( longitude, latitude, scale_factor, num_cols, num_rows, <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>-1, egm96GeoidList, &delta_height );
442
<a name="l00417"></a>00417 *geoidHeight = ellipsoidHeight - delta_height;
443
<a name="l00418"></a>00418 }
444
<a name="l00419"></a>00419
445
<a name="l00420"></a>00420
446
<a name="l00421"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a75723dea06d1708f8be9670e68fdfa4f">00421</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a75723dea06d1708f8be9670e68fdfa4f">GeoidLibrary::convertEllipsoidToEGM84TenDegBilinearHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> ellipsoidHeight, <span class="keywordtype">double</span> *geoidHeight )
447
<a name="l00422"></a>00422 {
448
<a name="l00423"></a>00423 <span class="comment">/*</span>
449
<a name="l00424"></a>00424 <span class="comment"> * The function convertEllipsoidToEGM84TenDegBilinearHeight converts the specified WGS84</span>
450
<a name="l00425"></a>00425 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
451
<a name="l00426"></a>00426 <span class="comment"> * geoid height, using the EGM84 gravity model and the bilinear interpolation method..</span>
452
<a name="l00427"></a>00427 <span class="comment"> *</span>
453
<a name="l00428"></a>00428 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
454
<a name="l00429"></a>00429 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
455
<a name="l00430"></a>00430 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
456
<a name="l00431"></a>00431 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
457
<a name="l00432"></a>00432 <span class="comment"> *</span>
458
<a name="l00433"></a>00433 <span class="comment"> */</span>
422
<a name="l00397"></a>00397 <span class="keywordtype">double</span> delta_height;
423
<a name="l00398"></a>00398
424
<a name="l00399"></a>00399 bilinearInterpolate(
425
<a name="l00400"></a>00400 longitude, latitude,
426
<a name="l00401"></a>00401 <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>, <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>, egm96GeoidList,
427
<a name="l00402"></a>00402 &delta_height );
428
<a name="l00403"></a>00403
429
<a name="l00404"></a>00404 *geoidHeight = ellipsoidHeight - delta_height;
430
<a name="l00405"></a>00405 }
431
<a name="l00406"></a>00406
432
<a name="l00407"></a>00407
433
<a name="l00408"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ada3967bed409cf6c185849da0fba7464">00408</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#ada3967bed409cf6c185849da0fba7464">GeoidLibrary::convertEllipsoidToEGM96VariableNaturalSplineHeight</a>(
434
<a name="l00409"></a>00409 <span class="keywordtype">double</span> longitude,
435
<a name="l00410"></a>00410 <span class="keywordtype">double</span> latitude,
436
<a name="l00411"></a>00411 <span class="keywordtype">double</span> ellipsoidHeight,
437
<a name="l00412"></a>00412 <span class="keywordtype">double</span> *geoidHeight )
438
<a name="l00413"></a>00413 {
439
<a name="l00414"></a>00414 <span class="comment">/*</span>
440
<a name="l00415"></a>00415 <span class="comment"> * The function convertEllipsoidToEGM96VariableNaturalSplineHeight converts the specified WGS84</span>
441
<a name="l00416"></a>00416 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
442
<a name="l00417"></a>00417 <span class="comment"> * geoid height, using the EGM96 gravity model and the natural spline interpolation method..</span>
443
<a name="l00418"></a>00418 <span class="comment"> *</span>
444
<a name="l00419"></a>00419 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
445
<a name="l00420"></a>00420 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
446
<a name="l00421"></a>00421 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
447
<a name="l00422"></a>00422 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
448
<a name="l00423"></a>00423 <span class="comment"> *</span>
449
<a name="l00424"></a>00424 <span class="comment"> */</span>
450
<a name="l00425"></a>00425
451
<a name="l00426"></a>00426 <span class="keywordtype">int</span> i = 0;
452
<a name="l00427"></a>00427 <span class="keywordtype">int</span> num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
453
<a name="l00428"></a>00428 <span class="keywordtype">int</span> num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
454
<a name="l00429"></a>00429 <span class="keywordtype">double</span> latitude_degrees = latitude * <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a>;
455
<a name="l00430"></a>00430 <span class="keywordtype">double</span> longitude_degrees = longitude * _180_OVER_PI;
456
<a name="l00431"></a>00431 <span class="keywordtype">double</span> scale_factor = <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>;
457
<a name="l00432"></a>00432 <span class="keywordtype">double</span> delta_height;
458
<a name="l00433"></a>00433 <span class="keywordtype">long</span> found = 0;
459
459
<a name="l00434"></a>00434
460
<a name="l00435"></a>00435 <span class="keywordtype">double</span> delta_height;
461
<a name="l00436"></a>00436
462
<a name="l00437"></a>00437 bilinearInterpolate( longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>, <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, egm84GeoidList, &delta_height );
463
<a name="l00438"></a>00438 *geoidHeight = ellipsoidHeight - delta_height;
464
<a name="l00439"></a>00439 }
465
<a name="l00440"></a>00440
466
<a name="l00441"></a>00441
467
<a name="l00442"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a655d40e5be6d0373ac5c473f745eaa16">00442</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a655d40e5be6d0373ac5c473f745eaa16">GeoidLibrary::convertEllipsoidToEGM84TenDegNaturalSplineHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> ellipsoidHeight, <span class="keywordtype">double</span> *geoidHeight )
468
<a name="l00443"></a>00443 {
469
<a name="l00444"></a>00444 <span class="comment">/*</span>
470
<a name="l00445"></a>00445 <span class="comment"> * The function convertEllipsoidToEGM84TenDegNaturalSplineHeight converts the specified WGS84</span>
471
<a name="l00446"></a>00446 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
472
<a name="l00447"></a>00447 <span class="comment"> * geoid height, using the EGM84 gravity model and the natural spline interpolation method..</span>
473
<a name="l00448"></a>00448 <span class="comment"> *</span>
474
<a name="l00449"></a>00449 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
475
<a name="l00450"></a>00450 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
476
<a name="l00451"></a>00451 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
477
<a name="l00452"></a>00452 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
478
<a name="l00453"></a>00453 <span class="comment"> *</span>
479
<a name="l00454"></a>00454 <span class="comment"> */</span>
480
<a name="l00455"></a>00455
481
<a name="l00456"></a>00456 <span class="keywordtype">double</span> delta_height;
482
<a name="l00457"></a>00457
483
<a name="l00458"></a>00458 naturalSplineInterpolate( longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>, <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>-1, egm84GeoidList, &delta_height );
484
<a name="l00459"></a>00459 *geoidHeight = ellipsoidHeight - delta_height;
485
<a name="l00460"></a>00460 }
486
<a name="l00461"></a>00461
487
<a name="l00462"></a>00462
488
<a name="l00463"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a319274ea37f279fab53f78ce8a609493">00463</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a319274ea37f279fab53f78ce8a609493">GeoidLibrary::convertEGM96FifteenMinBilinearGeoidToEllipsoidHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> geoidHeight, <span class="keywordtype">double</span> *ellipsoidHeight )
489
<a name="l00464"></a>00464 {
490
<a name="l00465"></a>00465 <span class="comment">/*</span>
491
<a name="l00466"></a>00466 <span class="comment"> * The function convertEGM96FifteenMinBilinearGeoidToEllipsoidHeight converts the specified WGS84</span>
492
<a name="l00467"></a>00467 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
493
<a name="l00468"></a>00468 <span class="comment"> * ellipsoid height, using the EGM96 gravity model and the bilinear interpolation method.</span>
494
<a name="l00469"></a>00469 <span class="comment"> *</span>
495
<a name="l00470"></a>00470 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
496
<a name="l00471"></a>00471 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
497
<a name="l00472"></a>00472 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
498
<a name="l00473"></a>00473 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
499
<a name="l00474"></a>00474 <span class="comment"> *</span>
500
<a name="l00475"></a>00475 <span class="comment"> */</span>
501
<a name="l00476"></a>00476
502
<a name="l00477"></a>00477 <span class="keywordtype">double</span> delta_height;
460
<a name="l00435"></a>00435 <span class="keywordflow">if</span>( longitude_degrees < 0.0 )
461
<a name="l00436"></a>00436 longitude_degrees += 360.0;
462
<a name="l00437"></a>00437
463
<a name="l00438"></a>00438 <span class="keywordflow">while</span>( !found && i < <a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a> )
464
<a name="l00439"></a>00439 {
465
<a name="l00440"></a>00440 <span class="keywordflow">if</span>(( latitude_degrees >= EGM96_Variable_Grid_Table[i].min_lat ) &&
466
<a name="l00441"></a>00441 ( longitude_degrees >= EGM96_Variable_Grid_Table[i].min_lon ) &&
467
<a name="l00442"></a>00442 ( latitude_degrees < EGM96_Variable_Grid_Table[i].max_lat ) &&
468
<a name="l00443"></a>00443 ( longitude_degrees < EGM96_Variable_Grid_Table[i].max_lon ) )
469
<a name="l00444"></a>00444 {
470
<a name="l00445"></a>00445 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a>; <span class="comment">// use 30 minute by 30 minute grid</span>
471
<a name="l00446"></a>00446 num_cols = 721;
472
<a name="l00447"></a>00447 num_rows = 361;
473
<a name="l00448"></a>00448 found = 1;
474
<a name="l00449"></a>00449 }
475
<a name="l00450"></a>00450
476
<a name="l00451"></a>00451 i++;
477
<a name="l00452"></a>00452 }
478
<a name="l00453"></a>00453
479
<a name="l00454"></a>00454 <span class="keywordflow">if</span>( !found )
480
<a name="l00455"></a>00455 {
481
<a name="l00456"></a>00456 <span class="keywordflow">if</span>( latitude_degrees >= -60.0 && latitude_degrees < 60.0 )
482
<a name="l00457"></a>00457 {
483
<a name="l00458"></a>00458 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a>; <span class="comment">// use 1 degree by 1 degree grid</span>
484
<a name="l00459"></a>00459 num_cols = 361;
485
<a name="l00460"></a>00460 num_rows = 181;
486
<a name="l00461"></a>00461 }
487
<a name="l00462"></a>00462 <span class="keywordflow">else</span>
488
<a name="l00463"></a>00463 {
489
<a name="l00464"></a>00464 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a>; <span class="comment">// use 2 degree by 2 degree grid</span>
490
<a name="l00465"></a>00465 num_cols = 181;
491
<a name="l00466"></a>00466 num_rows = 91;
492
<a name="l00467"></a>00467 }
493
<a name="l00468"></a>00468 }
494
<a name="l00469"></a>00469
495
<a name="l00470"></a>00470 naturalSplineInterpolate(
496
<a name="l00471"></a>00471 longitude, latitude,
497
<a name="l00472"></a>00472 scale_factor, num_cols, num_rows, <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>-1, egm96GeoidList,
498
<a name="l00473"></a>00473 &delta_height );
499
<a name="l00474"></a>00474
500
<a name="l00475"></a>00475 *geoidHeight = ellipsoidHeight - delta_height;
501
<a name="l00476"></a>00476 }
502
<a name="l00477"></a>00477
503
503
<a name="l00478"></a>00478
504
<a name="l00479"></a>00479 bilinearInterpolate( longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>, <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>, egm96GeoidList, &delta_height );
505
<a name="l00480"></a>00480 *ellipsoidHeight = geoidHeight + delta_height;
506
<a name="l00481"></a>00481 }
507
<a name="l00482"></a>00482
508
<a name="l00483"></a>00483
509
<a name="l00484"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#aa3f8de2392ade95991fd1519a7f869fa">00484</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#aa3f8de2392ade95991fd1519a7f869fa">GeoidLibrary::convertEGM96VariableNaturalSplineToEllipsoidHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> geoidHeight, <span class="keywordtype">double</span> *ellipsoidHeight )
510
<a name="l00485"></a>00485 {
511
<a name="l00486"></a>00486 <span class="comment">/*</span>
512
<a name="l00487"></a>00487 <span class="comment"> * The function convertEGM96VariableNaturalSplineToEllipsoidHeight converts the specified WGS84</span>
513
<a name="l00488"></a>00488 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
514
<a name="l00489"></a>00489 <span class="comment"> * ellipsoid height, using the EGM96 gravity model and the natural spline interpolation method.</span>
515
<a name="l00490"></a>00490 <span class="comment"> *</span>
516
<a name="l00491"></a>00491 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
517
<a name="l00492"></a>00492 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
518
<a name="l00493"></a>00493 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
519
<a name="l00494"></a>00494 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
520
<a name="l00495"></a>00495 <span class="comment"> *</span>
521
<a name="l00496"></a>00496 <span class="comment"> */</span>
522
<a name="l00497"></a>00497
523
<a name="l00498"></a>00498 <span class="keywordtype">int</span> i = 0;
524
<a name="l00499"></a>00499 <span class="keywordtype">int</span> num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
525
<a name="l00500"></a>00500 <span class="keywordtype">int</span> num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
526
<a name="l00501"></a>00501 <span class="keywordtype">double</span> latitude_degrees = latitude * <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a>;
527
<a name="l00502"></a>00502 <span class="keywordtype">double</span> longitude_degrees = longitude * _180_OVER_PI;
528
<a name="l00503"></a>00503 <span class="keywordtype">double</span> scale_factor = <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>;
529
<a name="l00504"></a>00504 <span class="keywordtype">double</span> delta_height;
530
<a name="l00505"></a>00505 <span class="keywordtype">long</span> found = 0;
504
<a name="l00479"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a75723dea06d1708f8be9670e68fdfa4f">00479</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a75723dea06d1708f8be9670e68fdfa4f">GeoidLibrary::convertEllipsoidToEGM84TenDegBilinearHeight</a>(
505
<a name="l00480"></a>00480 <span class="keywordtype">double</span> longitude,
506
<a name="l00481"></a>00481 <span class="keywordtype">double</span> latitude,
507
<a name="l00482"></a>00482 <span class="keywordtype">double</span> ellipsoidHeight,
508
<a name="l00483"></a>00483 <span class="keywordtype">double</span> *geoidHeight )
509
<a name="l00484"></a>00484 {
510
<a name="l00485"></a>00485 <span class="comment">/*</span>
511
<a name="l00486"></a>00486 <span class="comment"> * The function convertEllipsoidToEGM84TenDegBilinearHeight converts the specified WGS84</span>
512
<a name="l00487"></a>00487 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
513
<a name="l00488"></a>00488 <span class="comment"> * geoid height, using the EGM84 gravity model and the bilinear interpolation method.</span>
514
<a name="l00489"></a>00489 <span class="comment"> *</span>
515
<a name="l00490"></a>00490 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
516
<a name="l00491"></a>00491 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
517
<a name="l00492"></a>00492 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
518
<a name="l00493"></a>00493 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
519
<a name="l00494"></a>00494 <span class="comment"> *</span>
520
<a name="l00495"></a>00495 <span class="comment"> */</span>
521
<a name="l00496"></a>00496
522
<a name="l00497"></a>00497 <span class="keywordtype">double</span> delta_height;
523
<a name="l00498"></a>00498
524
<a name="l00499"></a>00499 bilinearInterpolate(
525
<a name="l00500"></a>00500 longitude, latitude,
526
<a name="l00501"></a>00501 <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>, <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, egm84GeoidList,
527
<a name="l00502"></a>00502 &delta_height );
528
<a name="l00503"></a>00503
529
<a name="l00504"></a>00504 *geoidHeight = ellipsoidHeight - delta_height;
530
<a name="l00505"></a>00505 }
531
531
<a name="l00506"></a>00506
532
<a name="l00507"></a>00507 <span class="keywordflow">if</span>( longitude_degrees < 0.0 )
533
<a name="l00508"></a>00508 longitude_degrees += 360.0;
534
<a name="l00509"></a>00509
535
<a name="l00510"></a>00510 <span class="keywordflow">while</span>( !found && i < <a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a> )
536
<a name="l00511"></a>00511 {
537
<a name="l00512"></a>00512 <span class="keywordflow">if</span>( ( latitude_degrees >= EGM96_Variable_Grid_Table[i].min_lat ) && ( longitude_degrees >= EGM96_Variable_Grid_Table[i].min_lon ) &&
538
<a name="l00513"></a>00513 ( latitude_degrees < EGM96_Variable_Grid_Table[i].max_lat ) && ( longitude_degrees < EGM96_Variable_Grid_Table[i].max_lon ) )
539
<a name="l00514"></a>00514 {
540
<a name="l00515"></a>00515 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a>; <span class="comment">// use 30 minute by 30 minute grid</span>
541
<a name="l00516"></a>00516 num_cols = 721;
542
<a name="l00517"></a>00517 num_rows = 361;
543
<a name="l00518"></a>00518 found = 1;
544
<a name="l00519"></a>00519 }
545
<a name="l00520"></a>00520
546
<a name="l00521"></a>00521 i++;
547
<a name="l00522"></a>00522 }
548
<a name="l00523"></a>00523
549
<a name="l00524"></a>00524 <span class="keywordflow">if</span>( !found )
550
<a name="l00525"></a>00525 {
551
<a name="l00526"></a>00526 <span class="keywordflow">if</span>( latitude_degrees >= -60.0 && latitude_degrees < 60.0 )
552
<a name="l00527"></a>00527 {
553
<a name="l00528"></a>00528 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a>; <span class="comment">// use 1 degree by 1 degree grid</span>
554
<a name="l00529"></a>00529 num_cols = 361;
555
<a name="l00530"></a>00530 num_rows = 181;
556
<a name="l00531"></a>00531 }
557
<a name="l00532"></a>00532 <span class="keywordflow">else</span>
558
<a name="l00533"></a>00533 {
559
<a name="l00534"></a>00534 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a>; <span class="comment">// use 2 degree by 2 degree grid</span>
560
<a name="l00535"></a>00535 num_cols = 181;
561
<a name="l00536"></a>00536 num_rows = 91;
562
<a name="l00537"></a>00537 }
563
<a name="l00538"></a>00538 }
564
<a name="l00539"></a>00539
565
<a name="l00540"></a>00540 naturalSplineInterpolate( longitude, latitude, scale_factor, num_cols, num_rows, <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>-1, egm96GeoidList, &delta_height );
566
<a name="l00541"></a>00541 *ellipsoidHeight = geoidHeight + delta_height;
567
<a name="l00542"></a>00542 }
568
<a name="l00543"></a>00543
569
<a name="l00544"></a>00544
570
<a name="l00545"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a68378c47cf83b97e59db68c7f287c4bc">00545</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a68378c47cf83b97e59db68c7f287c4bc">GeoidLibrary::convertEGM84TenDegBilinearToEllipsoidHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> geoidHeight, <span class="keywordtype">double</span> *ellipsoidHeight )
571
<a name="l00546"></a>00546 {
572
<a name="l00547"></a>00547 <span class="comment">/*</span>
573
<a name="l00548"></a>00548 <span class="comment"> * The function convertEGM84TenDegBilinearToEllipsoidHeight converts the specified WGS84</span>
574
<a name="l00549"></a>00549 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
575
<a name="l00550"></a>00550 <span class="comment"> * ellipsoid height, using the EGM84 gravity model and the bilinear interpolation method.</span>
532
<a name="l00507"></a>00507
533
<a name="l00508"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a655d40e5be6d0373ac5c473f745eaa16">00508</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a655d40e5be6d0373ac5c473f745eaa16">GeoidLibrary::convertEllipsoidToEGM84TenDegNaturalSplineHeight</a>(
534
<a name="l00509"></a>00509 <span class="keywordtype">double</span> longitude,
535
<a name="l00510"></a>00510 <span class="keywordtype">double</span> latitude,
536
<a name="l00511"></a>00511 <span class="keywordtype">double</span> ellipsoidHeight,
537
<a name="l00512"></a>00512 <span class="keywordtype">double</span> *geoidHeight )
538
<a name="l00513"></a>00513 {
539
<a name="l00514"></a>00514 <span class="comment">/*</span>
540
<a name="l00515"></a>00515 <span class="comment"> * The function convertEllipsoidToEGM84TenDegNaturalSplineHeight converts the specified WGS84</span>
541
<a name="l00516"></a>00516 <span class="comment"> * ellipsoid height at the specified geodetic coordinates to the equivalent</span>
542
<a name="l00517"></a>00517 <span class="comment"> * geoid height, using the EGM84 gravity model and the natural spline interpolation method..</span>
543
<a name="l00518"></a>00518 <span class="comment"> *</span>
544
<a name="l00519"></a>00519 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
545
<a name="l00520"></a>00520 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
546
<a name="l00521"></a>00521 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
547
<a name="l00522"></a>00522 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
548
<a name="l00523"></a>00523 <span class="comment"> *</span>
549
<a name="l00524"></a>00524 <span class="comment"> */</span>
550
<a name="l00525"></a>00525
551
<a name="l00526"></a>00526 <span class="keywordtype">double</span> delta_height;
552
<a name="l00527"></a>00527
553
<a name="l00528"></a>00528 naturalSplineInterpolate(
554
<a name="l00529"></a>00529 longitude, latitude,
555
<a name="l00530"></a>00530 <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>, <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>-1,
556
<a name="l00531"></a>00531 egm84GeoidList, &delta_height );
557
<a name="l00532"></a>00532
558
<a name="l00533"></a>00533 *geoidHeight = ellipsoidHeight - delta_height;
559
<a name="l00534"></a>00534 }
560
<a name="l00535"></a>00535
561
<a name="l00536"></a>00536
562
<a name="l00537"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#aae2680e03a7a00e9f046ee080cc63137">00537</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#aae2680e03a7a00e9f046ee080cc63137">GeoidLibrary::convertEllipsoidToEGM84ThirtyMinBiLinearHeight</a>(
563
<a name="l00538"></a>00538 <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> ellipsoidHeight,
564
<a name="l00539"></a>00539 <span class="keywordtype">double</span> *geoidHeight )
565
<a name="l00540"></a>00540 {
566
<a name="l00541"></a>00541 <span class="comment">/*</span>
567
<a name="l00542"></a>00542 <span class="comment"> * The function convertEllipsoidToEGM84ThirtyMinBiLinearHeight converts the </span>
568
<a name="l00543"></a>00543 <span class="comment"> * specified WGS84 ellipsoid height at the specified geodetic coordinates to the</span>
569
<a name="l00544"></a>00544 <span class="comment"> * equivalent geoid height, using the EGM84 gravity model and the bilinear </span>
570
<a name="l00545"></a>00545 <span class="comment"> * interpolation method..</span>
571
<a name="l00546"></a>00546 <span class="comment"> *</span>
572
<a name="l00547"></a>00547 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
573
<a name="l00548"></a>00548 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
574
<a name="l00549"></a>00549 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (input)</span>
575
<a name="l00550"></a>00550 <span class="comment"> * geoidHeight : Geoid height, in meters. (output)</span>
576
576
<a name="l00551"></a>00551 <span class="comment"> *</span>
577
<a name="l00552"></a>00552 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
578
<a name="l00553"></a>00553 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
579
<a name="l00554"></a>00554 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
580
<a name="l00555"></a>00555 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
581
<a name="l00556"></a>00556 <span class="comment"> *</span>
582
<a name="l00557"></a>00557 <span class="comment"> */</span>
583
<a name="l00558"></a>00558
584
<a name="l00559"></a>00559 <span class="keywordtype">double</span> delta_height;
577
<a name="l00552"></a>00552 <span class="comment"> */</span>
578
<a name="l00553"></a>00553
579
<a name="l00554"></a>00554 <span class="keywordtype">double</span> delta_height;
580
<a name="l00555"></a>00555
581
<a name="l00556"></a>00556 bilinearInterpolateDoubleHeights(
582
<a name="l00557"></a>00557 longitude, latitude,
583
<a name="l00558"></a>00558 <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a>, <a class="code" href="_geoid_library_8cpp.html#a580782fa058e200088161c1680921220">EGM84_30_MIN_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#a1ea8fba0f5bc5035bc21da731b6d8799">EGM84_30_MIN_ROWS</a>,
584
<a name="l00559"></a>00559 egm84ThirtyMinGeoidList, &delta_height );
585
585
<a name="l00560"></a>00560
586
<a name="l00561"></a>00561 bilinearInterpolate( longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>, <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, egm84GeoidList, &delta_height );
587
<a name="l00562"></a>00562 *ellipsoidHeight = geoidHeight + delta_height;
588
<a name="l00563"></a>00563 }
586
<a name="l00561"></a>00561 *geoidHeight = ellipsoidHeight - delta_height;
587
<a name="l00562"></a>00562 }
588
<a name="l00563"></a>00563
589
589
<a name="l00564"></a>00564
590
<a name="l00565"></a>00565
591
<a name="l00566"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a8f085fab24f9cd29bc1ad8d75c2ec798">00566</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a8f085fab24f9cd29bc1ad8d75c2ec798">GeoidLibrary::convertEGM84TenDegNaturalSplineToEllipsoidHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> geoidHeight, <span class="keywordtype">double</span> *ellipsoidHeight )
592
<a name="l00567"></a>00567 {
593
<a name="l00568"></a>00568 <span class="comment">/*</span>
594
<a name="l00569"></a>00569 <span class="comment"> * The function convertEGM84TenDegNaturalSplineToEllipsoidHeight converts the specified WGS84</span>
595
<a name="l00570"></a>00570 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
596
<a name="l00571"></a>00571 <span class="comment"> * ellipsoid height, using the EGM84 gravity model and the natural spline interpolation method.</span>
597
<a name="l00572"></a>00572 <span class="comment"> *</span>
598
<a name="l00573"></a>00573 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
599
<a name="l00574"></a>00574 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
600
<a name="l00575"></a>00575 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
601
<a name="l00576"></a>00576 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
602
<a name="l00577"></a>00577 <span class="comment"> *</span>
603
<a name="l00578"></a>00578 <span class="comment"> */</span>
604
<a name="l00579"></a>00579
605
<a name="l00580"></a>00580 <span class="keywordtype">double</span> delta_height;
606
<a name="l00581"></a>00581
607
<a name="l00582"></a>00582 naturalSplineInterpolate( longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>, <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>-1, egm84GeoidList, &delta_height );
608
<a name="l00583"></a>00583 *ellipsoidHeight = geoidHeight + delta_height;
609
<a name="l00584"></a>00584 }
590
<a name="l00565"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a319274ea37f279fab53f78ce8a609493">00565</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a319274ea37f279fab53f78ce8a609493">GeoidLibrary::convertEGM96FifteenMinBilinearGeoidToEllipsoidHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> geoidHeight, <span class="keywordtype">double</span> *ellipsoidHeight )
591
<a name="l00566"></a>00566 {
592
<a name="l00567"></a>00567 <span class="comment">/*</span>
593
<a name="l00568"></a>00568 <span class="comment"> * The function convertEGM96FifteenMinBilinearGeoidToEllipsoidHeight converts the specified WGS84</span>
594
<a name="l00569"></a>00569 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
595
<a name="l00570"></a>00570 <span class="comment"> * ellipsoid height, using the EGM96 gravity model and the bilinear interpolation method.</span>
596
<a name="l00571"></a>00571 <span class="comment"> *</span>
597
<a name="l00572"></a>00572 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
598
<a name="l00573"></a>00573 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
599
<a name="l00574"></a>00574 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
600
<a name="l00575"></a>00575 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
601
<a name="l00576"></a>00576 <span class="comment"> *</span>
602
<a name="l00577"></a>00577 <span class="comment"> */</span>
603
<a name="l00578"></a>00578
604
<a name="l00579"></a>00579 <span class="keywordtype">double</span> delta_height;
605
<a name="l00580"></a>00580
606
<a name="l00581"></a>00581 bilinearInterpolate(
607
<a name="l00582"></a>00582 longitude, latitude,
608
<a name="l00583"></a>00583 <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>, <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>, egm96GeoidList,
609
<a name="l00584"></a>00584 &delta_height );
610
610
<a name="l00585"></a>00585
611
<a name="l00586"></a>00586
612
<a name="l00587"></a>00587 <span class="comment">/************************************************************************/</span>
613
<a name="l00588"></a>00588 <span class="comment">/* PRIVATE FUNCTIONS </span>
614
<a name="l00589"></a>00589 <span class="comment"> *</span>
615
<a name="l00590"></a>00590 <span class="comment"> */</span>
616
<a name="l00591"></a>00591
617
<a name="l00592"></a>00592 <span class="keywordtype">void</span> GeoidLibrary::initializeEGM96Geoid()
618
<a name="l00593"></a>00593 {
619
<a name="l00594"></a>00594 <span class="comment">/*</span>
620
<a name="l00595"></a>00595 <span class="comment"> * The function initializeEGM96Geoid reads geoid separation data from the egm96.grd</span>
621
<a name="l00596"></a>00596 <span class="comment"> * file in the current directory and builds a geoid separation table from it.</span>
622
<a name="l00597"></a>00597 <span class="comment"> * If the separation file can not be found or accessed, an error code of</span>
623
<a name="l00598"></a>00598 <span class="comment"> * GEOID_FILE_OPEN_ERROR is returned, If the separation file is incomplete</span>
624
<a name="l00599"></a>00599 <span class="comment"> * or improperly formatted, an error code of GEOID_INITIALIZE_ERROR is returned,</span>
625
<a name="l00600"></a>00600 <span class="comment"> * otherwise GEOID_NO_ERROR is returned.</span>
626
<a name="l00601"></a>00601 <span class="comment"> */</span>
627
<a name="l00602"></a>00602
628
<a name="l00603"></a>00603 <span class="keywordtype">int</span> items_read = 0;
629
<a name="l00604"></a>00604 <span class="keywordtype">char</span>* file_name = 0;
630
<a name="l00605"></a>00605 <span class="keywordtype">char</span>* path_name = getenv( <span class="stringliteral">"MSPCCS_DATA"</span> );
631
<a name="l00606"></a>00606 <span class="keywordtype">long</span> elevations_read = 0;
632
<a name="l00607"></a>00607 <span class="keywordtype">long</span> items_discarded = 0;
633
<a name="l00608"></a>00608 <span class="keywordtype">long</span> num = 0;
634
<a name="l00609"></a>00609 FILE* geoid_height_file;
635
<a name="l00610"></a>00610
636
<a name="l00611"></a>00611 <a class="code" href="threads_8h.html#a662f8c1664a8ad6cf2bee661d61637b1">Thread_Mutex</a> geoid_96_mutex;
637
<a name="l00612"></a>00612 <span class="keywordtype">long</span> mutex_error = <a class="code" href="threads_8cpp.html#a12b23b7654d3bfb8595926789a6d5cb4">Threads_Create_Mutex</a>( &geoid_96_mutex );
638
<a name="l00613"></a>00613 <span class="keywordflow">if</span>( !mutex_error )
639
<a name="l00614"></a>00614 mutex_error = <a class="code" href="threads_8cpp.html#a50af35318da3e29aef5ae77632972927">Threads_Lock_Mutex</a>( geoid_96_mutex );
640
<a name="l00615"></a>00615
641
<a name="l00616"></a>00616 <span class="comment">/* Check the environment for a user provided path, else current directory; */</span>
642
<a name="l00617"></a>00617 <span class="comment">/* Build a File Name, including specified or default path: */</span>
643
<a name="l00618"></a>00618
644
<a name="l00619"></a>00619 <span class="keywordflow">if</span> (path_name != NULL)
645
<a name="l00620"></a>00620 {
646
<a name="l00621"></a>00621 file_name = <span class="keyword">new</span> <span class="keywordtype">char</span>[ strlen( path_name ) + 11 ];
647
<a name="l00622"></a>00622 strcpy( file_name, path_name );
648
<a name="l00623"></a>00623 strcat( file_name, <span class="stringliteral">"/"</span> );
649
<a name="l00624"></a>00624 }
650
<a name="l00625"></a>00625 <span class="keywordflow">else</span>
651
<a name="l00626"></a>00626 {
652
<a name="l00627"></a>00627 file_name = <span class="keyword">new</span> <span class="keywordtype">char</span>[ 18 ];
653
<a name="l00628"></a>00628 strcpy( file_name, <span class="stringliteral">"../../data/"</span> );
654
<a name="l00629"></a>00629 }
655
<a name="l00630"></a>00630 strcat( file_name, <span class="stringliteral">"egm96.grd"</span> );
656
<a name="l00631"></a>00631
657
<a name="l00632"></a>00632 <span class="comment">/* Open the File READONLY, or Return Error Condition: */</span>
658
<a name="l00633"></a>00633
659
<a name="l00634"></a>00634 <span class="keywordflow">if</span> ( ( geoid_height_file = fopen( file_name, <span class="stringliteral">"rb"</span> ) ) == NULL )
660
<a name="l00635"></a>00635 {
661
<a name="l00636"></a>00636 <span class="keyword">delete</span> [] file_name;
662
<a name="l00637"></a>00637 file_name = 0;
663
<a name="l00638"></a>00638
664
<a name="l00639"></a>00639 <span class="keywordflow">if</span>( !mutex_error )
665
<a name="l00640"></a>00640 mutex_error = <a class="code" href="threads_8cpp.html#af5eb9be225bb5c55f4592d1fd5fc527c">Threads_Unlock_Mutex</a>( geoid_96_mutex );
666
<a name="l00641"></a>00641 <span class="keywordflow">if</span>( !mutex_error )
667
<a name="l00642"></a>00642 <a class="code" href="threads_8cpp.html#a5fa63bacc9297e526ec2ee3e9388eac7">Threads_Destroy_Mutex</a>( geoid_96_mutex );
668
<a name="l00643"></a>00643
669
<a name="l00644"></a>00644 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
670
<a name="l00645"></a>00645 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a40a2ec1d77ae260338773a023518d449">ErrorMessages::geoidFileOpenError</a> );
671
<a name="l00646"></a>00646 strcat( message, <span class="stringliteral">": egm96.grd\n"</span> );
672
<a name="l00647"></a>00647 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
673
<a name="l00648"></a>00648 }
674
<a name="l00649"></a>00649
675
<a name="l00650"></a>00650 <span class="comment">/* Skip the Header Line: */</span>
611
<a name="l00586"></a>00586 *ellipsoidHeight = geoidHeight + delta_height;
612
<a name="l00587"></a>00587 }
613
<a name="l00588"></a>00588
614
<a name="l00589"></a>00589
615
<a name="l00590"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#aa3f8de2392ade95991fd1519a7f869fa">00590</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#aa3f8de2392ade95991fd1519a7f869fa">GeoidLibrary::convertEGM96VariableNaturalSplineToEllipsoidHeight</a>(
616
<a name="l00591"></a>00591 <span class="keywordtype">double</span> longitude,
617
<a name="l00592"></a>00592 <span class="keywordtype">double</span> latitude,
618
<a name="l00593"></a>00593 <span class="keywordtype">double</span> geoidHeight,
619
<a name="l00594"></a>00594 <span class="keywordtype">double</span> *ellipsoidHeight )
620
<a name="l00595"></a>00595 {
621
<a name="l00596"></a>00596 <span class="comment">/*</span>
622
<a name="l00597"></a>00597 <span class="comment"> * The function convertEGM96VariableNaturalSplineToEllipsoidHeight converts the specified WGS84</span>
623
<a name="l00598"></a>00598 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
624
<a name="l00599"></a>00599 <span class="comment"> * ellipsoid height, using the EGM96 gravity model and the natural spline interpolation method.</span>
625
<a name="l00600"></a>00600 <span class="comment"> *</span>
626
<a name="l00601"></a>00601 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
627
<a name="l00602"></a>00602 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
628
<a name="l00603"></a>00603 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
629
<a name="l00604"></a>00604 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
630
<a name="l00605"></a>00605 <span class="comment"> *</span>
631
<a name="l00606"></a>00606 <span class="comment"> */</span>
632
<a name="l00607"></a>00607
633
<a name="l00608"></a>00608 <span class="keywordtype">int</span> i = 0;
634
<a name="l00609"></a>00609 <span class="keywordtype">int</span> num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
635
<a name="l00610"></a>00610 <span class="keywordtype">int</span> num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
636
<a name="l00611"></a>00611 <span class="keywordtype">double</span> latitude_degrees = latitude * <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a>;
637
<a name="l00612"></a>00612 <span class="keywordtype">double</span> longitude_degrees = longitude * _180_OVER_PI;
638
<a name="l00613"></a>00613 <span class="keywordtype">double</span> scale_factor = <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a>;
639
<a name="l00614"></a>00614 <span class="keywordtype">double</span> delta_height;
640
<a name="l00615"></a>00615 <span class="keywordtype">long</span> found = 0;
641
<a name="l00616"></a>00616
642
<a name="l00617"></a>00617 <span class="keywordflow">if</span>( longitude_degrees < 0.0 )
643
<a name="l00618"></a>00618 longitude_degrees += 360.0;
644
<a name="l00619"></a>00619
645
<a name="l00620"></a>00620 <span class="keywordflow">while</span>( !found && i < <a class="code" href="_geoid_library_8cpp.html#adc700a7cefac4f0369fb4db08cdf20a4">EGM96_INSET_AREAS</a> )
646
<a name="l00621"></a>00621 {
647
<a name="l00622"></a>00622 <span class="keywordflow">if</span>(( latitude_degrees >= EGM96_Variable_Grid_Table[i].min_lat ) &&
648
<a name="l00623"></a>00623 ( longitude_degrees >= EGM96_Variable_Grid_Table[i].min_lon ) &&
649
<a name="l00624"></a>00624 ( latitude_degrees < EGM96_Variable_Grid_Table[i].max_lat ) &&
650
<a name="l00625"></a>00625 ( longitude_degrees < EGM96_Variable_Grid_Table[i].max_lon ) )
651
<a name="l00626"></a>00626 {
652
<a name="l00627"></a>00627 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a>; <span class="comment">// use 30 minute by 30 minute grid</span>
653
<a name="l00628"></a>00628 num_cols = 721;
654
<a name="l00629"></a>00629 num_rows = 361;
655
<a name="l00630"></a>00630 found = 1;
656
<a name="l00631"></a>00631 }
657
<a name="l00632"></a>00632
658
<a name="l00633"></a>00633 i++;
659
<a name="l00634"></a>00634 }
660
<a name="l00635"></a>00635
661
<a name="l00636"></a>00636 <span class="keywordflow">if</span>( !found )
662
<a name="l00637"></a>00637 {
663
<a name="l00638"></a>00638 <span class="keywordflow">if</span>( latitude_degrees >= -60.0 && latitude_degrees < 60.0 )
664
<a name="l00639"></a>00639 {
665
<a name="l00640"></a>00640 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a>; <span class="comment">// use 1 degree by 1 degree grid</span>
666
<a name="l00641"></a>00641 num_cols = 361;
667
<a name="l00642"></a>00642 num_rows = 181;
668
<a name="l00643"></a>00643 }
669
<a name="l00644"></a>00644 <span class="keywordflow">else</span>
670
<a name="l00645"></a>00645 {
671
<a name="l00646"></a>00646 scale_factor = <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a>; <span class="comment">// use 2 degree by 2 degree grid</span>
672
<a name="l00647"></a>00647 num_cols = 181;
673
<a name="l00648"></a>00648 num_rows = 91;
674
<a name="l00649"></a>00649 }
675
<a name="l00650"></a>00650 }
676
676
<a name="l00651"></a>00651
677
<a name="l00652"></a>00652 <span class="keywordflow">while</span> ( num < <a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">EGM96_HEADER_ITEMS</a> )
678
<a name="l00653"></a>00653 {
679
<a name="l00654"></a>00654 <span class="keywordflow">if</span> ( feof( geoid_height_file ) ) <span class="keywordflow">break</span>;
680
<a name="l00655"></a>00655 <span class="keywordflow">if</span> ( ferror( geoid_height_file ) ) <span class="keywordflow">break</span>;
681
<a name="l00656"></a>00656 egm96GeoidList.push_back( <a class="code" href="_geoid_library_8cpp.html#ac0f3a476600c87d7296394be3c5b047d">readFloat</a>( &items_read, geoid_height_file ) );
682
<a name="l00657"></a>00657 items_discarded += items_read;
683
<a name="l00658"></a>00658 num++;
684
<a name="l00659"></a>00659 }
677
<a name="l00652"></a>00652 naturalSplineInterpolate(
678
<a name="l00653"></a>00653 longitude, latitude,
679
<a name="l00654"></a>00654 scale_factor, num_cols, num_rows, <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>-1, egm96GeoidList,
680
<a name="l00655"></a>00655 &delta_height );
681
<a name="l00656"></a>00656
682
<a name="l00657"></a>00657 *ellipsoidHeight = geoidHeight + delta_height;
683
<a name="l00658"></a>00658 }
684
<a name="l00659"></a>00659
685
685
<a name="l00660"></a>00660
686
<a name="l00661"></a>00661 <span class="comment">/* Determine if header read properly, or NOT: */</span>
687
<a name="l00662"></a>00662
688
<a name="l00663"></a>00663 <span class="keywordflow">if</span>( egm96GeoidList[0] != -90.0 ||
689
<a name="l00664"></a>00664 egm96GeoidList[1] != 90.0 ||
690
<a name="l00665"></a>00665 egm96GeoidList[2] != 0.0 ||
691
<a name="l00666"></a>00666 egm96GeoidList[3] != 360.0 ||
692
<a name="l00667"></a>00667 egm96GeoidList[4] != <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a> ||
693
<a name="l00668"></a>00668 egm96GeoidList[5] != <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a> ||
694
<a name="l00669"></a>00669 items_discarded != <a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">EGM96_HEADER_ITEMS</a> )
695
<a name="l00670"></a>00670 {
696
<a name="l00671"></a>00671 fclose( geoid_height_file );
697
<a name="l00672"></a>00672 <span class="keyword">delete</span> [] file_name;
698
<a name="l00673"></a>00673 file_name = 0;
686
<a name="l00661"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a68378c47cf83b97e59db68c7f287c4bc">00661</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a68378c47cf83b97e59db68c7f287c4bc">GeoidLibrary::convertEGM84TenDegBilinearToEllipsoidHeight</a>( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> geoidHeight, <span class="keywordtype">double</span> *ellipsoidHeight )
687
<a name="l00662"></a>00662 {
688
<a name="l00663"></a>00663 <span class="comment">/*</span>
689
<a name="l00664"></a>00664 <span class="comment"> * The function convertEGM84TenDegBilinearToEllipsoidHeight converts the specified WGS84</span>
690
<a name="l00665"></a>00665 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
691
<a name="l00666"></a>00666 <span class="comment"> * ellipsoid height, using the EGM84 gravity model and the bilinear interpolation method.</span>
692
<a name="l00667"></a>00667 <span class="comment"> *</span>
693
<a name="l00668"></a>00668 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
694
<a name="l00669"></a>00669 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
695
<a name="l00670"></a>00670 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
696
<a name="l00671"></a>00671 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
697
<a name="l00672"></a>00672 <span class="comment"> *</span>
698
<a name="l00673"></a>00673 <span class="comment"> */</span>
699
699
<a name="l00674"></a>00674
700
<a name="l00675"></a>00675 <span class="keywordflow">if</span>( !mutex_error )
701
<a name="l00676"></a>00676 mutex_error = <a class="code" href="threads_8cpp.html#af5eb9be225bb5c55f4592d1fd5fc527c">Threads_Unlock_Mutex</a>( geoid_96_mutex );
702
<a name="l00677"></a>00677 <span class="keywordflow">if</span>( !mutex_error )
703
<a name="l00678"></a>00678 <a class="code" href="threads_8cpp.html#a5fa63bacc9297e526ec2ee3e9388eac7">Threads_Destroy_Mutex</a>( geoid_96_mutex );
704
<a name="l00679"></a>00679
705
<a name="l00680"></a>00680 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
706
<a name="l00681"></a>00681 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#af79e7384c2f3f25e9a45c4a63c14bfcb">ErrorMessages::geoidFileParseError</a> );
707
<a name="l00682"></a>00682 strcat( message, <span class="stringliteral">": egm96.grd\n"</span> );
708
<a name="l00683"></a>00683 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
709
<a name="l00684"></a>00684 }
700
<a name="l00675"></a>00675 <span class="keywordtype">double</span> delta_height;
701
<a name="l00676"></a>00676
702
<a name="l00677"></a>00677 bilinearInterpolate(
703
<a name="l00678"></a>00678 longitude, latitude,
704
<a name="l00679"></a>00679 <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>, <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, egm84GeoidList,
705
<a name="l00680"></a>00680 &delta_height );
706
<a name="l00681"></a>00681
707
<a name="l00682"></a>00682 *ellipsoidHeight = geoidHeight + delta_height;
708
<a name="l00683"></a>00683 }
709
<a name="l00684"></a>00684
710
710
<a name="l00685"></a>00685
711
<a name="l00686"></a>00686 <span class="comment">/* Extract elements from the file: */</span>
712
<a name="l00687"></a>00687 egm96GeoidList.clear();
713
<a name="l00688"></a>00688
714
<a name="l00689"></a>00689 num = 0;
715
<a name="l00690"></a>00690 <span class="keywordflow">while</span> ( num < <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a> )
716
<a name="l00691"></a>00691 {
717
<a name="l00692"></a>00692 <span class="keywordflow">if</span> ( feof( geoid_height_file ) ) <span class="keywordflow">break</span>;
718
<a name="l00693"></a>00693 <span class="keywordflow">if</span> ( ferror( geoid_height_file ) ) <span class="keywordflow">break</span>;
719
<a name="l00694"></a>00694 egm96GeoidList.push_back( <a class="code" href="_geoid_library_8cpp.html#ac0f3a476600c87d7296394be3c5b047d">readFloat</a> ( &items_read, geoid_height_file ) );
720
<a name="l00695"></a>00695 elevations_read += items_read;
721
<a name="l00696"></a>00696 num++;
722
<a name="l00697"></a>00697 }
723
<a name="l00698"></a>00698
724
<a name="l00699"></a>00699 <span class="comment">/* Determine if all elevations of file read properly, or NOT: */</span>
725
<a name="l00700"></a>00700
726
<a name="l00701"></a>00701 <span class="keywordflow">if</span> ( elevations_read != <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a> )
727
<a name="l00702"></a>00702 {
728
<a name="l00703"></a>00703 fclose(geoid_height_file);
729
<a name="l00704"></a>00704 <span class="keyword">delete</span> [] file_name;
730
<a name="l00705"></a>00705 file_name = 0;
731
<a name="l00706"></a>00706
732
<a name="l00707"></a>00707 <span class="keywordflow">if</span>( !mutex_error )
733
<a name="l00708"></a>00708 mutex_error = <a class="code" href="threads_8cpp.html#af5eb9be225bb5c55f4592d1fd5fc527c">Threads_Unlock_Mutex</a>( geoid_96_mutex );
734
<a name="l00709"></a>00709 <span class="keywordflow">if</span>( !mutex_error )
735
<a name="l00710"></a>00710 <a class="code" href="threads_8cpp.html#a5fa63bacc9297e526ec2ee3e9388eac7">Threads_Destroy_Mutex</a>( geoid_96_mutex );
736
<a name="l00711"></a>00711
737
<a name="l00712"></a>00712 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
738
<a name="l00713"></a>00713 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#af79e7384c2f3f25e9a45c4a63c14bfcb">ErrorMessages::geoidFileParseError</a> );
739
<a name="l00714"></a>00714 strcat( message, <span class="stringliteral">": egm96.grd\n"</span> );
740
<a name="l00715"></a>00715 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
741
<a name="l00716"></a>00716 }
742
<a name="l00717"></a>00717
743
<a name="l00718"></a>00718 fclose( geoid_height_file );
744
<a name="l00719"></a>00719 geoid_height_file = 0;
745
<a name="l00720"></a>00720
746
<a name="l00721"></a>00721 <span class="keywordflow">if</span>( !mutex_error )
747
<a name="l00722"></a>00722 mutex_error = <a class="code" href="threads_8cpp.html#af5eb9be225bb5c55f4592d1fd5fc527c">Threads_Unlock_Mutex</a>( geoid_96_mutex );
748
<a name="l00723"></a>00723 <span class="keywordflow">if</span>( !mutex_error )
749
<a name="l00724"></a>00724 <a class="code" href="threads_8cpp.html#a5fa63bacc9297e526ec2ee3e9388eac7">Threads_Destroy_Mutex</a>( geoid_96_mutex );
750
<a name="l00725"></a>00725
751
<a name="l00726"></a>00726 <span class="keyword">delete</span> [] file_name;
752
<a name="l00727"></a>00727 file_name = 0;
753
<a name="l00728"></a>00728 }
754
<a name="l00729"></a>00729
755
<a name="l00730"></a>00730
756
<a name="l00731"></a>00731 <span class="keywordtype">void</span> GeoidLibrary::initializeEGM84Geoid()
757
<a name="l00732"></a>00732 {
758
<a name="l00733"></a>00733 <span class="comment">/*</span>
759
<a name="l00734"></a>00734 <span class="comment"> * The function initializeEGM84Geoid reads geoid separation data from a file in</span>
760
<a name="l00735"></a>00735 <span class="comment"> * the current directory and builds the geoid separation table from it.</span>
761
<a name="l00736"></a>00736 <span class="comment"> * If the separation file can not be found or accessed, an error code of</span>
762
<a name="l00737"></a>00737 <span class="comment"> * GEOID_FILE_OPEN_ERROR is returned, If the separation file is incomplete</span>
763
<a name="l00738"></a>00738 <span class="comment"> * or improperly formatted, an error code of GEOID_INITIALIZE_ERROR is returned,</span>
764
<a name="l00739"></a>00739 <span class="comment"> * otherwise GEOID_NO_ERROR is returned.</span>
765
<a name="l00740"></a>00740 <span class="comment"> */</span>
766
<a name="l00741"></a>00741
767
<a name="l00742"></a>00742 <span class="keywordtype">int</span> items_read = 0;
768
<a name="l00743"></a>00743 <span class="keywordtype">char</span>* file_name = 0;
769
<a name="l00744"></a>00744 <span class="keywordtype">char</span>* path_name = getenv( <span class="stringliteral">"MSPCCS_DATA"</span> );
770
<a name="l00745"></a>00745 <span class="keywordtype">long</span> elevations_read = 0;
771
<a name="l00746"></a>00746 <span class="keywordtype">long</span> num = 0;
772
<a name="l00747"></a>00747 FILE* geoid_height_file;
773
<a name="l00748"></a>00748
774
<a name="l00749"></a>00749 <a class="code" href="threads_8h.html#a662f8c1664a8ad6cf2bee661d61637b1">Thread_Mutex</a> geoid_84_mutex;
775
<a name="l00750"></a>00750 <span class="keywordtype">long</span> mutex_error = <a class="code" href="threads_8cpp.html#a12b23b7654d3bfb8595926789a6d5cb4">Threads_Create_Mutex</a>( &geoid_84_mutex );
776
<a name="l00751"></a>00751 <span class="keywordflow">if</span>( !mutex_error )
777
<a name="l00752"></a>00752 mutex_error = <a class="code" href="threads_8cpp.html#a50af35318da3e29aef5ae77632972927">Threads_Lock_Mutex</a>( geoid_84_mutex );
778
<a name="l00753"></a>00753
779
<a name="l00754"></a>00754 <span class="comment">/* Check the environment for a user provided path, else current directory; */</span>
780
<a name="l00755"></a>00755 <span class="comment">/* Build a File Name, including specified or default path: */</span>
711
<a name="l00686"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a8f085fab24f9cd29bc1ad8d75c2ec798">00686</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a8f085fab24f9cd29bc1ad8d75c2ec798">GeoidLibrary::convertEGM84TenDegNaturalSplineToEllipsoidHeight</a>(
712
<a name="l00687"></a>00687 <span class="keywordtype">double</span> longitude,
713
<a name="l00688"></a>00688 <span class="keywordtype">double</span> latitude,
714
<a name="l00689"></a>00689 <span class="keywordtype">double</span> geoidHeight,
715
<a name="l00690"></a>00690 <span class="keywordtype">double</span> *ellipsoidHeight )
716
<a name="l00691"></a>00691 {
717
<a name="l00692"></a>00692 <span class="comment">/*</span>
718
<a name="l00693"></a>00693 <span class="comment"> * The function convertEGM84TenDegNaturalSplineToEllipsoidHeight converts the specified WGS84</span>
719
<a name="l00694"></a>00694 <span class="comment"> * geoid height at the specified geodetic coordinates to the equivalent</span>
720
<a name="l00695"></a>00695 <span class="comment"> * ellipsoid height, using the EGM84 gravity model and the natural spline interpolation method.</span>
721
<a name="l00696"></a>00696 <span class="comment"> *</span>
722
<a name="l00697"></a>00697 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
723
<a name="l00698"></a>00698 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
724
<a name="l00699"></a>00699 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
725
<a name="l00700"></a>00700 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
726
<a name="l00701"></a>00701 <span class="comment"> *</span>
727
<a name="l00702"></a>00702 <span class="comment"> */</span>
728
<a name="l00703"></a>00703
729
<a name="l00704"></a>00704 <span class="keywordtype">double</span> delta_height;
730
<a name="l00705"></a>00705
731
<a name="l00706"></a>00706 naturalSplineInterpolate(
732
<a name="l00707"></a>00707 longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#a68b1126f892291cb2f894021ba1068c4">SCALE_FACTOR_10_DEGREES</a>,
733
<a name="l00708"></a>00708 <a class="code" href="_geoid_library_8cpp.html#a99ccb09f1fef806451f84d9710ba26ca">EGM84_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#aaeb71f36b136d1e1b54901920f01f8a0">EGM84_ROWS</a>, <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>-1,
734
<a name="l00709"></a>00709 egm84GeoidList, &delta_height );
735
<a name="l00710"></a>00710
736
<a name="l00711"></a>00711 *ellipsoidHeight = geoidHeight + delta_height;
737
<a name="l00712"></a>00712 }
738
<a name="l00713"></a>00713
739
<a name="l00714"></a>00714
740
<a name="l00715"></a><a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a5c27cd4bc3de0f2b5525a44dba466f86">00715</a> <span class="keywordtype">void</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_geoid_library.html#a5c27cd4bc3de0f2b5525a44dba466f86">GeoidLibrary::convertEGM84ThirtyMinBiLinearToEllipsoidHeight</a>(
741
<a name="l00716"></a>00716 <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> geoidHeight,
742
<a name="l00717"></a>00717 <span class="keywordtype">double</span> *ellipsoidHeight )
743
<a name="l00718"></a>00718 {
744
<a name="l00719"></a>00719 <span class="comment">/*</span>
745
<a name="l00720"></a>00720 <span class="comment"> * The function convertEGM84ThirtyMinBiLinearToEllipsoidHeight converts the </span>
746
<a name="l00721"></a>00721 <span class="comment"> * specified WGS84 geoid height at the specified geodetic coordinates to the </span>
747
<a name="l00722"></a>00722 <span class="comment"> * equivalent ellipsoid height, using the EGM84 gravity model and the bilinear </span>
748
<a name="l00723"></a>00723 <span class="comment"> * interpolation method.</span>
749
<a name="l00724"></a>00724 <span class="comment"> *</span>
750
<a name="l00725"></a>00725 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
751
<a name="l00726"></a>00726 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
752
<a name="l00727"></a>00727 <span class="comment"> * geoidHeight : Geoid height, in meters. (input)</span>
753
<a name="l00728"></a>00728 <span class="comment"> * ellipsoidHeight : Ellipsoid height, in meters (output)</span>
754
<a name="l00729"></a>00729 <span class="comment"> *</span>
755
<a name="l00730"></a>00730 <span class="comment"> */</span>
756
<a name="l00731"></a>00731
757
<a name="l00732"></a>00732 <span class="keywordtype">double</span> delta_height;
758
<a name="l00733"></a>00733
759
<a name="l00734"></a>00734 bilinearInterpolateDoubleHeights( longitude, latitude, <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a>,
760
<a name="l00735"></a>00735 <a class="code" href="_geoid_library_8cpp.html#a580782fa058e200088161c1680921220">EGM84_30_MIN_COLS</a>, <a class="code" href="_geoid_library_8cpp.html#a1ea8fba0f5bc5035bc21da731b6d8799">EGM84_30_MIN_ROWS</a>, egm84ThirtyMinGeoidList,
761
<a name="l00736"></a>00736 &delta_height );
762
<a name="l00737"></a>00737 *ellipsoidHeight = geoidHeight + delta_height;
763
<a name="l00738"></a>00738 }
764
<a name="l00739"></a>00739
765
<a name="l00740"></a>00740
766
<a name="l00741"></a>00741 <span class="comment">/************************************************************************/</span>
767
<a name="l00742"></a>00742 <span class="comment">/* PRIVATE FUNCTIONS </span>
768
<a name="l00743"></a>00743 <span class="comment"> *</span>
769
<a name="l00744"></a>00744 <span class="comment"> */</span>
770
<a name="l00745"></a>00745
771
<a name="l00746"></a>00746 <span class="keywordtype">void</span> GeoidLibrary::initializeEGM96Geoid()
772
<a name="l00747"></a>00747 {
773
<a name="l00748"></a>00748 <span class="comment">/*</span>
774
<a name="l00749"></a>00749 <span class="comment"> * The function initializeEGM96Geoid reads geoid separation data from the egm96.grd</span>
775
<a name="l00750"></a>00750 <span class="comment"> * file in the current directory and builds a geoid separation table from it.</span>
776
<a name="l00751"></a>00751 <span class="comment"> * If the separation file can not be found or accessed, an error code of</span>
777
<a name="l00752"></a>00752 <span class="comment"> * GEOID_FILE_OPEN_ERROR is returned, If the separation file is incomplete</span>
778
<a name="l00753"></a>00753 <span class="comment"> * or improperly formatted, an error code of GEOID_INITIALIZE_ERROR is returned,</span>
779
<a name="l00754"></a>00754 <span class="comment"> * otherwise GEOID_NO_ERROR is returned.</span>
780
<a name="l00755"></a>00755 <span class="comment"> */</span>
781
781
<a name="l00756"></a>00756
782
<a name="l00757"></a>00757 <span class="keywordflow">if</span> (path_name != NULL)
783
<a name="l00758"></a>00758 {
784
<a name="l00759"></a>00759 file_name = <span class="keyword">new</span> <span class="keywordtype">char</span>[ strlen( path_name ) + 11 ];
785
<a name="l00760"></a>00760 strcpy( file_name, path_name );
786
<a name="l00761"></a>00761 strcat( file_name, <span class="stringliteral">"/"</span> );
787
<a name="l00762"></a>00762 }
788
<a name="l00763"></a>00763 <span class="keywordflow">else</span>
789
<a name="l00764"></a>00764 {
790
<a name="l00765"></a>00765 file_name =<span class="keyword">new</span> <span class="keywordtype">char</span>[ 18 ];
791
<a name="l00766"></a>00766 strcpy( file_name, <span class="stringliteral">"../../data/"</span> );
792
<a name="l00767"></a>00767 }
793
<a name="l00768"></a>00768 strcat( file_name, <span class="stringliteral">"egm84.grd"</span> );
782
<a name="l00757"></a>00757 <span class="keywordtype">int</span> items_read = 0;
783
<a name="l00758"></a>00758 <span class="keywordtype">char</span>* file_name = 0;
784
<a name="l00759"></a>00759 <span class="keywordtype">char</span>* path_name = getenv( <span class="stringliteral">"MSPCCS_DATA"</span> );
785
<a name="l00760"></a>00760 <span class="keywordtype">long</span> elevations_read = 0;
786
<a name="l00761"></a>00761 <span class="keywordtype">long</span> items_discarded = 0;
787
<a name="l00762"></a>00762 <span class="keywordtype">long</span> num = 0;
788
<a name="l00763"></a>00763 FILE* geoid_height_file;
789
<a name="l00764"></a>00764
790
<a name="l00765"></a>00765 <a class="code" href="class_m_s_p_1_1_c_c_s_thread_lock.html">CCSThreadLock</a> lock(&mutex);
791
<a name="l00766"></a>00766
792
<a name="l00767"></a>00767 <span class="comment">/* Check the environment for a user provided path, else current directory; */</span>
793
<a name="l00768"></a>00768 <span class="comment">/* Build a File Name, including specified or default path: */</span>
794
794
<a name="l00769"></a>00769
795
<a name="l00770"></a>00770 <span class="comment">/* Open the File READONLY, or Return Error Condition: */</span>
796
<a name="l00771"></a>00771
797
<a name="l00772"></a>00772 <span class="keywordflow">if</span>( ( geoid_height_file = fopen( file_name, <span class="stringliteral">"rb"</span> ) ) == NULL )
798
<a name="l00773"></a>00773 {
799
<a name="l00774"></a>00774 <span class="keyword">delete</span> [] file_name;
800
<a name="l00775"></a>00775 file_name = 0;
801
<a name="l00776"></a>00776
802
<a name="l00777"></a>00777 <span class="keywordflow">if</span>( !mutex_error )
803
<a name="l00778"></a>00778 mutex_error = <a class="code" href="threads_8cpp.html#af5eb9be225bb5c55f4592d1fd5fc527c">Threads_Unlock_Mutex</a>( geoid_84_mutex );
804
<a name="l00779"></a>00779 <span class="keywordflow">if</span>( !mutex_error )
805
<a name="l00780"></a>00780 <a class="code" href="threads_8cpp.html#a5fa63bacc9297e526ec2ee3e9388eac7">Threads_Destroy_Mutex</a>( geoid_84_mutex );
806
<a name="l00781"></a>00781
807
<a name="l00782"></a>00782 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
808
<a name="l00783"></a>00783 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a40a2ec1d77ae260338773a023518d449">ErrorMessages::geoidFileOpenError</a> );
809
<a name="l00784"></a>00784 strcat( message, <span class="stringliteral">": egm84.grd\n"</span> );
810
<a name="l00785"></a>00785 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
811
<a name="l00786"></a>00786 }
812
<a name="l00787"></a>00787
813
<a name="l00788"></a>00788
814
<a name="l00789"></a>00789 <span class="comment">/* Extract elements from the file: */</span>
815
<a name="l00790"></a>00790
816
<a name="l00791"></a>00791 num = 0;
817
<a name="l00792"></a>00792 <span class="keywordflow">while</span>( num < <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a> )
818
<a name="l00793"></a>00793 {
819
<a name="l00794"></a>00794 <span class="keywordflow">if</span> ( feof( geoid_height_file ) ) <span class="keywordflow">break</span>;
820
<a name="l00795"></a>00795 <span class="keywordflow">if</span> ( ferror( geoid_height_file ) ) <span class="keywordflow">break</span>;
821
<a name="l00796"></a>00796 egm84GeoidList.push_back( <a class="code" href="_geoid_library_8cpp.html#ac0f3a476600c87d7296394be3c5b047d">readFloat</a>( &items_read, geoid_height_file ) );
795
<a name="l00770"></a>00770 <span class="keywordflow">if</span> (path_name != NULL)
796
<a name="l00771"></a>00771 {
797
<a name="l00772"></a>00772 file_name = <span class="keyword">new</span> <span class="keywordtype">char</span>[ strlen( path_name ) + 11 ];
798
<a name="l00773"></a>00773 strcpy( file_name, path_name );
799
<a name="l00774"></a>00774 strcat( file_name, <span class="stringliteral">"/"</span> );
800
<a name="l00775"></a>00775 }
801
<a name="l00776"></a>00776 <span class="keywordflow">else</span>
802
<a name="l00777"></a>00777 {
803
<a name="l00778"></a>00778 file_name = <span class="keyword">new</span> <span class="keywordtype">char</span>[ 21 ];
804
<a name="l00779"></a>00779 strcpy( file_name, <span class="stringliteral">"../../data/"</span> );
805
<a name="l00780"></a>00780 }
806
<a name="l00781"></a>00781 strcat( file_name, <span class="stringliteral">"egm96.grd"</span> );
807
<a name="l00782"></a>00782
808
<a name="l00783"></a>00783 <span class="comment">/* Open the File READONLY, or Return Error Condition: */</span>
809
<a name="l00784"></a>00784
810
<a name="l00785"></a>00785 <span class="keywordflow">if</span> ( ( geoid_height_file = fopen( file_name, <span class="stringliteral">"rb"</span> ) ) == NULL )
811
<a name="l00786"></a>00786 {
812
<a name="l00787"></a>00787 <span class="keyword">delete</span> [] file_name;
813
<a name="l00788"></a>00788 file_name = 0;
814
<a name="l00789"></a>00789
815
<a name="l00790"></a>00790 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
816
<a name="l00791"></a>00791 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a40a2ec1d77ae260338773a023518d449">ErrorMessages::geoidFileOpenError</a> );
817
<a name="l00792"></a>00792 strcat( message, <span class="stringliteral">": egm96.grd\n"</span> );
818
<a name="l00793"></a>00793 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
819
<a name="l00794"></a>00794 }
820
<a name="l00795"></a>00795
821
<a name="l00796"></a>00796 <span class="comment">/* Skip the Header Line: */</span>
822
822
<a name="l00797"></a>00797
823
<a name="l00798"></a>00798 elevations_read += items_read;
824
<a name="l00799"></a>00799 num++;
825
<a name="l00800"></a>00800 }
826
<a name="l00801"></a>00801
827
<a name="l00802"></a>00802
828
<a name="l00803"></a>00803 <span class="comment">/* Determine if all elevations of file read properly, or NOT: */</span>
829
<a name="l00804"></a>00804
830
<a name="l00805"></a>00805 <span class="keywordflow">if</span>( ( elevations_read ) != <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a> )
831
<a name="l00806"></a>00806 {
832
<a name="l00807"></a>00807 fclose( geoid_height_file );
833
<a name="l00808"></a>00808 <span class="keyword">delete</span> [] file_name;
834
<a name="l00809"></a>00809 file_name = 0;
835
<a name="l00810"></a>00810
836
<a name="l00811"></a>00811 <span class="keywordflow">if</span>( !mutex_error )
837
<a name="l00812"></a>00812 mutex_error = <a class="code" href="threads_8cpp.html#af5eb9be225bb5c55f4592d1fd5fc527c">Threads_Unlock_Mutex</a>( geoid_84_mutex );
838
<a name="l00813"></a>00813 <span class="keywordflow">if</span>( !mutex_error )
839
<a name="l00814"></a>00814 <a class="code" href="threads_8cpp.html#a5fa63bacc9297e526ec2ee3e9388eac7">Threads_Destroy_Mutex</a>( geoid_84_mutex );
823
<a name="l00798"></a>00798 <span class="keywordtype">float</span> egm96GeoidHeaderList[<a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">EGM96_HEADER_ITEMS</a>];
824
<a name="l00799"></a>00799 items_discarded = <a class="code" href="_geoid_library_8cpp.html#a86b4187eeab007da13367ff3ac238250">readBinary</a>(
825
<a name="l00800"></a>00800 egm96GeoidHeaderList, 4, <a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">EGM96_HEADER_ITEMS</a>, geoid_height_file );
826
<a name="l00801"></a>00801
827
<a name="l00802"></a>00802 <span class="comment">/* Determine if header read properly, or NOT: */</span>
828
<a name="l00803"></a>00803
829
<a name="l00804"></a>00804 <span class="keywordflow">if</span>( egm96GeoidHeaderList[0] != -90.0 ||
830
<a name="l00805"></a>00805 egm96GeoidHeaderList[1] != 90.0 ||
831
<a name="l00806"></a>00806 egm96GeoidHeaderList[2] != 0.0 ||
832
<a name="l00807"></a>00807 egm96GeoidHeaderList[3] != 360.0 ||
833
<a name="l00808"></a>00808 egm96GeoidHeaderList[4] != <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a> ||
834
<a name="l00809"></a>00809 egm96GeoidHeaderList[5] != <a class="code" href="_geoid_library_8cpp.html#addc030e9b571c0bb88aeb0b4434b4cdb">SCALE_FACTOR_15_MINUTES</a> ||
835
<a name="l00810"></a>00810 items_discarded != <a class="code" href="_geoid_library_8cpp.html#a55024ccb2a1dc0f138cd29b9cb243d79">EGM96_HEADER_ITEMS</a> )
836
<a name="l00811"></a>00811 {
837
<a name="l00812"></a>00812 fclose( geoid_height_file );
838
<a name="l00813"></a>00813 <span class="keyword">delete</span> [] file_name;
839
<a name="l00814"></a>00814 file_name = 0;
840
840
<a name="l00815"></a>00815
841
841
<a name="l00816"></a>00816 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
842
842
<a name="l00817"></a>00817 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#af79e7384c2f3f25e9a45c4a63c14bfcb">ErrorMessages::geoidFileParseError</a> );
843
<a name="l00818"></a>00818 strcat( message, <span class="stringliteral">": egm84.grd\n"</span> );
843
<a name="l00818"></a>00818 strcat( message, <span class="stringliteral">": egm96.grd\n"</span> );
844
844
<a name="l00819"></a>00819 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
845
845
<a name="l00820"></a>00820 }
846
846
<a name="l00821"></a>00821
847
<a name="l00822"></a>00822 fclose( geoid_height_file );
848
<a name="l00823"></a>00823
849
<a name="l00824"></a>00824 <span class="keywordflow">if</span>( !mutex_error )
850
<a name="l00825"></a>00825 mutex_error = <a class="code" href="threads_8cpp.html#af5eb9be225bb5c55f4592d1fd5fc527c">Threads_Unlock_Mutex</a>( geoid_84_mutex );
851
<a name="l00826"></a>00826 <span class="keywordflow">if</span>( !mutex_error )
852
<a name="l00827"></a>00827 <a class="code" href="threads_8cpp.html#a5fa63bacc9297e526ec2ee3e9388eac7">Threads_Destroy_Mutex</a>( geoid_84_mutex );
853
<a name="l00828"></a>00828
854
<a name="l00829"></a>00829 <span class="keyword">delete</span> [] file_name;
855
<a name="l00830"></a>00830 file_name = 0;
856
<a name="l00831"></a>00831 }
857
<a name="l00832"></a>00832
847
<a name="l00822"></a>00822 <span class="comment">/* Extract elements from the file: */</span>
848
<a name="l00823"></a>00823 egm96GeoidList = <span class="keyword">new</span> <span class="keywordtype">float</span>[<a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>];
849
<a name="l00824"></a>00824 elevations_read = <a class="code" href="_geoid_library_8cpp.html#a86b4187eeab007da13367ff3ac238250">readBinary</a>(
850
<a name="l00825"></a>00825 egm96GeoidList, 4, <a class="code" href="_geoid_library_8cpp.html#a0be21fc9f3990edce742deda448ad66f">EGM96_ELEVATIONS</a>, geoid_height_file );
851
<a name="l00826"></a>00826
852
<a name="l00827"></a>00827 fclose( geoid_height_file );
853
<a name="l00828"></a>00828 geoid_height_file = 0;
854
<a name="l00829"></a>00829
855
<a name="l00830"></a>00830 <span class="keyword">delete</span> [] file_name;
856
<a name="l00831"></a>00831 file_name = 0;
857
<a name="l00832"></a>00832 }
858
858
<a name="l00833"></a>00833
859
<a name="l00834"></a>00834 <span class="keywordtype">void</span> GeoidLibrary::bilinearInterpolate( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> scale_factor, <span class="keywordtype">int</span> num_cols, <span class="keywordtype">int</span> num_rows,
860
<a name="l00835"></a>00835 std::vector<float>& height_buffer, <span class="keywordtype">double</span> *delta_height )
859
<a name="l00834"></a>00834
860
<a name="l00835"></a>00835 <span class="keywordtype">void</span> GeoidLibrary::initializeEGM84Geoid()
861
861
<a name="l00836"></a>00836 {
862
862
<a name="l00837"></a>00837 <span class="comment">/*</span>
863
<a name="l00838"></a>00838 <span class="comment"> * The private function bilinearInterpolate returns the height of the</span>
864
<a name="l00839"></a>00839 <span class="comment"> * WGS84 geoid above or below the WGS84 ellipsoid, at the</span>
865
<a name="l00840"></a>00840 <span class="comment"> * specified geodetic coordinates, using a grid of height adjustments</span>
866
<a name="l00841"></a>00841 <span class="comment"> * and the bilinear interpolation method.</span>
867
<a name="l00842"></a>00842 <span class="comment"> *</span>
868
<a name="l00843"></a>00843 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
869
<a name="l00844"></a>00844 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
870
<a name="l00845"></a>00845 <span class="comment"> * scale_factor : Grid scale factor (input)</span>
871
<a name="l00846"></a>00846 <span class="comment"> * num_cols : Number of columns in grid (input)</span>
872
<a name="l00847"></a>00847 <span class="comment"> * num_rows : Number of rows in grid (input)</span>
873
<a name="l00848"></a>00848 <span class="comment"> * height_buffer : Grid of height adjustments (input)</span>
874
<a name="l00849"></a>00849 <span class="comment"> * delta_height : Height Adjustment, in meters. (output)</span>
875
<a name="l00850"></a>00850 <span class="comment"> *</span>
876
<a name="l00851"></a>00851 <span class="comment"> */</span>
863
<a name="l00838"></a>00838 <span class="comment"> * The function initializeEGM84Geoid reads geoid separation data from a file in</span>
864
<a name="l00839"></a>00839 <span class="comment"> * the current directory and builds the geoid separation table from it.</span>
865
<a name="l00840"></a>00840 <span class="comment"> * If the separation file can not be found or accessed, an error code of</span>
866
<a name="l00841"></a>00841 <span class="comment"> * GEOID_FILE_OPEN_ERROR is returned, If the separation file is incomplete</span>
867
<a name="l00842"></a>00842 <span class="comment"> * or improperly formatted, an error code of GEOID_INITIALIZE_ERROR is returned,</span>
868
<a name="l00843"></a>00843 <span class="comment"> * otherwise GEOID_NO_ERROR is returned.</span>
869
<a name="l00844"></a>00844 <span class="comment"> */</span>
870
<a name="l00845"></a>00845
871
<a name="l00846"></a>00846 <span class="keywordtype">int</span> items_read = 0;
872
<a name="l00847"></a>00847 <span class="keywordtype">char</span>* file_name = 0;
873
<a name="l00848"></a>00848 <span class="keywordtype">char</span>* path_name = getenv( <span class="stringliteral">"MSPCCS_DATA"</span> );
874
<a name="l00849"></a>00849 <span class="keywordtype">long</span> elevations_read = 0;
875
<a name="l00850"></a>00850 <span class="keywordtype">long</span> num = 0;
876
<a name="l00851"></a>00851 FILE* geoid_height_file;
877
877
<a name="l00852"></a>00852
878
<a name="l00853"></a>00853 <span class="keywordtype">int</span> index;
879
<a name="l00854"></a>00854 <span class="keywordtype">int</span> post_x, post_y;
880
<a name="l00855"></a>00855 <span class="keywordtype">double</span> offset_x, offset_y;
881
<a name="l00856"></a>00856 <span class="keywordtype">double</span> delta_x, delta_y;
882
<a name="l00857"></a>00857 <span class="keywordtype">double</span> _1_minus_delta_x, _1_minus_delta_y;
883
<a name="l00858"></a>00858 <span class="keywordtype">double</span> latitude_dd, longitude_dd;
884
<a name="l00859"></a>00859 <span class="keywordtype">double</span> height_se, height_sw, height_ne, height_nw;
885
<a name="l00860"></a>00860 <span class="keywordtype">double</span> w_sw, w_se, w_ne, w_nw;
886
<a name="l00861"></a>00861 <span class="keywordtype">double</span> south_lat, west_lon;
887
<a name="l00862"></a>00862 <span class="keywordtype">int</span> end_index = 0;
888
<a name="l00863"></a>00863 <span class="keywordtype">int</span> max_index = num_rows * num_cols - 1;
889
<a name="l00864"></a>00864 <span class="keywordtype">char</span> errorStatus[50] = <span class="stringliteral">""</span>;
890
<a name="l00865"></a>00865
891
<a name="l00866"></a>00866 <span class="keywordflow">if</span>( ( latitude < -<a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) || ( latitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) )
892
<a name="l00867"></a>00867 { <span class="comment">/* Latitude out of range */</span>
893
<a name="l00868"></a>00868 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a7b9a89219866d2577632c4a2a29f21a7">ErrorMessages::latitude</a> );
894
<a name="l00869"></a>00869 }
895
<a name="l00870"></a>00870 <span class="keywordflow">if</span>( ( longitude < -<a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> ) || ( longitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#a820597b124334bf7fb85214114f4876d">TWO_PI</a> ) )
896
<a name="l00871"></a>00871 { <span class="comment">/* Longitude out of range */</span>
897
<a name="l00872"></a>00872 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a148e03a19b3a95cae62e48e757245d12">ErrorMessages::longitude</a> );
898
<a name="l00873"></a>00873 }
899
<a name="l00874"></a>00874
900
<a name="l00875"></a>00875 <span class="keywordflow">if</span>( strlen( errorStatus ) > 0)
901
<a name="l00876"></a>00876 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( errorStatus );
878
<a name="l00853"></a>00853 <a class="code" href="class_m_s_p_1_1_c_c_s_thread_lock.html">CCSThreadLock</a> lock(&mutex);
879
<a name="l00854"></a>00854
880
<a name="l00855"></a>00855 <span class="comment">/* Check the environment for a user provided path, else current directory; */</span>
881
<a name="l00856"></a>00856 <span class="comment">/* Build a File Name, including specified or default path: */</span>
882
<a name="l00857"></a>00857
883
<a name="l00858"></a>00858 <span class="keywordflow">if</span> (path_name != NULL)
884
<a name="l00859"></a>00859 {
885
<a name="l00860"></a>00860 file_name = <span class="keyword">new</span> <span class="keywordtype">char</span>[ strlen( path_name ) + 11 ];
886
<a name="l00861"></a>00861 strcpy( file_name, path_name );
887
<a name="l00862"></a>00862 strcat( file_name, <span class="stringliteral">"/"</span> );
888
<a name="l00863"></a>00863 }
889
<a name="l00864"></a>00864 <span class="keywordflow">else</span>
890
<a name="l00865"></a>00865 {
891
<a name="l00866"></a>00866 file_name =<span class="keyword">new</span> <span class="keywordtype">char</span>[ 21 ];
892
<a name="l00867"></a>00867 strcpy( file_name, <span class="stringliteral">"../../data/"</span> );
893
<a name="l00868"></a>00868 }
894
<a name="l00869"></a>00869 strcat( file_name, <span class="stringliteral">"egm84.grd"</span> );
895
<a name="l00870"></a>00870
896
<a name="l00871"></a>00871 <span class="comment">/* Open the File READONLY, or Return Error Condition: */</span>
897
<a name="l00872"></a>00872
898
<a name="l00873"></a>00873 <span class="keywordflow">if</span>( ( geoid_height_file = fopen( file_name, <span class="stringliteral">"rb"</span> ) ) == NULL )
899
<a name="l00874"></a>00874 {
900
<a name="l00875"></a>00875 <span class="keyword">delete</span> [] file_name;
901
<a name="l00876"></a>00876 file_name = 0;
902
902
<a name="l00877"></a>00877
903
<a name="l00878"></a>00878 latitude_dd = latitude * <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a>;
904
<a name="l00879"></a>00879 longitude_dd = longitude * _180_OVER_PI;
905
<a name="l00880"></a>00880
906
<a name="l00881"></a>00881 <span class="comment">/* Compute X and Y Offsets into Geoid Height Array: */</span>
907
<a name="l00882"></a>00882
908
<a name="l00883"></a>00883 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
909
<a name="l00884"></a>00884 {
910
<a name="l00885"></a>00885 offset_x = ( longitude_dd + 360.0 ) / scale_factor;
911
<a name="l00886"></a>00886 }
912
<a name="l00887"></a>00887 <span class="keywordflow">else</span>
913
<a name="l00888"></a>00888 {
914
<a name="l00889"></a>00889 offset_x = longitude_dd / scale_factor;
915
<a name="l00890"></a>00890 }
916
<a name="l00891"></a>00891 offset_y = ( 90 - latitude_dd ) / scale_factor;
917
<a name="l00892"></a>00892
918
<a name="l00893"></a>00893 <span class="comment">/* Find Four Nearest Geoid Height Cells for specified latitude, longitude; */</span>
919
<a name="l00894"></a>00894 <span class="comment">/* Assumes that (0,0) of Geoid Height Array is at Northwest corner: */</span>
903
<a name="l00878"></a>00878 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
904
<a name="l00879"></a>00879 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a40a2ec1d77ae260338773a023518d449">ErrorMessages::geoidFileOpenError</a> );
905
<a name="l00880"></a>00880 strcat( message, <span class="stringliteral">": egm84.grd\n"</span> );
906
<a name="l00881"></a>00881 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
907
<a name="l00882"></a>00882 }
908
<a name="l00883"></a>00883
909
<a name="l00884"></a>00884
910
<a name="l00885"></a>00885 <span class="comment">/* Extract elements from the file: */</span>
911
<a name="l00886"></a>00886 egm84GeoidList = <span class="keyword">new</span> <span class="keywordtype">float</span>[<a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>];
912
<a name="l00887"></a>00887 elevations_read = <a class="code" href="_geoid_library_8cpp.html#a86b4187eeab007da13367ff3ac238250">readBinary</a>(
913
<a name="l00888"></a>00888 egm84GeoidList, 4, <a class="code" href="_geoid_library_8cpp.html#affd9297851046ceb405e6fa44c81ddb9">EGM84_ELEVATIONS</a>, geoid_height_file );
914
<a name="l00889"></a>00889
915
<a name="l00890"></a>00890 fclose( geoid_height_file );
916
<a name="l00891"></a>00891
917
<a name="l00892"></a>00892 <span class="keyword">delete</span> [] file_name;
918
<a name="l00893"></a>00893 file_name = 0;
919
<a name="l00894"></a>00894 }
920
920
<a name="l00895"></a>00895
921
<a name="l00896"></a>00896 post_x = ( int )( offset_x );
922
<a name="l00897"></a>00897 <span class="keywordflow">if</span>( ( post_x + 1 ) == num_cols )
923
<a name="l00898"></a>00898 post_x--;
924
<a name="l00899"></a>00899 post_y = ( int )( offset_y + 1.0e-11 );
925
<a name="l00900"></a>00900 <span class="keywordflow">if</span>( ( post_y + 1 ) == num_rows )
926
<a name="l00901"></a>00901 post_y--;
927
<a name="l00902"></a>00902
928
<a name="l00903"></a>00903 <span class="comment">// NW Height</span>
929
<a name="l00904"></a>00904 index = post_y * num_cols + post_x;
930
<a name="l00905"></a>00905 <span class="keywordflow">if</span>( index < 0 )
931
<a name="l00906"></a>00906 height_nw = height_buffer[ 0 ];
932
<a name="l00907"></a>00907 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
933
<a name="l00908"></a>00908 height_nw = height_buffer[ max_index ];
934
<a name="l00909"></a>00909 <span class="keywordflow">else</span>
935
<a name="l00910"></a>00910 height_nw = height_buffer[ index ];
936
<a name="l00911"></a>00911 <span class="comment">// NE Height</span>
937
<a name="l00912"></a>00912 end_index = index + 1;
938
<a name="l00913"></a>00913 <span class="keywordflow">if</span>( end_index > max_index )
939
<a name="l00914"></a>00914 height_ne = height_buffer[ max_index ];
940
<a name="l00915"></a>00915 <span class="keywordflow">else</span>
941
<a name="l00916"></a>00916 height_ne = height_buffer[ end_index ];
942
<a name="l00917"></a>00917
943
<a name="l00918"></a>00918 <span class="comment">// SW Height</span>
944
<a name="l00919"></a>00919 index = ( post_y + 1 ) * num_cols + post_x;
945
<a name="l00920"></a>00920 <span class="keywordflow">if</span>( index < 0 )
946
<a name="l00921"></a>00921 height_sw = height_buffer[ 0 ];
947
<a name="l00922"></a>00922 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
948
<a name="l00923"></a>00923 height_sw = height_buffer[ max_index ];
949
<a name="l00924"></a>00924 <span class="keywordflow">else</span>
950
<a name="l00925"></a>00925 height_sw = height_buffer[ index ];
951
<a name="l00926"></a>00926 <span class="comment">// SE Height</span>
952
<a name="l00927"></a>00927 end_index = index + 1;
953
<a name="l00928"></a>00928 <span class="keywordflow">if</span>( end_index > max_index )
954
<a name="l00929"></a>00929 height_se = height_buffer[ max_index ];
955
<a name="l00930"></a>00930 <span class="keywordflow">else</span>
956
<a name="l00931"></a>00931 height_se = height_buffer[ end_index ];
957
<a name="l00932"></a>00932
958
<a name="l00933"></a>00933 west_lon = post_x * scale_factor;
959
<a name="l00934"></a>00934
960
<a name="l00935"></a>00935 <span class="comment">// North latitude - scale_factor</span>
961
<a name="l00936"></a>00936 south_lat = ( 90 - ( post_y * scale_factor ) ) - scale_factor;
962
<a name="l00937"></a>00937
963
<a name="l00938"></a>00938 <span class="comment">/* Perform Bi-Linear Interpolation to compute Height above Ellipsoid: */</span>
964
<a name="l00939"></a>00939
965
<a name="l00940"></a>00940 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
966
<a name="l00941"></a>00941 delta_x = ( longitude_dd + 360.0 - west_lon ) / scale_factor;
967
<a name="l00942"></a>00942 <span class="keywordflow">else</span>
968
<a name="l00943"></a>00943 delta_x = ( longitude_dd - west_lon ) / scale_factor;
969
<a name="l00944"></a>00944 delta_y = ( latitude_dd - south_lat ) / scale_factor;
921
<a name="l00896"></a>00896 <span class="keywordtype">void</span> GeoidLibrary::initializeEGM84ThirtyMinGeoid()
922
<a name="l00897"></a>00897 {
923
<a name="l00898"></a>00898 <span class="comment">/*</span>
924
<a name="l00899"></a>00899 <span class="comment"> * The function initializeEGM84ThirtyMinGeoid reads geoid separation data from </span>
925
<a name="l00900"></a>00900 <span class="comment"> * a file in the current directory and builds the geoid separation table from it.</span>
926
<a name="l00901"></a>00901 <span class="comment"> * If the separation file can not be found or accessed, an error code of</span>
927
<a name="l00902"></a>00902 <span class="comment"> * GEOID_FILE_OPEN_ERROR is returned, If the separation file is incomplete</span>
928
<a name="l00903"></a>00903 <span class="comment"> * or improperly formatted, an error code of GEOID_INITIALIZE_ERROR is returned,</span>
929
<a name="l00904"></a>00904 <span class="comment"> * otherwise GEOID_NO_ERROR is returned.</span>
930
<a name="l00905"></a>00905 <span class="comment"> */</span>
931
<a name="l00906"></a>00906
932
<a name="l00907"></a>00907 <span class="keywordtype">int</span> items_read = 0;
933
<a name="l00908"></a>00908 <span class="keywordtype">char</span>* file_name = 0;
934
<a name="l00909"></a>00909 <span class="keywordtype">char</span>* path_name = getenv( <span class="stringliteral">"MSPCCS_DATA"</span> );
935
<a name="l00910"></a>00910 <span class="keywordtype">long</span> elevations_read = 0;
936
<a name="l00911"></a>00911 <span class="keywordtype">long</span> num = 0;
937
<a name="l00912"></a>00912 FILE* geoid_height_file;
938
<a name="l00913"></a>00913
939
<a name="l00914"></a>00914 <a class="code" href="class_m_s_p_1_1_c_c_s_thread_lock.html">CCSThreadLock</a> lock(&mutex);
940
<a name="l00915"></a>00915
941
<a name="l00916"></a>00916 <span class="comment">/* Check the environment for a user provided path, else current directory; */</span>
942
<a name="l00917"></a>00917 <span class="comment">/* Build a File Name, including specified or default path: */</span>
943
<a name="l00918"></a>00918
944
<a name="l00919"></a>00919 <span class="keywordflow">if</span> (path_name != NULL)
945
<a name="l00920"></a>00920 {
946
<a name="l00921"></a>00921 file_name = <span class="keyword">new</span> <span class="keywordtype">char</span>[ strlen( path_name ) + 12 ];
947
<a name="l00922"></a>00922 strcpy( file_name, path_name );
948
<a name="l00923"></a>00923 strcat( file_name, <span class="stringliteral">"/"</span> );
949
<a name="l00924"></a>00924 }
950
<a name="l00925"></a>00925 <span class="keywordflow">else</span>
951
<a name="l00926"></a>00926 {
952
<a name="l00927"></a>00927 file_name =<span class="keyword">new</span> <span class="keywordtype">char</span>[ 22 ];
953
<a name="l00928"></a>00928 strcpy( file_name, <span class="stringliteral">"../../data/"</span> );
954
<a name="l00929"></a>00929 }
955
<a name="l00930"></a>00930 strcat( file_name, <span class="stringliteral">"wwgrid.bin"</span> );
956
<a name="l00931"></a>00931
957
<a name="l00932"></a>00932 <span class="comment">/* Open the File READONLY, or Return Error Condition: */</span>
958
<a name="l00933"></a>00933
959
<a name="l00934"></a>00934 <span class="keywordflow">if</span>( ( geoid_height_file = fopen( file_name, <span class="stringliteral">"rb"</span> ) ) == NULL )
960
<a name="l00935"></a>00935 {
961
<a name="l00936"></a>00936 <span class="keyword">delete</span> [] file_name;
962
<a name="l00937"></a>00937 file_name = 0;
963
<a name="l00938"></a>00938
964
<a name="l00939"></a>00939 <span class="keywordtype">char</span> message[256] = <span class="stringliteral">""</span>;
965
<a name="l00940"></a>00940 strcpy( message, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a40a2ec1d77ae260338773a023518d449">ErrorMessages::geoidFileOpenError</a> );
966
<a name="l00941"></a>00941 strcat( message, <span class="stringliteral">": wwgrid.bin\n"</span> );
967
<a name="l00942"></a>00942 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( message );
968
<a name="l00943"></a>00943 }
969
<a name="l00944"></a>00944
970
970
<a name="l00945"></a>00945
971
<a name="l00946"></a>00946 _1_minus_delta_x = 1 - delta_x;
972
<a name="l00947"></a>00947 _1_minus_delta_y = 1 - delta_y;
973
<a name="l00948"></a>00948
974
<a name="l00949"></a>00949 w_sw = _1_minus_delta_x * _1_minus_delta_y;
975
<a name="l00950"></a>00950 w_se = delta_x * _1_minus_delta_y;
976
<a name="l00951"></a>00951 w_ne = delta_x * delta_y;
977
<a name="l00952"></a>00952 w_nw = _1_minus_delta_x * delta_y;
971
<a name="l00946"></a>00946 <span class="comment">/* Extract elements from the file: */</span>
972
<a name="l00947"></a>00947
973
<a name="l00948"></a>00948 egm84ThirtyMinGeoidList = <span class="keyword">new</span> <span class="keywordtype">double</span>[<a class="code" href="_geoid_library_8cpp.html#aa3581e7c236f1d903b4b3ccf242d1088">EGM84_30_MIN_ELEVATIONS</a>];
974
<a name="l00949"></a>00949 elevations_read = <a class="code" href="_geoid_library_8cpp.html#a86b4187eeab007da13367ff3ac238250">readBinary</a>(
975
<a name="l00950"></a>00950 egm84ThirtyMinGeoidList, 8, <a class="code" href="_geoid_library_8cpp.html#aa3581e7c236f1d903b4b3ccf242d1088">EGM84_30_MIN_ELEVATIONS</a>, geoid_height_file );
976
<a name="l00951"></a>00951
977
<a name="l00952"></a>00952 fclose( geoid_height_file );
978
978
<a name="l00953"></a>00953
979
<a name="l00954"></a>00954 *delta_height = height_sw * w_sw + height_se * w_se + height_ne * w_ne + height_nw * w_nw;
980
<a name="l00955"></a>00955 }
981
<a name="l00956"></a>00956
979
<a name="l00954"></a>00954 <span class="keyword">delete</span> [] file_name;
980
<a name="l00955"></a>00955 file_name = 0;
981
<a name="l00956"></a>00956 }
982
982
<a name="l00957"></a>00957
983
<a name="l00958"></a>00958 <span class="keywordtype">void</span> GeoidLibrary::naturalSplineInterpolate( <span class="keywordtype">double</span> longitude, <span class="keywordtype">double</span> latitude, <span class="keywordtype">double</span> scale_factor, <span class="keywordtype">int</span> num_cols, <span class="keywordtype">int</span> num_rows,
984
<a name="l00959"></a>00959 <span class="keywordtype">int</span> max_index, std::vector<float>& height_buffer, <span class="keywordtype">double</span> *delta_height )
985
<a name="l00960"></a>00960 {
986
<a name="l00961"></a>00961 <span class="comment">/*</span>
987
<a name="l00962"></a>00962 <span class="comment"> * The private function naturalSplineInterpolate returns the height of the</span>
988
<a name="l00963"></a>00963 <span class="comment"> * WGS84 geoid above or below the WGS84 ellipsoid, at the</span>
989
<a name="l00964"></a>00964 <span class="comment"> * specified geodetic coordinates, using a grid of height adjustments</span>
990
<a name="l00965"></a>00965 <span class="comment"> * and the natural spline interpolation method.</span>
991
<a name="l00966"></a>00966 <span class="comment"> *</span>
992
<a name="l00967"></a>00967 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
993
<a name="l00968"></a>00968 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
994
<a name="l00969"></a>00969 <span class="comment"> * scale_factor : Grid scale factor (input)</span>
995
<a name="l00970"></a>00970 <span class="comment"> * num_cols : Number of columns in grid (input)</span>
996
<a name="l00971"></a>00971 <span class="comment"> * num_rows : Number of rows in grid (input)</span>
997
<a name="l00972"></a>00972 <span class="comment"> * height_buffer : Grid of height adjustments (input)</span>
998
<a name="l00973"></a>00973 <span class="comment"> * DeltaHeight : Height Adjustment, in meters. (output)</span>
999
<a name="l00974"></a>00974 <span class="comment"> *</span>
1000
<a name="l00975"></a>00975 <span class="comment"> */</span>
1001
<a name="l00976"></a>00976
1002
<a name="l00977"></a>00977 <span class="keywordtype">int</span> index;
1003
<a name="l00978"></a>00978 <span class="keywordtype">int</span> post_x, post_y;
1004
<a name="l00979"></a>00979 <span class="keywordtype">int</span> temp_offset_x, temp_offset_y;
1005
<a name="l00980"></a>00980 <span class="keywordtype">double</span> offset_x, offset_y;
1006
<a name="l00981"></a>00981 <span class="keywordtype">double</span> delta_x, delta_y;
1007
<a name="l00982"></a>00982 <span class="keywordtype">double</span> delta_x2, delta_y2;
1008
<a name="l00983"></a>00983 <span class="keywordtype">double</span> _1_minus_delta_x, _1_minus_delta_y;
1009
<a name="l00984"></a>00984 <span class="keywordtype">double</span> _1_minus_delta_x2, _1_minus_delta_y2;
1010
<a name="l00985"></a>00985 <span class="keywordtype">double</span> _3_minus_2_times_1_minus_delta_x, _3_minus_2_times_1_minus_delta_y;
1011
<a name="l00986"></a>00986 <span class="keywordtype">double</span> _3_minus_2_times_delta_x, _3_minus_2_times_delta_y;
1012
<a name="l00987"></a>00987 <span class="keywordtype">double</span> latitude_dd, longitude_dd;
1013
<a name="l00988"></a>00988 <span class="keywordtype">double</span> height_se, height_sw, height_ne, height_nw;
1014
<a name="l00989"></a>00989 <span class="keywordtype">double</span> w_sw, w_se, w_ne, w_nw;
1015
<a name="l00990"></a>00990 <span class="keywordtype">double</span> south_lat, west_lon;
1016
<a name="l00991"></a>00991 <span class="keywordtype">int</span> end_index = 0;
1017
<a name="l00992"></a>00992 <span class="keywordtype">double</span> skip_factor = 1.0;
1018
<a name="l00993"></a>00993 <span class="keywordtype">char</span> errorStatus[50] = <span class="stringliteral">""</span>;
1019
<a name="l00994"></a>00994
1020
<a name="l00995"></a>00995 <span class="keywordflow">if</span>( ( latitude < -<a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) || ( latitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) )
1021
<a name="l00996"></a>00996 { <span class="comment">/* latitude out of range */</span>
1022
<a name="l00997"></a>00997 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a7b9a89219866d2577632c4a2a29f21a7">ErrorMessages::latitude</a> );
1023
<a name="l00998"></a>00998 }
1024
<a name="l00999"></a>00999 <span class="keywordflow">if</span>( ( longitude < -<a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> ) || ( longitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#a820597b124334bf7fb85214114f4876d">TWO_PI</a> ) )
1025
<a name="l01000"></a>01000 { <span class="comment">/* longitude out of range */</span>
1026
<a name="l01001"></a>01001 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a148e03a19b3a95cae62e48e757245d12">ErrorMessages::longitude</a> );
1027
<a name="l01002"></a>01002 }
1028
<a name="l01003"></a>01003
1029
<a name="l01004"></a>01004 <span class="keywordflow">if</span>( strlen( errorStatus ) > 0)
1030
<a name="l01005"></a>01005 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( errorStatus );
1031
<a name="l01006"></a>01006
1032
<a name="l01007"></a>01007 latitude_dd = latitude * _180_OVER_PI;
1033
<a name="l01008"></a>01008 longitude_dd = longitude * _180_OVER_PI;
1034
<a name="l01009"></a>01009
1035
<a name="l01010"></a>01010 <span class="comment">/* Compute X and Y Offsets into Geoid Height Array: */</span>
1036
<a name="l01011"></a>01011
1037
<a name="l01012"></a>01012 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1038
<a name="l01013"></a>01013 {
1039
<a name="l01014"></a>01014 offset_x = ( longitude_dd + 360.0 ) / scale_factor;
1040
<a name="l01015"></a>01015 }
1041
<a name="l01016"></a>01016 <span class="keywordflow">else</span>
1042
<a name="l01017"></a>01017 {
1043
<a name="l01018"></a>01018 offset_x = longitude_dd / scale_factor;
1044
<a name="l01019"></a>01019 }
1045
<a name="l01020"></a>01020 offset_y = ( 90.0 - latitude_dd ) / scale_factor;
1046
<a name="l01021"></a>01021
1047
<a name="l01022"></a>01022 <span class="comment">/* Find Four Nearest Geoid Height Cells for specified latitude, longitude; */</span>
1048
<a name="l01023"></a>01023 <span class="comment">/* Assumes that (0,0) of Geoid Height Array is at Northwest corner: */</span>
1049
<a name="l01024"></a>01024
1050
<a name="l01025"></a>01025 post_x = ( int ) offset_x;
1051
<a name="l01026"></a>01026 <span class="keywordflow">if</span> ( ( post_x + 1 ) == num_cols)
1052
<a name="l01027"></a>01027 post_x--;
1053
<a name="l01028"></a>01028 post_y = ( int )( offset_y + 1.0e-11 );
1054
<a name="l01029"></a>01029 <span class="keywordflow">if</span> ( ( post_y + 1 ) == num_rows)
1055
<a name="l01030"></a>01030 post_y--;
1056
<a name="l01031"></a>01031
1057
<a name="l01032"></a>01032 <span class="keywordflow">if</span>( scale_factor == <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a> )
1058
<a name="l01033"></a>01033 {
1059
<a name="l01034"></a>01034 skip_factor = 2.0;
1060
<a name="l01035"></a>01035 num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
1061
<a name="l01036"></a>01036 num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
1062
<a name="l01037"></a>01037 }
1063
<a name="l01038"></a>01038 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( scale_factor == <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a> )
1064
<a name="l01039"></a>01039 {
1065
<a name="l01040"></a>01040 skip_factor = 4.0;
1066
<a name="l01041"></a>01041 num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
1067
<a name="l01042"></a>01042 num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
1068
<a name="l01043"></a>01043 }
1069
<a name="l01044"></a>01044 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( scale_factor == <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a> )
1070
<a name="l01045"></a>01045 {
1071
<a name="l01046"></a>01046 skip_factor = 8.0;
1072
<a name="l01047"></a>01047 num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
1073
<a name="l01048"></a>01048 num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
1074
<a name="l01049"></a>01049 }
1075
<a name="l01050"></a>01050
1076
<a name="l01051"></a>01051 temp_offset_x = ( int )( post_x * skip_factor );
1077
<a name="l01052"></a>01052 temp_offset_y = ( int )( post_y * skip_factor + 1.0e-11 );
1078
<a name="l01053"></a>01053 <span class="keywordflow">if</span> ( ( temp_offset_x + 1 ) == num_cols )
1079
<a name="l01054"></a>01054 temp_offset_x--;
1080
<a name="l01055"></a>01055 <span class="keywordflow">if</span> ( ( temp_offset_y + 1 ) == num_rows )
1081
<a name="l01056"></a>01056 temp_offset_y--;
1082
<a name="l01057"></a>01057
1083
<a name="l01058"></a>01058 <span class="comment">// NW Height</span>
1084
<a name="l01059"></a>01059 index = ( int )( temp_offset_y * num_cols + temp_offset_x );
1085
<a name="l01060"></a>01060 <span class="keywordflow">if</span>( index < 0 )
1086
<a name="l01061"></a>01061 height_nw = height_buffer[ 0 ];
1087
<a name="l01062"></a>01062 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1088
<a name="l01063"></a>01063 height_nw = height_buffer[ max_index ];
1089
<a name="l01064"></a>01064 <span class="keywordflow">else</span>
1090
<a name="l01065"></a>01065 height_nw = height_buffer[ index ];
1091
<a name="l01066"></a>01066 <span class="comment">// NE Height</span>
1092
<a name="l01067"></a>01067 end_index = index + (int)skip_factor;
1093
<a name="l01068"></a>01068 <span class="keywordflow">if</span>( end_index < 0 )
1094
<a name="l01069"></a>01069 height_ne = height_buffer[ 0 ];
1095
<a name="l01070"></a>01070 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( end_index > max_index )
1096
<a name="l01071"></a>01071 height_ne = height_buffer[ max_index ];
1097
<a name="l01072"></a>01072 <span class="keywordflow">else</span>
1098
<a name="l01073"></a>01073 height_ne = height_buffer[ end_index ];
1099
<a name="l01074"></a>01074
1100
<a name="l01075"></a>01075 <span class="comment">// SW Height</span>
1101
<a name="l01076"></a>01076 index = ( int )( ( temp_offset_y + skip_factor ) * num_cols + temp_offset_x );
1102
<a name="l01077"></a>01077 <span class="keywordflow">if</span>( index < 0 )
1103
<a name="l01078"></a>01078 height_sw = height_buffer[ 0 ];
1104
<a name="l01079"></a>01079 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1105
<a name="l01080"></a>01080 height_sw = height_buffer[ max_index ];
1106
<a name="l01081"></a>01081 <span class="keywordflow">else</span>
1107
<a name="l01082"></a>01082 height_sw = height_buffer[ index ];
1108
<a name="l01083"></a>01083 <span class="comment">// SE Height</span>
1109
<a name="l01084"></a>01084 end_index = index + (int)skip_factor;
1110
<a name="l01085"></a>01085 <span class="keywordflow">if</span>( end_index < 0 )
1111
<a name="l01086"></a>01086 height_se = height_buffer[ 0 ];
1112
<a name="l01087"></a>01087 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( end_index > max_index )
1113
<a name="l01088"></a>01088 height_se = height_buffer[ max_index ];
1114
<a name="l01089"></a>01089 <span class="keywordflow">else</span>
1115
<a name="l01090"></a>01090 height_se = height_buffer[ end_index ];
1116
<a name="l01091"></a>01091
1117
<a name="l01092"></a>01092 west_lon = post_x * scale_factor;
1118
<a name="l01093"></a>01093
1119
<a name="l01094"></a>01094 <span class="comment">// North latitude - scale_factor</span>
1120
<a name="l01095"></a>01095 south_lat = ( 90 - ( post_y * scale_factor ) ) - scale_factor;
1121
<a name="l01096"></a>01096
1122
<a name="l01097"></a>01097 <span class="comment">/* Perform Non-Linear Interpolation to compute Height above Ellipsoid: */</span>
1123
<a name="l01098"></a>01098
1124
<a name="l01099"></a>01099 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1125
<a name="l01100"></a>01100 delta_x = ( longitude_dd + 360.0 - west_lon ) / scale_factor;
1126
<a name="l01101"></a>01101 <span class="keywordflow">else</span>
1127
<a name="l01102"></a>01102 delta_x = ( longitude_dd - west_lon ) / scale_factor;
1128
<a name="l01103"></a>01103 delta_y = ( latitude_dd - south_lat ) / scale_factor;
1129
<a name="l01104"></a>01104
1130
<a name="l01105"></a>01105 delta_x2 = delta_x * delta_x;
1131
<a name="l01106"></a>01106 delta_y2 = delta_y * delta_y;
1132
<a name="l01107"></a>01107
1133
<a name="l01108"></a>01108 _1_minus_delta_x = 1 - delta_x;
1134
<a name="l01109"></a>01109 _1_minus_delta_y = 1 - delta_y;
1135
<a name="l01110"></a>01110
1136
<a name="l01111"></a>01111 _1_minus_delta_x2 = _1_minus_delta_x * _1_minus_delta_x;
1137
<a name="l01112"></a>01112 _1_minus_delta_y2 = _1_minus_delta_y * _1_minus_delta_y;
1138
<a name="l01113"></a>01113
1139
<a name="l01114"></a>01114 _3_minus_2_times_1_minus_delta_x = 3 - 2 * _1_minus_delta_x;
1140
<a name="l01115"></a>01115 _3_minus_2_times_1_minus_delta_y = 3 - 2 * _1_minus_delta_y;
1141
<a name="l01116"></a>01116
1142
<a name="l01117"></a>01117 _3_minus_2_times_delta_x = 3 - 2 * delta_x;
1143
<a name="l01118"></a>01118 _3_minus_2_times_delta_y = 3 - 2 * delta_y;
1144
<a name="l01119"></a>01119
1145
<a name="l01120"></a>01120 w_sw = _1_minus_delta_x2 * _1_minus_delta_y2 * ( _3_minus_2_times_1_minus_delta_x * _3_minus_2_times_1_minus_delta_y );
1146
<a name="l01121"></a>01121 w_se = delta_x2 * _1_minus_delta_y2 * ( _3_minus_2_times_delta_x * _3_minus_2_times_1_minus_delta_y );
1147
<a name="l01122"></a>01122 w_ne = delta_x2 * delta_y2 * ( _3_minus_2_times_delta_x * _3_minus_2_times_delta_y );
1148
<a name="l01123"></a>01123 w_nw = _1_minus_delta_x2 * delta_y2 * ( _3_minus_2_times_1_minus_delta_x * _3_minus_2_times_delta_y );
1149
<a name="l01124"></a>01124
1150
<a name="l01125"></a>01125 *delta_height = height_sw * w_sw + height_se * w_se + height_ne * w_ne + height_nw * w_nw;
1151
<a name="l01126"></a>01126 }
983
<a name="l00958"></a>00958 <span class="keywordtype">void</span> GeoidLibrary::bilinearInterpolateDoubleHeights(
984
<a name="l00959"></a>00959 <span class="keywordtype">double</span> longitude,
985
<a name="l00960"></a>00960 <span class="keywordtype">double</span> latitude,
986
<a name="l00961"></a>00961 <span class="keywordtype">double</span> scale_factor,
987
<a name="l00962"></a>00962 <span class="keywordtype">int</span> num_cols,
988
<a name="l00963"></a>00963 <span class="keywordtype">int</span> num_rows,
989
<a name="l00964"></a>00964 <span class="keywordtype">double</span> *height_buffer,
990
<a name="l00965"></a>00965 <span class="keywordtype">double</span> *delta_height )
991
<a name="l00966"></a>00966 {
992
<a name="l00967"></a>00967 <span class="comment">/*</span>
993
<a name="l00968"></a>00968 <span class="comment"> * The private function bilinearInterpolateDoubleHeights returns the height of the</span>
994
<a name="l00969"></a>00969 <span class="comment"> * WGS84 geoid above or below the WGS84 ellipsoid, at the</span>
995
<a name="l00970"></a>00970 <span class="comment"> * specified geodetic coordinates, using a grid of height adjustments</span>
996
<a name="l00971"></a>00971 <span class="comment"> * and the bilinear interpolation method.</span>
997
<a name="l00972"></a>00972 <span class="comment"> *</span>
998
<a name="l00973"></a>00973 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
999
<a name="l00974"></a>00974 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
1000
<a name="l00975"></a>00975 <span class="comment"> * scale_factor : Grid scale factor (input)</span>
1001
<a name="l00976"></a>00976 <span class="comment"> * num_cols : Number of columns in grid (input)</span>
1002
<a name="l00977"></a>00977 <span class="comment"> * num_rows : Number of rows in grid (input)</span>
1003
<a name="l00978"></a>00978 <span class="comment"> * height_buffer : Grid of height adjustments, doubles (input)</span>
1004
<a name="l00979"></a>00979 <span class="comment"> * delta_height : Height Adjustment, in meters. (output)</span>
1005
<a name="l00980"></a>00980 <span class="comment"> *</span>
1006
<a name="l00981"></a>00981 <span class="comment"> */</span>
1007
<a name="l00982"></a>00982
1008
<a name="l00983"></a>00983 <span class="keywordtype">int</span> index;
1009
<a name="l00984"></a>00984 <span class="keywordtype">int</span> post_x, post_y;
1010
<a name="l00985"></a>00985 <span class="keywordtype">double</span> offset_x, offset_y;
1011
<a name="l00986"></a>00986 <span class="keywordtype">double</span> delta_x, delta_y;
1012
<a name="l00987"></a>00987 <span class="keywordtype">double</span> _1_minus_delta_x, _1_minus_delta_y;
1013
<a name="l00988"></a>00988 <span class="keywordtype">double</span> latitude_dd, longitude_dd;
1014
<a name="l00989"></a>00989 <span class="keywordtype">double</span> height_se, height_sw, height_ne, height_nw;
1015
<a name="l00990"></a>00990 <span class="keywordtype">double</span> w_sw, w_se, w_ne, w_nw;
1016
<a name="l00991"></a>00991 <span class="keywordtype">double</span> south_lat, west_lon;
1017
<a name="l00992"></a>00992 <span class="keywordtype">int</span> end_index = 0;
1018
<a name="l00993"></a>00993 <span class="keywordtype">int</span> max_index = num_rows * num_cols - 1;
1019
<a name="l00994"></a>00994 <span class="keywordtype">char</span> errorStatus[50] = <span class="stringliteral">""</span>;
1020
<a name="l00995"></a>00995
1021
<a name="l00996"></a>00996 <span class="keywordflow">if</span>( ( latitude < -<a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) || ( latitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) )
1022
<a name="l00997"></a>00997 { <span class="comment">/* Latitude out of range */</span>
1023
<a name="l00998"></a>00998 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a7b9a89219866d2577632c4a2a29f21a7">ErrorMessages::latitude</a> );
1024
<a name="l00999"></a>00999 }
1025
<a name="l01000"></a>01000 <span class="keywordflow">if</span>( ( longitude < -<a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> ) || ( longitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#a820597b124334bf7fb85214114f4876d">TWO_PI</a> ) )
1026
<a name="l01001"></a>01001 { <span class="comment">/* Longitude out of range */</span>
1027
<a name="l01002"></a>01002 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a148e03a19b3a95cae62e48e757245d12">ErrorMessages::longitude</a> );
1028
<a name="l01003"></a>01003 }
1029
<a name="l01004"></a>01004
1030
<a name="l01005"></a>01005 <span class="keywordflow">if</span>( strlen( errorStatus ) > 0)
1031
<a name="l01006"></a>01006 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( errorStatus );
1032
<a name="l01007"></a>01007
1033
<a name="l01008"></a>01008 latitude_dd = latitude * <a class="code" href="_datum_library_implementation_8cpp.html#a673bac0558e3a6ca72b12a4f7f84bef2">_180_OVER_PI</a>;
1034
<a name="l01009"></a>01009 longitude_dd = longitude * _180_OVER_PI;
1035
<a name="l01010"></a>01010
1036
<a name="l01011"></a>01011 <span class="comment">/* Compute X and Y Offsets into Geoid Height Array: */</span>
1037
<a name="l01012"></a>01012
1038
<a name="l01013"></a>01013 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1039
<a name="l01014"></a>01014 {
1040
<a name="l01015"></a>01015 offset_x = ( longitude_dd + 360.0 ) / scale_factor;
1041
<a name="l01016"></a>01016 }
1042
<a name="l01017"></a>01017 <span class="keywordflow">else</span>
1043
<a name="l01018"></a>01018 {
1044
<a name="l01019"></a>01019 offset_x = longitude_dd / scale_factor;
1045
<a name="l01020"></a>01020 }
1046
<a name="l01021"></a>01021 offset_y = ( 90 - latitude_dd ) / scale_factor;
1047
<a name="l01022"></a>01022
1048
<a name="l01023"></a>01023 <span class="comment">/* Find Four Nearest Geoid Height Cells for specified latitude, longitude; */</span>
1049
<a name="l01024"></a>01024 <span class="comment">/* Assumes that (0,0) of Geoid Height Array is at Northwest corner: */</span>
1050
<a name="l01025"></a>01025
1051
<a name="l01026"></a>01026 post_x = ( int )( offset_x );
1052
<a name="l01027"></a>01027 <span class="keywordflow">if</span>( ( post_x + 1 ) == num_cols )
1053
<a name="l01028"></a>01028 post_x--;
1054
<a name="l01029"></a>01029 post_y = ( int )( offset_y + 1.0e-11 );
1055
<a name="l01030"></a>01030 <span class="keywordflow">if</span>( ( post_y + 1 ) == num_rows )
1056
<a name="l01031"></a>01031 post_y--;
1057
<a name="l01032"></a>01032
1058
<a name="l01033"></a>01033 <span class="comment">// NW Height</span>
1059
<a name="l01034"></a>01034 index = post_y * num_cols + post_x;
1060
<a name="l01035"></a>01035 <span class="keywordflow">if</span>( index < 0 )
1061
<a name="l01036"></a>01036 height_nw = height_buffer[ 0 ];
1062
<a name="l01037"></a>01037 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1063
<a name="l01038"></a>01038 height_nw = height_buffer[ max_index ];
1064
<a name="l01039"></a>01039 <span class="keywordflow">else</span>
1065
<a name="l01040"></a>01040 height_nw = height_buffer[ index ];
1066
<a name="l01041"></a>01041 <span class="comment">// NE Height</span>
1067
<a name="l01042"></a>01042 end_index = index + 1;
1068
<a name="l01043"></a>01043 <span class="keywordflow">if</span>( end_index > max_index )
1069
<a name="l01044"></a>01044 height_ne = height_buffer[ max_index ];
1070
<a name="l01045"></a>01045 <span class="keywordflow">else</span>
1071
<a name="l01046"></a>01046 height_ne = height_buffer[ end_index ];
1072
<a name="l01047"></a>01047
1073
<a name="l01048"></a>01048 <span class="comment">// SW Height</span>
1074
<a name="l01049"></a>01049 index = ( post_y + 1 ) * num_cols + post_x;
1075
<a name="l01050"></a>01050 <span class="keywordflow">if</span>( index < 0 )
1076
<a name="l01051"></a>01051 height_sw = height_buffer[ 0 ];
1077
<a name="l01052"></a>01052 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1078
<a name="l01053"></a>01053 height_sw = height_buffer[ max_index ];
1079
<a name="l01054"></a>01054 <span class="keywordflow">else</span>
1080
<a name="l01055"></a>01055 height_sw = height_buffer[ index ];
1081
<a name="l01056"></a>01056
1082
<a name="l01057"></a>01057 <span class="comment">// SE Height</span>
1083
<a name="l01058"></a>01058 end_index = index + 1;
1084
<a name="l01059"></a>01059 <span class="keywordflow">if</span>( end_index > max_index )
1085
<a name="l01060"></a>01060 height_se = height_buffer[ max_index ];
1086
<a name="l01061"></a>01061 <span class="keywordflow">else</span>
1087
<a name="l01062"></a>01062 height_se = height_buffer[ end_index ];
1088
<a name="l01063"></a>01063
1089
<a name="l01064"></a>01064 west_lon = post_x * scale_factor;
1090
<a name="l01065"></a>01065
1091
<a name="l01066"></a>01066 <span class="comment">// North latitude - scale_factor</span>
1092
<a name="l01067"></a>01067 south_lat = ( 90 - ( post_y * scale_factor ) ) - scale_factor;
1093
<a name="l01068"></a>01068
1094
<a name="l01069"></a>01069 <span class="comment">/* Perform Bi-Linear Interpolation to compute Height above Ellipsoid: */</span>
1095
<a name="l01070"></a>01070
1096
<a name="l01071"></a>01071 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1097
<a name="l01072"></a>01072 delta_x = ( longitude_dd + 360.0 - west_lon ) / scale_factor;
1098
<a name="l01073"></a>01073 <span class="keywordflow">else</span>
1099
<a name="l01074"></a>01074 delta_x = ( longitude_dd - west_lon ) / scale_factor;
1100
<a name="l01075"></a>01075 delta_y = ( latitude_dd - south_lat ) / scale_factor;
1101
<a name="l01076"></a>01076
1102
<a name="l01077"></a>01077 _1_minus_delta_x = 1 - delta_x;
1103
<a name="l01078"></a>01078 _1_minus_delta_y = 1 - delta_y;
1104
<a name="l01079"></a>01079
1105
<a name="l01080"></a>01080 w_sw = _1_minus_delta_x * _1_minus_delta_y;
1106
<a name="l01081"></a>01081 w_se = delta_x * _1_minus_delta_y;
1107
<a name="l01082"></a>01082 w_ne = delta_x * delta_y;
1108
<a name="l01083"></a>01083 w_nw = _1_minus_delta_x * delta_y;
1109
<a name="l01084"></a>01084
1110
<a name="l01085"></a>01085 *delta_height =
1111
<a name="l01086"></a>01086 height_sw * w_sw + height_se * w_se + height_ne * w_ne + height_nw * w_nw;
1112
<a name="l01087"></a>01087 }
1113
<a name="l01088"></a>01088
1114
<a name="l01089"></a>01089
1115
<a name="l01090"></a>01090 <span class="keywordtype">void</span> GeoidLibrary::bilinearInterpolate(
1116
<a name="l01091"></a>01091 <span class="keywordtype">double</span> longitude,
1117
<a name="l01092"></a>01092 <span class="keywordtype">double</span> latitude,
1118
<a name="l01093"></a>01093 <span class="keywordtype">double</span> scale_factor,
1119
<a name="l01094"></a>01094 <span class="keywordtype">int</span> num_cols,
1120
<a name="l01095"></a>01095 <span class="keywordtype">int</span> num_rows,
1121
<a name="l01096"></a>01096 <span class="keywordtype">float</span> *height_buffer,
1122
<a name="l01097"></a>01097 <span class="keywordtype">double</span> *delta_height )
1123
<a name="l01098"></a>01098 {
1124
<a name="l01099"></a>01099 <span class="comment">/*</span>
1125
<a name="l01100"></a>01100 <span class="comment"> * The private function bilinearInterpolate returns the height of the</span>
1126
<a name="l01101"></a>01101 <span class="comment"> * WGS84 geoid above or below the WGS84 ellipsoid, at the</span>
1127
<a name="l01102"></a>01102 <span class="comment"> * specified geodetic coordinates, using a grid of height adjustments</span>
1128
<a name="l01103"></a>01103 <span class="comment"> * and the bilinear interpolation method.</span>
1129
<a name="l01104"></a>01104 <span class="comment"> *</span>
1130
<a name="l01105"></a>01105 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
1131
<a name="l01106"></a>01106 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
1132
<a name="l01107"></a>01107 <span class="comment"> * scale_factor : Grid scale factor (input)</span>
1133
<a name="l01108"></a>01108 <span class="comment"> * num_cols : Number of columns in grid (input)</span>
1134
<a name="l01109"></a>01109 <span class="comment"> * num_rows : Number of rows in grid (input)</span>
1135
<a name="l01110"></a>01110 <span class="comment"> * height_buffer : Grid of height adjustments, floats (input)</span>
1136
<a name="l01111"></a>01111 <span class="comment"> * delta_height : Height Adjustment, in meters. (output)</span>
1137
<a name="l01112"></a>01112 <span class="comment"> *</span>
1138
<a name="l01113"></a>01113 <span class="comment"> */</span>
1139
<a name="l01114"></a>01114
1140
<a name="l01115"></a>01115 <span class="keywordtype">int</span> index;
1141
<a name="l01116"></a>01116 <span class="keywordtype">int</span> post_x, post_y;
1142
<a name="l01117"></a>01117 <span class="keywordtype">double</span> offset_x, offset_y;
1143
<a name="l01118"></a>01118 <span class="keywordtype">double</span> delta_x, delta_y;
1144
<a name="l01119"></a>01119 <span class="keywordtype">double</span> _1_minus_delta_x, _1_minus_delta_y;
1145
<a name="l01120"></a>01120 <span class="keywordtype">double</span> latitude_dd, longitude_dd;
1146
<a name="l01121"></a>01121 <span class="keywordtype">double</span> height_se, height_sw, height_ne, height_nw;
1147
<a name="l01122"></a>01122 <span class="keywordtype">double</span> w_sw, w_se, w_ne, w_nw;
1148
<a name="l01123"></a>01123 <span class="keywordtype">double</span> south_lat, west_lon;
1149
<a name="l01124"></a>01124 <span class="keywordtype">int</span> end_index = 0;
1150
<a name="l01125"></a>01125 <span class="keywordtype">int</span> max_index = num_rows * num_cols - 1;
1151
<a name="l01126"></a>01126 <span class="keywordtype">char</span> errorStatus[50] = <span class="stringliteral">""</span>;
1152
1152
<a name="l01127"></a>01127
1153
<a name="l01128"></a>01128 <span class="comment">// CLASSIFICATION: UNCLASSIFIED</span>
1153
<a name="l01128"></a>01128 <span class="keywordflow">if</span>( ( latitude < -<a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) || ( latitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) )
1154
<a name="l01129"></a>01129 { <span class="comment">/* Latitude out of range */</span>
1155
<a name="l01130"></a>01130 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a7b9a89219866d2577632c4a2a29f21a7">ErrorMessages::latitude</a> );
1156
<a name="l01131"></a>01131 }
1157
<a name="l01132"></a>01132 <span class="keywordflow">if</span>( ( longitude < -<a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> ) || ( longitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#a820597b124334bf7fb85214114f4876d">TWO_PI</a> ) )
1158
<a name="l01133"></a>01133 { <span class="comment">/* Longitude out of range */</span>
1159
<a name="l01134"></a>01134 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a148e03a19b3a95cae62e48e757245d12">ErrorMessages::longitude</a> );
1160
<a name="l01135"></a>01135 }
1161
<a name="l01136"></a>01136
1162
<a name="l01137"></a>01137 <span class="keywordflow">if</span>( strlen( errorStatus ) > 0)
1163
<a name="l01138"></a>01138 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( errorStatus );
1164
<a name="l01139"></a>01139
1165
<a name="l01140"></a>01140 latitude_dd = latitude * _180_OVER_PI;
1166
<a name="l01141"></a>01141 longitude_dd = longitude * _180_OVER_PI;
1167
<a name="l01142"></a>01142
1168
<a name="l01143"></a>01143 <span class="comment">/* Compute X and Y Offsets into Geoid Height Array: */</span>
1169
<a name="l01144"></a>01144
1170
<a name="l01145"></a>01145 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1171
<a name="l01146"></a>01146 {
1172
<a name="l01147"></a>01147 offset_x = ( longitude_dd + 360.0 ) / scale_factor;
1173
<a name="l01148"></a>01148 }
1174
<a name="l01149"></a>01149 <span class="keywordflow">else</span>
1175
<a name="l01150"></a>01150 {
1176
<a name="l01151"></a>01151 offset_x = longitude_dd / scale_factor;
1177
<a name="l01152"></a>01152 }
1178
<a name="l01153"></a>01153 offset_y = ( 90 - latitude_dd ) / scale_factor;
1179
<a name="l01154"></a>01154
1180
<a name="l01155"></a>01155 <span class="comment">/* Find Four Nearest Geoid Height Cells for specified latitude, longitude; */</span>
1181
<a name="l01156"></a>01156 <span class="comment">/* Assumes that (0,0) of Geoid Height Array is at Northwest corner: */</span>
1182
<a name="l01157"></a>01157
1183
<a name="l01158"></a>01158 post_x = ( int )( offset_x );
1184
<a name="l01159"></a>01159 <span class="keywordflow">if</span>( ( post_x + 1 ) == num_cols )
1185
<a name="l01160"></a>01160 post_x--;
1186
<a name="l01161"></a>01161 post_y = ( int )( offset_y + 1.0e-11 );
1187
<a name="l01162"></a>01162 <span class="keywordflow">if</span>( ( post_y + 1 ) == num_rows )
1188
<a name="l01163"></a>01163 post_y--;
1189
<a name="l01164"></a>01164
1190
<a name="l01165"></a>01165 <span class="comment">// NW Height</span>
1191
<a name="l01166"></a>01166 index = post_y * num_cols + post_x;
1192
<a name="l01167"></a>01167 <span class="keywordflow">if</span>( index < 0 )
1193
<a name="l01168"></a>01168 height_nw = height_buffer[ 0 ];
1194
<a name="l01169"></a>01169 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1195
<a name="l01170"></a>01170 height_nw = height_buffer[ max_index ];
1196
<a name="l01171"></a>01171 <span class="keywordflow">else</span>
1197
<a name="l01172"></a>01172 height_nw = height_buffer[ index ];
1198
<a name="l01173"></a>01173 <span class="comment">// NE Height</span>
1199
<a name="l01174"></a>01174 end_index = index + 1;
1200
<a name="l01175"></a>01175 <span class="keywordflow">if</span>( end_index > max_index )
1201
<a name="l01176"></a>01176 height_ne = height_buffer[ max_index ];
1202
<a name="l01177"></a>01177 <span class="keywordflow">else</span>
1203
<a name="l01178"></a>01178 height_ne = height_buffer[ end_index ];
1204
<a name="l01179"></a>01179
1205
<a name="l01180"></a>01180 <span class="comment">// SW Height</span>
1206
<a name="l01181"></a>01181 index = ( post_y + 1 ) * num_cols + post_x;
1207
<a name="l01182"></a>01182 <span class="keywordflow">if</span>( index < 0 )
1208
<a name="l01183"></a>01183 height_sw = height_buffer[ 0 ];
1209
<a name="l01184"></a>01184 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1210
<a name="l01185"></a>01185 height_sw = height_buffer[ max_index ];
1211
<a name="l01186"></a>01186 <span class="keywordflow">else</span>
1212
<a name="l01187"></a>01187 height_sw = height_buffer[ index ];
1213
<a name="l01188"></a>01188
1214
<a name="l01189"></a>01189 <span class="comment">// SE Height</span>
1215
<a name="l01190"></a>01190 end_index = index + 1;
1216
<a name="l01191"></a>01191 <span class="keywordflow">if</span>( end_index > max_index )
1217
<a name="l01192"></a>01192 height_se = height_buffer[ max_index ];
1218
<a name="l01193"></a>01193 <span class="keywordflow">else</span>
1219
<a name="l01194"></a>01194 height_se = height_buffer[ end_index ];
1220
<a name="l01195"></a>01195
1221
<a name="l01196"></a>01196 west_lon = post_x * scale_factor;
1222
<a name="l01197"></a>01197
1223
<a name="l01198"></a>01198 <span class="comment">// North latitude - scale_factor</span>
1224
<a name="l01199"></a>01199 south_lat = ( 90 - ( post_y * scale_factor ) ) - scale_factor;
1225
<a name="l01200"></a>01200
1226
<a name="l01201"></a>01201 <span class="comment">/* Perform Bi-Linear Interpolation to compute Height above Ellipsoid: */</span>
1227
<a name="l01202"></a>01202
1228
<a name="l01203"></a>01203 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1229
<a name="l01204"></a>01204 delta_x = ( longitude_dd + 360.0 - west_lon ) / scale_factor;
1230
<a name="l01205"></a>01205 <span class="keywordflow">else</span>
1231
<a name="l01206"></a>01206 delta_x = ( longitude_dd - west_lon ) / scale_factor;
1232
<a name="l01207"></a>01207 delta_y = ( latitude_dd - south_lat ) / scale_factor;
1233
<a name="l01208"></a>01208
1234
<a name="l01209"></a>01209 _1_minus_delta_x = 1 - delta_x;
1235
<a name="l01210"></a>01210 _1_minus_delta_y = 1 - delta_y;
1236
<a name="l01211"></a>01211
1237
<a name="l01212"></a>01212 w_sw = _1_minus_delta_x * _1_minus_delta_y;
1238
<a name="l01213"></a>01213 w_se = delta_x * _1_minus_delta_y;
1239
<a name="l01214"></a>01214 w_ne = delta_x * delta_y;
1240
<a name="l01215"></a>01215 w_nw = _1_minus_delta_x * delta_y;
1241
<a name="l01216"></a>01216
1242
<a name="l01217"></a>01217 *delta_height =
1243
<a name="l01218"></a>01218 height_sw * w_sw + height_se * w_se + height_ne * w_ne + height_nw * w_nw;
1244
<a name="l01219"></a>01219 }
1245
<a name="l01220"></a>01220
1246
<a name="l01221"></a>01221
1247
<a name="l01222"></a>01222 <span class="keywordtype">void</span> GeoidLibrary::naturalSplineInterpolate(
1248
<a name="l01223"></a>01223 <span class="keywordtype">double</span> longitude,
1249
<a name="l01224"></a>01224 <span class="keywordtype">double</span> latitude,
1250
<a name="l01225"></a>01225 <span class="keywordtype">double</span> scale_factor,
1251
<a name="l01226"></a>01226 <span class="keywordtype">int</span> num_cols,
1252
<a name="l01227"></a>01227 <span class="keywordtype">int</span> num_rows,
1253
<a name="l01228"></a>01228 <span class="keywordtype">int</span> max_index,
1254
<a name="l01229"></a>01229 <span class="keywordtype">float</span> *height_buffer,
1255
<a name="l01230"></a>01230 <span class="keywordtype">double</span> *delta_height )
1256
<a name="l01231"></a>01231 {
1257
<a name="l01232"></a>01232 <span class="comment">/*</span>
1258
<a name="l01233"></a>01233 <span class="comment"> * The private function naturalSplineInterpolate returns the height of the</span>
1259
<a name="l01234"></a>01234 <span class="comment"> * WGS84 geoid above or below the WGS84 ellipsoid, at the</span>
1260
<a name="l01235"></a>01235 <span class="comment"> * specified geodetic coordinates, using a grid of height adjustments</span>
1261
<a name="l01236"></a>01236 <span class="comment"> * and the natural spline interpolation method.</span>
1262
<a name="l01237"></a>01237 <span class="comment"> *</span>
1263
<a name="l01238"></a>01238 <span class="comment"> * longitude : Geodetic longitude in radians (input)</span>
1264
<a name="l01239"></a>01239 <span class="comment"> * latitude : Geodetic latitude in radians (input)</span>
1265
<a name="l01240"></a>01240 <span class="comment"> * scale_factor : Grid scale factor (input)</span>
1266
<a name="l01241"></a>01241 <span class="comment"> * num_cols : Number of columns in grid (input)</span>
1267
<a name="l01242"></a>01242 <span class="comment"> * num_rows : Number of rows in grid (input)</span>
1268
<a name="l01243"></a>01243 <span class="comment"> * height_buffer : Grid of height adjustments (input)</span>
1269
<a name="l01244"></a>01244 <span class="comment"> * DeltaHeight : Height Adjustment, in meters. (output)</span>
1270
<a name="l01245"></a>01245 <span class="comment"> *</span>
1271
<a name="l01246"></a>01246 <span class="comment"> */</span>
1272
<a name="l01247"></a>01247
1273
<a name="l01248"></a>01248 <span class="keywordtype">int</span> index;
1274
<a name="l01249"></a>01249 <span class="keywordtype">int</span> post_x, post_y;
1275
<a name="l01250"></a>01250 <span class="keywordtype">int</span> temp_offset_x, temp_offset_y;
1276
<a name="l01251"></a>01251 <span class="keywordtype">double</span> offset_x, offset_y;
1277
<a name="l01252"></a>01252 <span class="keywordtype">double</span> delta_x, delta_y;
1278
<a name="l01253"></a>01253 <span class="keywordtype">double</span> delta_x2, delta_y2;
1279
<a name="l01254"></a>01254 <span class="keywordtype">double</span> _1_minus_delta_x, _1_minus_delta_y;
1280
<a name="l01255"></a>01255 <span class="keywordtype">double</span> _1_minus_delta_x2, _1_minus_delta_y2;
1281
<a name="l01256"></a>01256 <span class="keywordtype">double</span> _3_minus_2_times_1_minus_delta_x, _3_minus_2_times_1_minus_delta_y;
1282
<a name="l01257"></a>01257 <span class="keywordtype">double</span> _3_minus_2_times_delta_x, _3_minus_2_times_delta_y;
1283
<a name="l01258"></a>01258 <span class="keywordtype">double</span> latitude_dd, longitude_dd;
1284
<a name="l01259"></a>01259 <span class="keywordtype">double</span> height_se, height_sw, height_ne, height_nw;
1285
<a name="l01260"></a>01260 <span class="keywordtype">double</span> w_sw, w_se, w_ne, w_nw;
1286
<a name="l01261"></a>01261 <span class="keywordtype">double</span> south_lat, west_lon;
1287
<a name="l01262"></a>01262 <span class="keywordtype">int</span> end_index = 0;
1288
<a name="l01263"></a>01263 <span class="keywordtype">double</span> skip_factor = 1.0;
1289
<a name="l01264"></a>01264 <span class="keywordtype">char</span> errorStatus[50] = <span class="stringliteral">""</span>;
1290
<a name="l01265"></a>01265
1291
<a name="l01266"></a>01266 <span class="keywordflow">if</span>( ( latitude < -<a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) || ( latitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#aa3c1d8e4403d757bf419703c040db87d">PI_OVER_2</a> ) )
1292
<a name="l01267"></a>01267 { <span class="comment">/* latitude out of range */</span>
1293
<a name="l01268"></a>01268 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a7b9a89219866d2577632c4a2a29f21a7">ErrorMessages::latitude</a> );
1294
<a name="l01269"></a>01269 }
1295
<a name="l01270"></a>01270 <span class="keywordflow">if</span>( ( longitude < -<a class="code" href="_coordinate_conversion_service_8cpp.html#a952eac791b596a61bba0a133a3bb439f">PI</a> ) || ( longitude > <a class="code" href="_albers_equal_area_conic_8cpp.html#a820597b124334bf7fb85214114f4876d">TWO_PI</a> ) )
1296
<a name="l01271"></a>01271 { <span class="comment">/* longitude out of range */</span>
1297
<a name="l01272"></a>01272 strcat( errorStatus, <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_error_messages.html#a148e03a19b3a95cae62e48e757245d12">ErrorMessages::longitude</a> );
1298
<a name="l01273"></a>01273 }
1299
<a name="l01274"></a>01274
1300
<a name="l01275"></a>01275 <span class="keywordflow">if</span>( strlen( errorStatus ) > 0)
1301
<a name="l01276"></a>01276 <span class="keywordflow">throw</span> <a class="code" href="class_m_s_p_1_1_c_c_s_1_1_coordinate_conversion_exception.html">CoordinateConversionException</a>( errorStatus );
1302
<a name="l01277"></a>01277
1303
<a name="l01278"></a>01278 latitude_dd = latitude * _180_OVER_PI;
1304
<a name="l01279"></a>01279 longitude_dd = longitude * _180_OVER_PI;
1305
<a name="l01280"></a>01280
1306
<a name="l01281"></a>01281 <span class="comment">/* Compute X and Y Offsets into Geoid Height Array: */</span>
1307
<a name="l01282"></a>01282
1308
<a name="l01283"></a>01283 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1309
<a name="l01284"></a>01284 {
1310
<a name="l01285"></a>01285 offset_x = ( longitude_dd + 360.0 ) / scale_factor;
1311
<a name="l01286"></a>01286 }
1312
<a name="l01287"></a>01287 <span class="keywordflow">else</span>
1313
<a name="l01288"></a>01288 {
1314
<a name="l01289"></a>01289 offset_x = longitude_dd / scale_factor;
1315
<a name="l01290"></a>01290 }
1316
<a name="l01291"></a>01291 offset_y = ( 90.0 - latitude_dd ) / scale_factor;
1317
<a name="l01292"></a>01292
1318
<a name="l01293"></a>01293 <span class="comment">/* Find Four Nearest Geoid Height Cells for specified latitude, longitude; */</span>
1319
<a name="l01294"></a>01294 <span class="comment">/* Assumes that (0,0) of Geoid Height Array is at Northwest corner: */</span>
1320
<a name="l01295"></a>01295
1321
<a name="l01296"></a>01296 post_x = ( int ) offset_x;
1322
<a name="l01297"></a>01297 <span class="keywordflow">if</span> ( ( post_x + 1 ) == num_cols)
1323
<a name="l01298"></a>01298 post_x--;
1324
<a name="l01299"></a>01299 post_y = ( int )( offset_y + 1.0e-11 );
1325
<a name="l01300"></a>01300 <span class="keywordflow">if</span> ( ( post_y + 1 ) == num_rows)
1326
<a name="l01301"></a>01301 post_y--;
1327
<a name="l01302"></a>01302
1328
<a name="l01303"></a>01303 <span class="keywordflow">if</span>( scale_factor == <a class="code" href="_geoid_library_8cpp.html#a8d54ee4d70e48db25d9b488202563d85">SCALE_FACTOR_30_MINUTES</a> )
1329
<a name="l01304"></a>01304 {
1330
<a name="l01305"></a>01305 skip_factor = 2.0;
1331
<a name="l01306"></a>01306 num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
1332
<a name="l01307"></a>01307 num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
1333
<a name="l01308"></a>01308 }
1334
<a name="l01309"></a>01309 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( scale_factor == <a class="code" href="_geoid_library_8cpp.html#a771e8a447a6e73dc00bf7f3066358470">SCALE_FACTOR_1_DEGREE</a> )
1335
<a name="l01310"></a>01310 {
1336
<a name="l01311"></a>01311 skip_factor = 4.0;
1337
<a name="l01312"></a>01312 num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
1338
<a name="l01313"></a>01313 num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
1339
<a name="l01314"></a>01314 }
1340
<a name="l01315"></a>01315 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( scale_factor == <a class="code" href="_geoid_library_8cpp.html#a56a3075ab7678a7a4a1fc99cf8c97195">SCALE_FACTOR_2_DEGREES</a> )
1341
<a name="l01316"></a>01316 {
1342
<a name="l01317"></a>01317 skip_factor = 8.0;
1343
<a name="l01318"></a>01318 num_rows = <a class="code" href="_geoid_library_8cpp.html#a6b9c7e7d8bf217110a863d62fd264c33">EGM96_ROWS</a>;
1344
<a name="l01319"></a>01319 num_cols = <a class="code" href="_geoid_library_8cpp.html#ab1dcd661211cae9b1ba54b541767b65d">EGM96_COLS</a>;
1345
<a name="l01320"></a>01320 }
1346
<a name="l01321"></a>01321
1347
<a name="l01322"></a>01322 temp_offset_x = ( int )( post_x * skip_factor );
1348
<a name="l01323"></a>01323 temp_offset_y = ( int )( post_y * skip_factor + 1.0e-11 );
1349
<a name="l01324"></a>01324 <span class="keywordflow">if</span> ( ( temp_offset_x + 1 ) == num_cols )
1350
<a name="l01325"></a>01325 temp_offset_x--;
1351
<a name="l01326"></a>01326 <span class="keywordflow">if</span> ( ( temp_offset_y + 1 ) == num_rows )
1352
<a name="l01327"></a>01327 temp_offset_y--;
1353
<a name="l01328"></a>01328
1354
<a name="l01329"></a>01329 <span class="comment">// NW Height</span>
1355
<a name="l01330"></a>01330 index = ( int )( temp_offset_y * num_cols + temp_offset_x );
1356
<a name="l01331"></a>01331 <span class="keywordflow">if</span>( index < 0 )
1357
<a name="l01332"></a>01332 height_nw = height_buffer[ 0 ];
1358
<a name="l01333"></a>01333 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1359
<a name="l01334"></a>01334 height_nw = height_buffer[ max_index ];
1360
<a name="l01335"></a>01335 <span class="keywordflow">else</span>
1361
<a name="l01336"></a>01336 height_nw = height_buffer[ index ];
1362
<a name="l01337"></a>01337 <span class="comment">// NE Height</span>
1363
<a name="l01338"></a>01338 end_index = index + (int)skip_factor;
1364
<a name="l01339"></a>01339 <span class="keywordflow">if</span>( end_index < 0 )
1365
<a name="l01340"></a>01340 height_ne = height_buffer[ 0 ];
1366
<a name="l01341"></a>01341 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( end_index > max_index )
1367
<a name="l01342"></a>01342 height_ne = height_buffer[ max_index ];
1368
<a name="l01343"></a>01343 <span class="keywordflow">else</span>
1369
<a name="l01344"></a>01344 height_ne = height_buffer[ end_index ];
1370
<a name="l01345"></a>01345
1371
<a name="l01346"></a>01346 <span class="comment">// SW Height</span>
1372
<a name="l01347"></a>01347 index = ( int )( ( temp_offset_y + skip_factor ) * num_cols + temp_offset_x );
1373
<a name="l01348"></a>01348 <span class="keywordflow">if</span>( index < 0 )
1374
<a name="l01349"></a>01349 height_sw = height_buffer[ 0 ];
1375
<a name="l01350"></a>01350 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( index > max_index )
1376
<a name="l01351"></a>01351 height_sw = height_buffer[ max_index ];
1377
<a name="l01352"></a>01352 <span class="keywordflow">else</span>
1378
<a name="l01353"></a>01353 height_sw = height_buffer[ index ];
1379
<a name="l01354"></a>01354 <span class="comment">// SE Height</span>
1380
<a name="l01355"></a>01355 end_index = index + (int)skip_factor;
1381
<a name="l01356"></a>01356 <span class="keywordflow">if</span>( end_index < 0 )
1382
<a name="l01357"></a>01357 height_se = height_buffer[ 0 ];
1383
<a name="l01358"></a>01358 <span class="keywordflow">else</span> <span class="keywordflow">if</span>( end_index > max_index )
1384
<a name="l01359"></a>01359 height_se = height_buffer[ max_index ];
1385
<a name="l01360"></a>01360 <span class="keywordflow">else</span>
1386
<a name="l01361"></a>01361 height_se = height_buffer[ end_index ];
1387
<a name="l01362"></a>01362
1388
<a name="l01363"></a>01363 west_lon = post_x * scale_factor;
1389
<a name="l01364"></a>01364
1390
<a name="l01365"></a>01365 <span class="comment">// North latitude - scale_factor</span>
1391
<a name="l01366"></a>01366 south_lat = ( 90 - ( post_y * scale_factor ) ) - scale_factor;
1392
<a name="l01367"></a>01367
1393
<a name="l01368"></a>01368 <span class="comment">/* Perform Non-Linear Interpolation to compute Height above Ellipsoid: */</span>
1394
<a name="l01369"></a>01369
1395
<a name="l01370"></a>01370 <span class="keywordflow">if</span>( longitude_dd < 0.0 )
1396
<a name="l01371"></a>01371 delta_x = ( longitude_dd + 360.0 - west_lon ) / scale_factor;
1397
<a name="l01372"></a>01372 <span class="keywordflow">else</span>
1398
<a name="l01373"></a>01373 delta_x = ( longitude_dd - west_lon ) / scale_factor;
1399
<a name="l01374"></a>01374 delta_y = ( latitude_dd - south_lat ) / scale_factor;
1400
<a name="l01375"></a>01375
1401
<a name="l01376"></a>01376 delta_x2 = delta_x * delta_x;
1402
<a name="l01377"></a>01377 delta_y2 = delta_y * delta_y;
1403
<a name="l01378"></a>01378
1404
<a name="l01379"></a>01379 _1_minus_delta_x = 1 - delta_x;
1405
<a name="l01380"></a>01380 _1_minus_delta_y = 1 - delta_y;
1406
<a name="l01381"></a>01381
1407
<a name="l01382"></a>01382 _1_minus_delta_x2 = _1_minus_delta_x * _1_minus_delta_x;
1408
<a name="l01383"></a>01383 _1_minus_delta_y2 = _1_minus_delta_y * _1_minus_delta_y;
1409
<a name="l01384"></a>01384
1410
<a name="l01385"></a>01385 _3_minus_2_times_1_minus_delta_x = 3 - 2 * _1_minus_delta_x;
1411
<a name="l01386"></a>01386 _3_minus_2_times_1_minus_delta_y = 3 - 2 * _1_minus_delta_y;
1412
<a name="l01387"></a>01387
1413
<a name="l01388"></a>01388 _3_minus_2_times_delta_x = 3 - 2 * delta_x;
1414
<a name="l01389"></a>01389 _3_minus_2_times_delta_y = 3 - 2 * delta_y;
1415
<a name="l01390"></a>01390
1416
<a name="l01391"></a>01391 w_sw = _1_minus_delta_x2 * _1_minus_delta_y2 *
1417
<a name="l01392"></a>01392 ( _3_minus_2_times_1_minus_delta_x * _3_minus_2_times_1_minus_delta_y );
1418
<a name="l01393"></a>01393
1419
<a name="l01394"></a>01394 w_se = delta_x2 * _1_minus_delta_y2 *
1420
<a name="l01395"></a>01395 ( _3_minus_2_times_delta_x * _3_minus_2_times_1_minus_delta_y );
1421
<a name="l01396"></a>01396
1422
<a name="l01397"></a>01397 w_ne = delta_x2 * delta_y2 *
1423
<a name="l01398"></a>01398 ( _3_minus_2_times_delta_x * _3_minus_2_times_delta_y );
1424
<a name="l01399"></a>01399
1425
<a name="l01400"></a>01400 w_nw = _1_minus_delta_x2 * delta_y2 *
1426
<a name="l01401"></a>01401 ( _3_minus_2_times_1_minus_delta_x * _3_minus_2_times_delta_y );
1427
<a name="l01402"></a>01402
1428
<a name="l01403"></a>01403 *delta_height =
1429
<a name="l01404"></a>01404 height_sw * w_sw + height_se * w_se + height_ne * w_ne + height_nw * w_nw;
1430
<a name="l01405"></a>01405 }
1431
<a name="l01406"></a>01406
1432
<a name="l01407"></a>01407 <span class="comment">// CLASSIFICATION: UNCLASSIFIED</span>
1154
1433
</pre></div></div>
1155
<hr size="1"/><address style="text-align: right;"><small>Generated on Wed Nov 25 10:53:09 2009 for MSP GEOTRANS by
1434
<hr size="1"/><address style="text-align: right;"><small>Generated on Tue Aug 3 10:29:15 2010 for MSP GEOTRANS by
1156
1435
<a href="http://www.doxygen.org/index.html">
1157
1436
<img class="footer" src="doxygen.png" alt="doxygen"/></a> 1.6.1 </small></address>
added
removed
CCS/docs/html/_geoid_library_8cpp_source.html
