Polar极坐标地图投影,主要用于714数字雷达、多普勒天气雷达图像的显示、处理,以及屏幕坐标(像素点)与经纬度坐标的转换。
Polar.java
001 /**
002 *
003 * Polar 投影(扫描方式,自正北方向顺时针)
004 *
005 * PACKAGE: cma.common.projection
006 * FILENAME: Polar.java
007 * LANGUAGE: Java2 v1.4
008 * ORIGINAL: 无
009 * DESCRIPTION: 极坐标投影(主要用于雷达图像处理)
010 * RELATED: cma.common.projection.Lambert(兰勃特投影)
011 * EDITOR: UltraEdit-32 v12.20a(Windows) NEdit(Linux)
012 * CREATE: 2007-05-06 20:08:23
013 * UPDATE: 2007-06-09
014 * AUTHOR: 刘泽军 (BJ0773@gmail.com)
015 * 广西气象减灾研究所
016 * Guangxi Institude of Meteorology and Disaster-reducing Research(GIMDR)
017 *
018 * Compile : javac Polar.java
019 *
020 * How to use Polar class:
021 *
022 * Polar polar = new Polar(new Point(240, 240), 109.24, 24.35, 1.5);//构造函数
023 * polar.setScale(1.0);//设置比例尺,1公里对应1个像素点
024 * ...
025 *
026 */
027
028 /**
029 *
030 * 扫描平面
031 * /
032 * /
033 * /
034 * /
035 * /
036 * / 仰角
037 * -------------------- 0度平面
038 *
039 * 如图所示:
040 * 扫描平面=>0度平面,需要乘以cos(仰角)
041 * 0度平面=>扫描平面,需要除以cos(仰角)
042 *
043 * 注意,日常显示的雷达图是扫描平面上的图。本类所说的屏幕指扫描平面。
044 *
045 */
046
047 /**
048 * 雷达扫描示意图
049 *
050 * 359 0
051 * | radius
052 * | /
053 * | /
054 * |angle/
055 * | /
056 * | ^ /
057 * | /
058 * | /
059 * |/
060 * 270 -----------------中心----------------- 90
061 * |
062 * |
063 * |
064 * |
065 * |
066 * |
067 * |
068 * |
069 * |
070 * 180
071 */
072
073 /**
074 * Polar.java的函数列表
075 * 1、计算球面距离
076 * double distanceOfSphere(double lon1, double lat1, double lon2, double lat2);
077 * 2、构造函数,不考虑仰角
078 * Polar(Point pos, double lon, double lat);
079 * 3、构造函数,考虑仰角
080 * Polar(Point pos, double lon, double lat, double elv);
081 * 4、获得极坐标中心点的屏幕位置
082 * getCenterPosition();
083 * 5、设置极坐标中心点的屏幕位置
084 * setCenterPosition(Point pos);
085 * 6、获得极坐标中心点的经度
086 * getCenterLongitude();
087 * 7、获得极坐标中心点的纬度
088 * getCenterLatitude();
089 * 8、设置极坐标中心点的经纬度
090 * setCenterCoordinate(double lon, double lat);
091 * 9、获得仰角
092 * getElevation();
093 * 10、设置仰角
094 * setElevation(double elv);
095 * 11、获得缩放比例尺
096 * getScale();
097 * 12、设置缩放比例尺
098 * setScale(double value);
099 * 13、获得经纬度对应的屏幕坐标(扫描平面),主要用于体扫数据显示、底图叠加等。
100 * getPosition(double lon, double lat);
101 * 14、获得极坐标对应的屏幕坐标(扫描平面),主要用于体扫数据显示、底图叠加等。
102 * getXY(double radius, double angle);
103 * 15、根据屏幕坐标获得极半径
104 * getRadius(int x, int y);
105 * 16、根据经纬度坐标获得
106 * getRadius(double lon, double lat);
107 * 17、根据屏幕坐标获得极角
108 * getAngle(int x, int y);
109 * 18、根据经纬度坐标获得极角
110 * getAngle(double lon, double lat);
111 * 19、根据屏幕坐标获得对应的经度值,主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
112 * getLongitude(int x, int y);
113 * 20、根据屏幕坐标获得对应的纬度值,主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
114 * getLatitude(int x, int y);
115 * 21、获得屏幕坐标对应的经纬度
116 * public Point2D.Double getCoordinate(int x, int y);
117 */
118
119 package cma.common.projection;
120
121 import java.awt.*;
122 import java.awt.geom.*;
123 import java.lang.Math.*;
124
125 public class Polar {
126
127 //静态常量,地球半径,来源:《大气科学常用公式》,P601,附录
128 public static double RADIUS = 6371.004 ; //地球平均半径,单位:公里(Km)。
129 public static double RADIUS_POLAR = 6356.755 ; //地球两极半径,单位:公里(Km)。
130 public static double RADIUS_EQUATOR = 6373.140 ; //地球赤道半径,单位:公里(Km)。
131
132 //私有成员
133 private Point centerPosition; //中心对应的屏幕位置
134 private double centerLongitude = 0.0 ; //中心经度
135 private double centerLatitude = 0.0 ; //中心纬度
136 private double perKilometer = 0.0 ; //比例尺:一公里对应的像素点数(扫描平面)
137 private double elevation = 0.0 ; //仰角
138 private double cosineElevation = 1.0 ; //仰角的余弦值
139 private double kmPerDegreeX = 0.0 ; //1经度对应的距离(公里),不同纬度数值不同
140 private double kmPerDegreeY = 0.0 ; //1纬度对应的距离(公里),不同纬度数值不同
141
142 /**
143 * 功能:计算球面上两点间的距离(单位:公里),原在edu.gimdr.Atmos.Meteorology类中写有,为避免import过多的类,故重写一份
144 * 参数:
145 * lon1,lat1 - 第1点的位置(经纬度)
146 * lon2,lat2 - 第2点的位置(经纬度)
147 * 返回值:
148 * 球面距离
149 */
150 public static double distanceOfSphere ( double lon1, double lat1, double lon2, double lat2 ) {
151 /* 公式:
152 A(x,y) B(a,b)
153 AB点的球面距离=R*{arccos[cos(b)*cos(y)*cos(a-x)+sin(b)*sin(y)]}, by Google
154 */
155
156 double rlon1 = Math.toRadians ( lon1 ) ;
157 double rlat1 = Math.toRadians ( lat1 ) ;
158 double rlon2 = Math.toRadians ( lon2 ) ;
159 double rlat2 = Math.toRadians ( lat2 ) ;
160
161 return ( Polar.RADIUS* ( Math.acos ( Math.cos ( rlat2 ) *Math.cos ( rlat1 ) *Math.cos ( rlon2-rlon1 ) +Math.sin ( rlat2 ) *Math.sin ( rlat1 )))) ;
162 }
163
164 /**
165 * 功能:构造函数
166 * 参数:
167 * pos - 中心对应的屏幕位置
168 * lon - 中心对应的经度坐标
169 * lat - 中心对应的纬度坐标
170 * 返回值:
171 * 无
172 */
173 public Polar ( Point pos, double lon, double lat ) {
174 elevation = 0.0 ; //无仰角
175 cosineElevation = 1.0 ;
176 setCenterPosition ( pos ) ;
177 setCenterCoordinate ( lon, lat ) ;
178 }
179
180 /**
181 * 功能:构造函数
182 * 参数:
183 * pos - 中心对应的屏幕位置
184 * lon - 中心对应的经度坐标
185 * lat - 中心对应的纬度坐标
186 * elv = 仰角
187 * 返回值:
188 * 无
189 */
190 public Polar ( Point pos, double lon, double lat, double elv ) {
191 elevation = Math.toRadians ( Math.IEEEremainder ( Math.abs ( elv ) , 90.0 )) ; //在0-90度之间,但不能为90度
192 cosineElevation = Math.cos ( elevation ) ; //仰角的余弦值
193 setCenterPosition ( pos ) ;
194 setCenterCoordinate ( lon, lat ) ;
195 }
196
197 /**
198 * 功能:获得极坐标中心位置
199 * 参数:
200 * 无
201 * 返回值:
202 * 极坐标中心对应的屏幕位置
203 */
204 public Point getCenterPosition () {
205 return ( centerPosition ) ;
206 }
207
208 /**
209 * 功能:设置极坐标的中心位置(屏幕坐标)
210 * 参数:
211 * pos - 新的中心位置(屏幕坐标)
212 * 返回值:
213 * 无
214 */
215 public void setCenterPosition ( Point pos ) {
216 centerPosition = pos;
217 }
218
219 /**
220 * 功能:获得极坐标中心对应的经度坐标
221 * 参数:
222 * 无
223 * 返回值:
224 * 极坐标中心对应的经度坐标
225 */
226 public double getCenterLongitude () {
227 return ( centerLongitude ) ;
228 }
229
230 /**
231 * 功能:获得极坐标中心对应的纬度坐标
232 * 参数:
233 * 无
234 * 返回值:
235 * 极坐标中心对应的纬度坐标
236 */
237 public double getCenterLatitude () {
238 return ( centerLatitude ) ;
239 }
240
241 /**
242 * 功能:设置极坐标的中心位置(经纬度坐标)
243 * 参数:
244 * lon - 新的中心位置(经度值)
245 * lat - 新的中心位置(纬度值)
246 * 返回值:
247 * 无
248 */
249 public void setCenterCoordinate ( double lon, double lat ) {
250 centerLongitude = lon;
251 centerLatitude = lat;
252 //中心经纬度或仰角发生改变,必须重新计算经向和纬向的1度对应的球面距离
253 kmPerDegreeX = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude+ 1.0 , centerLatitude ) / cosineElevation;
254 kmPerDegreeY = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude, centerLatitude+ 1.0 ) / cosineElevation;
255 }
256
257 /**
258 * 功能:获得仰角
259 * 参数:
260 * 无
261 * 返回值:
262 * 仰角的度数
263 */
264 public double getElevation () {
265 return ( Math.toDegrees ( elevation )) ;
266 }
267
268 /**
269 * 功能:设置仰角
270 * 参数:
271 * elv - 新的仰角
272 * 返回值:
273 * 无
274 */
275 public void setElevation ( double elv ) {
276 elevation = Math.toRadians ( Math.IEEEremainder ( Math.abs ( elv ) , 90.0 )) ; //在0-90度之间,但不能为90度
277 cosineElevation = Math.cos ( elevation ) ; //仰角的余弦值
278 //中心经纬度或仰角发生改变,必须重新计算经向和纬向的1度对应的球面距离
279 kmPerDegreeX = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude+ 1.0 , centerLatitude ) / cosineElevation;
280 kmPerDegreeY = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude, centerLatitude+ 1.0 ) / cosineElevation;
281 }
282
283 /**
284 * 功能:获得比例尺,即1公里对应的像素点数
285 * 参数:
286 * 无
287 * 返回值:
288 * 比例尺
289 */
290 public double getScale () {
291 return ( perKilometer ) ;
292 }
293
294 /**
295 * 功能:设置比例尺,即1公里对应的像素点数
296 * 参数:
297 * value - 比例尺数值
298 * 返回值:
299 * 无
300 */
301 public void setScale ( double value ) {
302 perKilometer = value;
303 }
304
305 /**
306 * 功能:获得经纬度对应的屏幕像素坐标,与雷达仰角有关,主要用于体扫数据显示、底图叠加等。
307 * 参数:
308 * lon - 经度
309 * lat - 纬度
310 * 返回值:
311 * 对应的屏幕坐标
312 */
313 public Point getPosition ( double lon, double lat ) {
314 double disX = distanceOfSphere ( lon, centerLatitude, centerLongitude, centerLatitude ) /cosineElevation;
315 double disY = distanceOfSphere ( centerLongitude, lat, centerLongitude, centerLatitude ) /cosineElevation;
316 return ( new Point ( centerPosition.x+ ( lon>centerLongitude? 1 :- 1 ) * ( int )( disX*perKilometer ) ,centerPosition.y- ( lat>centerLatitude? 1 :- 1 ) * ( int )( disY*perKilometer ))) ;
317 }
318
319 /**
320 * 功能:获得极坐标对应的屏幕像素坐标,与雷达仰角无关,主要用于体扫数据显示、底图叠加等。
321 * 参数:
322 * radius - 极半径
323 * angle - 角度(以正北方向顺时针)
324 * 返回值:
325 * 对应的屏幕坐标
326 */
327
328 public Point getXY ( double radius, double angle ) {
329 int x = ( new Double ( radius * Math.sin ( Math.toRadians ( angle )))) .intValue () ;
330 int y = ( new Double ( radius * Math.cos ( Math.toRadians ( angle )))) .intValue () ;
331 return ( new Point ( centerPosition.x+x, centerPosition.y-y )) ;
332 }
333
334 /**
335 * 功能:获得屏幕像素点位置的极坐标半径,由于是输入参数是扫描平面上的值,故与雷达仰角无关。
336 * 参数:
337 * x - 水平坐标
338 * y - 垂直坐标
339 * 返回值:
340 * 与极坐标中心的距离,即极半径
341 */
342 public double getRadius ( int x, int y ) {
343 return ( Math.sqrt ( 1.0 * ( x-centerPosition.x ) * ( x-centerPosition.x ) + 1.0 * ( y-centerPosition.y ) * ( y-centerPosition.y ))) ;
344 }
345
346 /**
347 * 功能:获得经纬度位置的极坐标半径,与雷达仰角有关。
348 * 参数:
349 * lon - 经度坐标
350 * lat - 纬度坐标
351 * 返回值:
352 * 与极坐标中心的距离(象素点),即极半径
353 */
354 public double getRadius ( double lon, double lat ) {
355 Point pos = getPosition ( lon, lat ) ; //此函数已经考虑了仰角的影响
356 return ( getRadius ( pos.x, pos.y )) ;
357 }
358
359 /**
360 * 功能:获得屏幕像素点位置的极坐标角度(扫描平面与0度平面均相同),与雷达仰角无关。
361 * 参数:
362 * x - 水平坐标
363 * y - 垂直坐标
364 * 返回值:
365 * 角度值,自正北方向顺时针
366 */
367 public double getAngle ( int x, int y ) {
368 double agl = 0.0 ;
369 if ( x == centerPosition.x && y == centerPosition.y ) { //重合
370 agl = 0.0 ;
371 }
372 else if ( x == centerPosition.x ) {
373 agl = y > centerPosition.y ? 180.0 : 360.0 ;
374 }
375 else if ( y == centerPosition.y ) {
376 agl = x > centerPosition.x ? 90.0 : 270.0 ;
377 }
378 else {
379 agl = Math.toDegrees ( Math.atan ( 1.0 *Math.abs ( x-centerPosition.x ) /Math.abs ( y-centerPosition.y ))) ;
380 agl =
381 x > centerPosition.x && y < centerPosition.y ? agl : //直角坐标的第一象限
382 x < centerPosition.x && y < centerPosition.y ? 180.0 - agl : //直角坐标的第二象限
383 x < centerPosition.x && y > centerPosition.y ? 180.0 + agl : //直角坐标的第三象限
384 x > centerPosition.x && y > centerPosition.y ? 360.0 - agl : //直角坐标的第四象限
385 agl;
386 }
387 System.out.println ( agl ) ;
388 return ( agl ) ;
389 }
390
391 /**
392 * 功能:获得经纬度位置的极坐标角度(扫描平面与0度平面均相同),与雷达仰角无关。
393 * 参数:
394 * lon - 水平坐标
395 * lat - 垂直坐标
396 * 返回值:
397 * 角度值,自正北方向顺时针
398 */
399 public double getAngle ( double lon, double lat ) {
400 /*
401 //若通过获得屏幕坐标来计算角度,精度比较差,特别是在极坐标中心附近
402 Point p = getPosition(lon, lat);
403 return(getAngle(p.x, p.y);
404 */
405 double agl = 0.0 ;
406 if ( lon == centerLongitude && lat == centerLatitude ) { //重合
407 agl = 0.0 ;
408 }
409 else if ( lon == centerLongitude ) {
410 agl = lat > centerLatitude ? 360.0 : 180.0 ;
411 }
412 else if ( lat == centerLatitude ) {
413 agl = lon > centerLongitude ? 90.0 : 270.0 ;
414 }
415 else {
416 //注:由于经向和纬向的球面距离不等(华南,经向>纬向),故点(1,1)与中心点(0,0)的极角不等45度,而应是略大于45度
417 agl = Math.toDegrees ( Math.atan (( Math.abs ( lon-centerLongitude ) *kmPerDegreeX ) / ( Math.abs ( lat-centerLatitude ) *kmPerDegreeY ))) ;
418 agl =
419 lon > centerLongitude && lat > centerLatitude ? agl : //第一象限
420 lon < centerLongitude && lat > centerLatitude ? 180.0 - agl : //第二象限
421 lon < centerLongitude && lat < centerLatitude ? 180.0 + agl : //第三象限
422 lon > centerLongitude && lat < centerLatitude ? 360.0 - agl : //第四象限
423 agl;
424 }
425 return ( agl ) ;
426 }
427
428 /**
429 * 功能:获得屏幕坐标对应的经度值(根据目标点的经向球面距离来计算,雷达南面和北面的值略有差别),与雷达仰角有关。
430 * 主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
431 * 参数:
432 * x - 横坐标,自西向东
433 * y - 纵坐标,自北向南
434 * 返回值:
435 * 对应的经度坐标
436 */
437 public double getLongitude ( int x, int y ) {
438 double lat = getLatitude ( x, y ) ;
439 double disX0 = distanceOfSphere ( centerLongitude, lat, centerLongitude+ 1.0 , lat ) ; //0度平面上1经度的球面距离
440 double disX = disX0 / cosineElevation; //扫描平面上1经度的距离
441 double perDegreeX = disX * perKilometer; //扫描平面上1经度的对应的像素点数
442 return ( centerLongitude + ( x - centerPosition.x ) / perDegreeX ) ;
443 }
444
445 /**
446 * 功能:获得屏幕坐标对应的纬度值(根据极坐标中心点的纬向球面距离来计算),与雷达仰角有关。
447 * 主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
448 * 参数:
449 * x - 横坐标,自西向东(未用到)
450 * y - 纵坐标,自北向南
451 * 返回值:
452 * 对应的纬度坐标
453 */
454 public double getLatitude ( int x, int y ) {
455 /*
456 目标点 A(X,Y) 弧度
457 中心点 B(A,B) 弧度
458 AB球面距离=R*{arccos[cos(B)*cos(Y)*cos(A-X)+sin(B)*sin(Y)]}, by Google
459 经度相同 => AB = R*{arccos[cos(B)*cos(Y)+sin(B)*sin(Y)]}
460 => AB = R*{arccos[cos(B-Y)]}
461 => AB = R * (B-Y)
462 => AB / R = B - Y
463 => Y = B - AB /R
464 => Y = B - (y-centerPosition.y)/perKilometer/R
465 */
466 return ( Math.toDegrees ( Math.toRadians ( centerLatitude ) + ( centerPosition.y-y ) *cosineElevation/perKilometer/Polar.RADIUS )) ;
467 }
468
469 /**
470 * 功能:
471 * 获得屏幕坐标对应的经纬度
472 * 参数:
473 * x - 屏幕水平坐标
474 * y - 屏幕垂直坐标
475 * 返回值:
476 * 对应的经纬度
477 */
478 public Point2D.Double getCoordinate ( int x, int y ) {
479 double lat = getLatitude ( x, y ) ;
480 double disX0 = distanceOfSphere ( centerLongitude, lat, centerLongitude+ 1.0 , lat ) ; //0度平面上1经度的球面距离
481 double disX = disX0 / cosineElevation; //扫描平面上1经度的距离
482 double perDegreeX = disX * perKilometer; //扫描平面上1经度的对应的像素点数
483 double lon = centerLongitude + ( x - centerPosition.x ) / perDegreeX;
484 return ( new Point2D.Double ( lon, lat )) ;
485 }
486
487 /**
488 * 功能:
489 * 画经线、纬线
490 * 参数:
491 * g - 图形设备
492 * f - 字体
493 * c - 画线颜色
494 * inc_lon - 经线间隔
495 * inc_lat - 纬线间隔
496 * 返回值:
497 * 无
498 */
499 public void drawGridLine ( Graphics2D g, Font f, Color c, int inc_lon, int inc_lat ) {
500 return ;
501 /*
502 Color saveColor = g.getColor();
503 Font saveFont = g.getFont();
504 g.setFont(saveFont);
505 g.setColor(saveColor);
506 */
507 }
508 }
Polar.java
001 /**
002 *
003 * Polar 投影(扫描方式,自正北方向顺时针)
004 *
005 * PACKAGE: cma.common.projection
006 * FILENAME: Polar.java
007 * LANGUAGE: Java2 v1.4
008 * ORIGINAL: 无
009 * DESCRIPTION: 极坐标投影(主要用于雷达图像处理)
010 * RELATED: cma.common.projection.Lambert(兰勃特投影)
011 * EDITOR: UltraEdit-32 v12.20a(Windows) NEdit(Linux)
012 * CREATE: 2007-05-06 20:08:23
013 * UPDATE: 2007-06-09
014 * AUTHOR: 刘泽军 (BJ0773@gmail.com)
015 * 广西气象减灾研究所
016 * Guangxi Institude of Meteorology and Disaster-reducing Research(GIMDR)
017 *
018 * Compile : javac Polar.java
019 *
020 * How to use Polar class:
021 *
022 * Polar polar = new Polar(new Point(240, 240), 109.24, 24.35, 1.5);//构造函数
023 * polar.setScale(1.0);//设置比例尺,1公里对应1个像素点
024 * ...
025 *
026 */
027
028 /**
029 *
030 * 扫描平面
031 * /
032 * /
033 * /
034 * /
035 * /
036 * / 仰角
037 * -------------------- 0度平面
038 *
039 * 如图所示:
040 * 扫描平面=>0度平面,需要乘以cos(仰角)
041 * 0度平面=>扫描平面,需要除以cos(仰角)
042 *
043 * 注意,日常显示的雷达图是扫描平面上的图。本类所说的屏幕指扫描平面。
044 *
045 */
046
047 /**
048 * 雷达扫描示意图
049 *
050 * 359 0
051 * | radius
052 * | /
053 * | /
054 * |angle/
055 * | /
056 * | ^ /
057 * | /
058 * | /
059 * |/
060 * 270 -----------------中心----------------- 90
061 * |
062 * |
063 * |
064 * |
065 * |
066 * |
067 * |
068 * |
069 * |
070 * 180
071 */
072
073 /**
074 * Polar.java的函数列表
075 * 1、计算球面距离
076 * double distanceOfSphere(double lon1, double lat1, double lon2, double lat2);
077 * 2、构造函数,不考虑仰角
078 * Polar(Point pos, double lon, double lat);
079 * 3、构造函数,考虑仰角
080 * Polar(Point pos, double lon, double lat, double elv);
081 * 4、获得极坐标中心点的屏幕位置
082 * getCenterPosition();
083 * 5、设置极坐标中心点的屏幕位置
084 * setCenterPosition(Point pos);
085 * 6、获得极坐标中心点的经度
086 * getCenterLongitude();
087 * 7、获得极坐标中心点的纬度
088 * getCenterLatitude();
089 * 8、设置极坐标中心点的经纬度
090 * setCenterCoordinate(double lon, double lat);
091 * 9、获得仰角
092 * getElevation();
093 * 10、设置仰角
094 * setElevation(double elv);
095 * 11、获得缩放比例尺
096 * getScale();
097 * 12、设置缩放比例尺
098 * setScale(double value);
099 * 13、获得经纬度对应的屏幕坐标(扫描平面),主要用于体扫数据显示、底图叠加等。
100 * getPosition(double lon, double lat);
101 * 14、获得极坐标对应的屏幕坐标(扫描平面),主要用于体扫数据显示、底图叠加等。
102 * getXY(double radius, double angle);
103 * 15、根据屏幕坐标获得极半径
104 * getRadius(int x, int y);
105 * 16、根据经纬度坐标获得
106 * getRadius(double lon, double lat);
107 * 17、根据屏幕坐标获得极角
108 * getAngle(int x, int y);
109 * 18、根据经纬度坐标获得极角
110 * getAngle(double lon, double lat);
111 * 19、根据屏幕坐标获得对应的经度值,主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
112 * getLongitude(int x, int y);
113 * 20、根据屏幕坐标获得对应的纬度值,主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
114 * getLatitude(int x, int y);
115 * 21、获得屏幕坐标对应的经纬度
116 * public Point2D.Double getCoordinate(int x, int y);
117 */
118
119 package cma.common.projection;
120
121 import java.awt.*;
122 import java.awt.geom.*;
123 import java.lang.Math.*;
124
125 public class Polar {
126
127 //静态常量,地球半径,来源:《大气科学常用公式》,P601,附录
128 public static double RADIUS = 6371.004 ; //地球平均半径,单位:公里(Km)。
129 public static double RADIUS_POLAR = 6356.755 ; //地球两极半径,单位:公里(Km)。
130 public static double RADIUS_EQUATOR = 6373.140 ; //地球赤道半径,单位:公里(Km)。
131
132 //私有成员
133 private Point centerPosition; //中心对应的屏幕位置
134 private double centerLongitude = 0.0 ; //中心经度
135 private double centerLatitude = 0.0 ; //中心纬度
136 private double perKilometer = 0.0 ; //比例尺:一公里对应的像素点数(扫描平面)
137 private double elevation = 0.0 ; //仰角
138 private double cosineElevation = 1.0 ; //仰角的余弦值
139 private double kmPerDegreeX = 0.0 ; //1经度对应的距离(公里),不同纬度数值不同
140 private double kmPerDegreeY = 0.0 ; //1纬度对应的距离(公里),不同纬度数值不同
141
142 /**
143 * 功能:计算球面上两点间的距离(单位:公里),原在edu.gimdr.Atmos.Meteorology类中写有,为避免import过多的类,故重写一份
144 * 参数:
145 * lon1,lat1 - 第1点的位置(经纬度)
146 * lon2,lat2 - 第2点的位置(经纬度)
147 * 返回值:
148 * 球面距离
149 */
150 public static double distanceOfSphere ( double lon1, double lat1, double lon2, double lat2 ) {
151 /* 公式:
152 A(x,y) B(a,b)
153 AB点的球面距离=R*{arccos[cos(b)*cos(y)*cos(a-x)+sin(b)*sin(y)]}, by Google
154 */
155
156 double rlon1 = Math.toRadians ( lon1 ) ;
157 double rlat1 = Math.toRadians ( lat1 ) ;
158 double rlon2 = Math.toRadians ( lon2 ) ;
159 double rlat2 = Math.toRadians ( lat2 ) ;
160
161 return ( Polar.RADIUS* ( Math.acos ( Math.cos ( rlat2 ) *Math.cos ( rlat1 ) *Math.cos ( rlon2-rlon1 ) +Math.sin ( rlat2 ) *Math.sin ( rlat1 )))) ;
162 }
163
164 /**
165 * 功能:构造函数
166 * 参数:
167 * pos - 中心对应的屏幕位置
168 * lon - 中心对应的经度坐标
169 * lat - 中心对应的纬度坐标
170 * 返回值:
171 * 无
172 */
173 public Polar ( Point pos, double lon, double lat ) {
174 elevation = 0.0 ; //无仰角
175 cosineElevation = 1.0 ;
176 setCenterPosition ( pos ) ;
177 setCenterCoordinate ( lon, lat ) ;
178 }
179
180 /**
181 * 功能:构造函数
182 * 参数:
183 * pos - 中心对应的屏幕位置
184 * lon - 中心对应的经度坐标
185 * lat - 中心对应的纬度坐标
186 * elv = 仰角
187 * 返回值:
188 * 无
189 */
190 public Polar ( Point pos, double lon, double lat, double elv ) {
191 elevation = Math.toRadians ( Math.IEEEremainder ( Math.abs ( elv ) , 90.0 )) ; //在0-90度之间,但不能为90度
192 cosineElevation = Math.cos ( elevation ) ; //仰角的余弦值
193 setCenterPosition ( pos ) ;
194 setCenterCoordinate ( lon, lat ) ;
195 }
196
197 /**
198 * 功能:获得极坐标中心位置
199 * 参数:
200 * 无
201 * 返回值:
202 * 极坐标中心对应的屏幕位置
203 */
204 public Point getCenterPosition () {
205 return ( centerPosition ) ;
206 }
207
208 /**
209 * 功能:设置极坐标的中心位置(屏幕坐标)
210 * 参数:
211 * pos - 新的中心位置(屏幕坐标)
212 * 返回值:
213 * 无
214 */
215 public void setCenterPosition ( Point pos ) {
216 centerPosition = pos;
217 }
218
219 /**
220 * 功能:获得极坐标中心对应的经度坐标
221 * 参数:
222 * 无
223 * 返回值:
224 * 极坐标中心对应的经度坐标
225 */
226 public double getCenterLongitude () {
227 return ( centerLongitude ) ;
228 }
229
230 /**
231 * 功能:获得极坐标中心对应的纬度坐标
232 * 参数:
233 * 无
234 * 返回值:
235 * 极坐标中心对应的纬度坐标
236 */
237 public double getCenterLatitude () {
238 return ( centerLatitude ) ;
239 }
240
241 /**
242 * 功能:设置极坐标的中心位置(经纬度坐标)
243 * 参数:
244 * lon - 新的中心位置(经度值)
245 * lat - 新的中心位置(纬度值)
246 * 返回值:
247 * 无
248 */
249 public void setCenterCoordinate ( double lon, double lat ) {
250 centerLongitude = lon;
251 centerLatitude = lat;
252 //中心经纬度或仰角发生改变,必须重新计算经向和纬向的1度对应的球面距离
253 kmPerDegreeX = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude+ 1.0 , centerLatitude ) / cosineElevation;
254 kmPerDegreeY = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude, centerLatitude+ 1.0 ) / cosineElevation;
255 }
256
257 /**
258 * 功能:获得仰角
259 * 参数:
260 * 无
261 * 返回值:
262 * 仰角的度数
263 */
264 public double getElevation () {
265 return ( Math.toDegrees ( elevation )) ;
266 }
267
268 /**
269 * 功能:设置仰角
270 * 参数:
271 * elv - 新的仰角
272 * 返回值:
273 * 无
274 */
275 public void setElevation ( double elv ) {
276 elevation = Math.toRadians ( Math.IEEEremainder ( Math.abs ( elv ) , 90.0 )) ; //在0-90度之间,但不能为90度
277 cosineElevation = Math.cos ( elevation ) ; //仰角的余弦值
278 //中心经纬度或仰角发生改变,必须重新计算经向和纬向的1度对应的球面距离
279 kmPerDegreeX = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude+ 1.0 , centerLatitude ) / cosineElevation;
280 kmPerDegreeY = distanceOfSphere ( centerLongitude, centerLatitude, centerLongitude, centerLatitude+ 1.0 ) / cosineElevation;
281 }
282
283 /**
284 * 功能:获得比例尺,即1公里对应的像素点数
285 * 参数:
286 * 无
287 * 返回值:
288 * 比例尺
289 */
290 public double getScale () {
291 return ( perKilometer ) ;
292 }
293
294 /**
295 * 功能:设置比例尺,即1公里对应的像素点数
296 * 参数:
297 * value - 比例尺数值
298 * 返回值:
299 * 无
300 */
301 public void setScale ( double value ) {
302 perKilometer = value;
303 }
304
305 /**
306 * 功能:获得经纬度对应的屏幕像素坐标,与雷达仰角有关,主要用于体扫数据显示、底图叠加等。
307 * 参数:
308 * lon - 经度
309 * lat - 纬度
310 * 返回值:
311 * 对应的屏幕坐标
312 */
313 public Point getPosition ( double lon, double lat ) {
314 double disX = distanceOfSphere ( lon, centerLatitude, centerLongitude, centerLatitude ) /cosineElevation;
315 double disY = distanceOfSphere ( centerLongitude, lat, centerLongitude, centerLatitude ) /cosineElevation;
316 return ( new Point ( centerPosition.x+ ( lon>centerLongitude? 1 :- 1 ) * ( int )( disX*perKilometer ) ,centerPosition.y- ( lat>centerLatitude? 1 :- 1 ) * ( int )( disY*perKilometer ))) ;
317 }
318
319 /**
320 * 功能:获得极坐标对应的屏幕像素坐标,与雷达仰角无关,主要用于体扫数据显示、底图叠加等。
321 * 参数:
322 * radius - 极半径
323 * angle - 角度(以正北方向顺时针)
324 * 返回值:
325 * 对应的屏幕坐标
326 */
327
328 public Point getXY ( double radius, double angle ) {
329 int x = ( new Double ( radius * Math.sin ( Math.toRadians ( angle )))) .intValue () ;
330 int y = ( new Double ( radius * Math.cos ( Math.toRadians ( angle )))) .intValue () ;
331 return ( new Point ( centerPosition.x+x, centerPosition.y-y )) ;
332 }
333
334 /**
335 * 功能:获得屏幕像素点位置的极坐标半径,由于是输入参数是扫描平面上的值,故与雷达仰角无关。
336 * 参数:
337 * x - 水平坐标
338 * y - 垂直坐标
339 * 返回值:
340 * 与极坐标中心的距离,即极半径
341 */
342 public double getRadius ( int x, int y ) {
343 return ( Math.sqrt ( 1.0 * ( x-centerPosition.x ) * ( x-centerPosition.x ) + 1.0 * ( y-centerPosition.y ) * ( y-centerPosition.y ))) ;
344 }
345
346 /**
347 * 功能:获得经纬度位置的极坐标半径,与雷达仰角有关。
348 * 参数:
349 * lon - 经度坐标
350 * lat - 纬度坐标
351 * 返回值:
352 * 与极坐标中心的距离(象素点),即极半径
353 */
354 public double getRadius ( double lon, double lat ) {
355 Point pos = getPosition ( lon, lat ) ; //此函数已经考虑了仰角的影响
356 return ( getRadius ( pos.x, pos.y )) ;
357 }
358
359 /**
360 * 功能:获得屏幕像素点位置的极坐标角度(扫描平面与0度平面均相同),与雷达仰角无关。
361 * 参数:
362 * x - 水平坐标
363 * y - 垂直坐标
364 * 返回值:
365 * 角度值,自正北方向顺时针
366 */
367 public double getAngle ( int x, int y ) {
368 double agl = 0.0 ;
369 if ( x == centerPosition.x && y == centerPosition.y ) { //重合
370 agl = 0.0 ;
371 }
372 else if ( x == centerPosition.x ) {
373 agl = y > centerPosition.y ? 180.0 : 360.0 ;
374 }
375 else if ( y == centerPosition.y ) {
376 agl = x > centerPosition.x ? 90.0 : 270.0 ;
377 }
378 else {
379 agl = Math.toDegrees ( Math.atan ( 1.0 *Math.abs ( x-centerPosition.x ) /Math.abs ( y-centerPosition.y ))) ;
380 agl =
381 x > centerPosition.x && y < centerPosition.y ? agl : //直角坐标的第一象限
382 x < centerPosition.x && y < centerPosition.y ? 180.0 - agl : //直角坐标的第二象限
383 x < centerPosition.x && y > centerPosition.y ? 180.0 + agl : //直角坐标的第三象限
384 x > centerPosition.x && y > centerPosition.y ? 360.0 - agl : //直角坐标的第四象限
385 agl;
386 }
387 System.out.println ( agl ) ;
388 return ( agl ) ;
389 }
390
391 /**
392 * 功能:获得经纬度位置的极坐标角度(扫描平面与0度平面均相同),与雷达仰角无关。
393 * 参数:
394 * lon - 水平坐标
395 * lat - 垂直坐标
396 * 返回值:
397 * 角度值,自正北方向顺时针
398 */
399 public double getAngle ( double lon, double lat ) {
400 /*
401 //若通过获得屏幕坐标来计算角度,精度比较差,特别是在极坐标中心附近
402 Point p = getPosition(lon, lat);
403 return(getAngle(p.x, p.y);
404 */
405 double agl = 0.0 ;
406 if ( lon == centerLongitude && lat == centerLatitude ) { //重合
407 agl = 0.0 ;
408 }
409 else if ( lon == centerLongitude ) {
410 agl = lat > centerLatitude ? 360.0 : 180.0 ;
411 }
412 else if ( lat == centerLatitude ) {
413 agl = lon > centerLongitude ? 90.0 : 270.0 ;
414 }
415 else {
416 //注:由于经向和纬向的球面距离不等(华南,经向>纬向),故点(1,1)与中心点(0,0)的极角不等45度,而应是略大于45度
417 agl = Math.toDegrees ( Math.atan (( Math.abs ( lon-centerLongitude ) *kmPerDegreeX ) / ( Math.abs ( lat-centerLatitude ) *kmPerDegreeY ))) ;
418 agl =
419 lon > centerLongitude && lat > centerLatitude ? agl : //第一象限
420 lon < centerLongitude && lat > centerLatitude ? 180.0 - agl : //第二象限
421 lon < centerLongitude && lat < centerLatitude ? 180.0 + agl : //第三象限
422 lon > centerLongitude && lat < centerLatitude ? 360.0 - agl : //第四象限
423 agl;
424 }
425 return ( agl ) ;
426 }
427
428 /**
429 * 功能:获得屏幕坐标对应的经度值(根据目标点的经向球面距离来计算,雷达南面和北面的值略有差别),与雷达仰角有关。
430 * 主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
431 * 参数:
432 * x - 横坐标,自西向东
433 * y - 纵坐标,自北向南
434 * 返回值:
435 * 对应的经度坐标
436 */
437 public double getLongitude ( int x, int y ) {
438 double lat = getLatitude ( x, y ) ;
439 double disX0 = distanceOfSphere ( centerLongitude, lat, centerLongitude+ 1.0 , lat ) ; //0度平面上1经度的球面距离
440 double disX = disX0 / cosineElevation; //扫描平面上1经度的距离
441 double perDegreeX = disX * perKilometer; //扫描平面上1经度的对应的像素点数
442 return ( centerLongitude + ( x - centerPosition.x ) / perDegreeX ) ;
443 }
444
445 /**
446 * 功能:获得屏幕坐标对应的纬度值(根据极坐标中心点的纬向球面距离来计算),与雷达仰角有关。
447 * 主要用于雷达产品的定位、底图叠加、转换为经纬度网格产品、拼图等。
448 * 参数:
449 * x - 横坐标,自西向东(未用到)
450 * y - 纵坐标,自北向南
451 * 返回值:
452 * 对应的纬度坐标
453 */
454 public double getLatitude ( int x, int y ) {
455 /*
456 目标点 A(X,Y) 弧度
457 中心点 B(A,B) 弧度
458 AB球面距离=R*{arccos[cos(B)*cos(Y)*cos(A-X)+sin(B)*sin(Y)]}, by Google
459 经度相同 => AB = R*{arccos[cos(B)*cos(Y)+sin(B)*sin(Y)]}
460 => AB = R*{arccos[cos(B-Y)]}
461 => AB = R * (B-Y)
462 => AB / R = B - Y
463 => Y = B - AB /R
464 => Y = B - (y-centerPosition.y)/perKilometer/R
465 */
466 return ( Math.toDegrees ( Math.toRadians ( centerLatitude ) + ( centerPosition.y-y ) *cosineElevation/perKilometer/Polar.RADIUS )) ;
467 }
468
469 /**
470 * 功能:
471 * 获得屏幕坐标对应的经纬度
472 * 参数:
473 * x - 屏幕水平坐标
474 * y - 屏幕垂直坐标
475 * 返回值:
476 * 对应的经纬度
477 */
478 public Point2D.Double getCoordinate ( int x, int y ) {
479 double lat = getLatitude ( x, y ) ;
480 double disX0 = distanceOfSphere ( centerLongitude, lat, centerLongitude+ 1.0 , lat ) ; //0度平面上1经度的球面距离
481 double disX = disX0 / cosineElevation; //扫描平面上1经度的距离
482 double perDegreeX = disX * perKilometer; //扫描平面上1经度的对应的像素点数
483 double lon = centerLongitude + ( x - centerPosition.x ) / perDegreeX;
484 return ( new Point2D.Double ( lon, lat )) ;
485 }
486
487 /**
488 * 功能:
489 * 画经线、纬线
490 * 参数:
491 * g - 图形设备
492 * f - 字体
493 * c - 画线颜色
494 * inc_lon - 经线间隔
495 * inc_lat - 纬线间隔
496 * 返回值:
497 * 无
498 */
499 public void drawGridLine ( Graphics2D g, Font f, Color c, int inc_lon, int inc_lat ) {
500 return ;
501 /*
502 Color saveColor = g.getColor();
503 Font saveFont = g.getFont();
504 g.setFont(saveFont);
505 g.setColor(saveColor);
506 */
507 }
508 }