本文介绍了一种称为极射赤面投影的空间映射方法，该方法主要用于气象领域中的地图绘制。文章详细解释了如何通过Java实现这一投影，并提供了具体的算法步骤，包括参数设置、坐标转换等关键环节。
Stereogram.java
001 /*
002
003 极射赤面投影(Stereogram projection)
004
005 PACKAGE: cma.common.projection
006 FILENAME: Coordinate.java
007 LANGUAGE: Java2 v1.4
008 ORIGINAL: none
009 DESCRIPTION:
010 CREATE: 2007-07-19 09:29:06
011 UPDATE: 2007-07-23
012 AUTHOR: 刘泽军 (BJ0773@gmail.com)
013 广西气象减灾研究所
014 Guangxi Institude of Meteorology and Disaster-reducing Research(GIMDR)
015 REFERENCE: http://mathworld.wolfram.com/StereographicProjection.html
016
017 Compile: javac Coordinate.java Stereogram.java
018
019 */
020
021 package cma.common.projection;
022
023 import java.io.*;
024 import java.awt.*;
025 import java.awt.geom.*;
026 import java.util.*;
027 import java.lang.Math.*;
028 import java.text.DecimalFormat;
029
030 import cma.common.algorithm.Geometry;//矢量的加法、减法、叉乘
031
032 public class Stereogram extends Coordinate {
033
034 /**
035 * 功能:
036 * 重置参数
037 * 参数:
038 * lon0,lat0 - 标准经纬度，lat0只能为90或-90，即BBQ或NBQ
039 * lon1,lat1 - 中心经纬度(用于定位，不一定是中心)
040 * x,y - 中心经纬度对应的屏幕坐标
041 * sc - 缩放系数
042 * 返回值:
043 * 无
044 */
045 private void reset(//Degree
046 double lon0, double lat0,
047 double lon1, double lat1,
048 int x, int y,
049 double sc) {
050 type = lat0 >= 0.0 ? Coordinate.BBQ : Coordinate.NBQ;
051 standard = new Point2D.Double(Math.IEEEremainder(Math.abs(lon0), 360.0), lat0>=0.0?90.0:-90.0);
052 center = new Point2D.Double(
053 Math.IEEEremainder(Math.abs(lon1), 360.0),
054 type == Coordinate.BBQ ?
055 (lat1 < -85.0 ? -85.0 : lat1 > 90.0 ? 90.0 : lat1) : //北半球可以以北极为定位点
056 (lat1 > 85.0 ? 85.0 : lat1 < -90.0 ? -90.0 : lat1) //南半球可以以南极为定位点
057 );
058 place = new Point(x, y);
059 scale = sc <= 0.0 ? 1.0 : sc;
060 double _lon0 = Math.toRadians(standard.x);
061 double _lat0 = Math.toRadians(standard.y);
062 double _lon1 = Math.toRadians(lon1);
063 double _lat1 = Math.toRadians(lat1);
064 //0.04149是与Micaps1.0匹配的参数
065 double k = scale * 0.04149 * 2.0*Coordinate.RADIUS/(1.0+Math.sin(_lat0)*Math.sin(_lat1)+Math.cos(_lat0)*Math.cos(_lat1)*Math.cos(_lon1-_lon0));
066 double x0 = k * Math.cos(_lat1)*Math.sin(_lon1-_lon0);
067 double y0 = k * (Math.cos(_lat0)*Math.sin(_lat1)-Math.sin(_lat0)*Math.cos(_lat1)*Math.cos(_lon1-_lon0));
068 offset = new Point((int)(0.5+place.x-x0), (int)(0.5+place.y-y0));
069 }
070
071 /**
072 * 功能:
073 * 构造函数
074 * 参数:
075 * 无(使用缺省值)
076 * 返回值:
077 * 无
078 */
079 public Stereogram() {
080 reset(110.0, 90.0, 109.40, 24.35, 640, 480, 1.0);
081 }
082
083 /**
084 * 功能:
085 * 构造函数
086 * 参数:
087 * lon,lat - 中心经纬度，
088 * px,py - 中心经纬度对应的屏幕坐标
089 * sc - 缩放系数
090 * 返回值:
091 * 无
092 */
093 public Stereogram(double lon, double lat, int px, int py, double sc) {
094 reset(110.0, 35.0, lon, lat, px, py, sc);
095 }
096
097 /**
098 * 功能:
099 * 构造函数
100 * 参数:
101 * lon,lat - 中心经纬度，
102 * px,py - 中心经纬度对应的屏幕坐标
103 * sx,sy - 缩放比例
104 * 返回值:
105 * 无
106 */
107 public Stereogram(double lon0, double lat0, double lon1, double lat1, int px, int py, double sc) {
108 reset(lon0, lat0, lon1, lat1, px, py, sc);
109 }
110
111 /**
112 * 功能:
113 * 获得屏幕坐标
114 * 参数:
115 * lon - 经度
116 * lat - 纬度
117 * 返回值:
118 * 屏幕坐标
119 */
120 public Point getPosition(double lon, double lat) {
121 /*
122 double _lon0 = Math.toRadians(standard.x);
123 double _lat0 = Math.toRadians(standard.y);
124 double _lon1 = Math.toRadians(lon);
125 double _lat1 = Math.toRadians(lat);
126 //0.04149是与Micaps1.0匹配的参数，
127 double k = scale * 0.04149 * 2.0*Coordinate.RADIUS/(1.0+Math.sin(_lat0)*Math.sin(_lat1)+Math.cos(_lat0)*Math.cos(_lat1)*Math.cos(_lon1-_lon0));
128 double x = k * Math.cos(_lat1)*Math.sin(_lon1-_lon0);
129 double y = k * (Math.cos(_lat0)*Math.sin(_lat1)-Math.sin(_lat0)*Math.cos(_lat1)*Math.cos(_lon1-_lon0));
130 return(new Point((int)(0.5+offset.x+x), (int)(0.5+offset.y-y)));
131 */
132 double deltaLon = Math.toRadians(lon-standard.x);
133 double radianLat = Math.toRadians(lat);
134 //sign=sin(+/-90)=+/-1, cos(90)=cos(-90)=0.0，优化，去掉与标准纬度(即两个极点)有关的三角函数
135 double sgn = standard.y == 90.0 ? 1.0 : -1.0;
136 //0.04149是与Micaps1.0匹配的参数，
137 double k = scale * 0.04149 * 2.0*Coordinate.RADIUS/(1.0+sgn*Math.sin(radianLat));
138 double x = k * Math.cos(radianLat)*Math.sin(deltaLon);
139 double y = k * (-sgn*Math.cos(radianLat)*Math.cos(deltaLon));
140 return(new Point((int)(0.5+offset.x+x), (int)(0.5+offset.y-y)));
141 }
142
143 /**
144 * 功能:
145 * 获得屏幕坐标对应的经纬度
146 * 参数:
147 * x - 屏幕水平坐标
148 * y - 屏幕垂直坐标
149 * 返回值:
150 * 对应的经纬度
151 */
152 public Point2D.Double getCoordinate(int x, int y) {
153 Point p0 = this.getPosition(0.0, standard.y);
154 if( p0.x == x && p0.y == y ) {
155 return(new Point2D.Double(0.0, standard.y));
156 }
157 Point p1 = this.getPosition(0.0, 0.0);
158 if( p1.x == x && p1.y == y ) {
159 return(new Point2D.Double(0.0, 0.0));
160 }
161 //求未知点与极点、(0.0E,0.0N)的夹角
162 double x1 = p0.x;
163 double y1 = p0.y;
164 double x2 = p1.x;
165 double y2 = p1.y;
166 double x3 = x;
167 double y3 = y;
168 //增加临时变量，但节约了计算量，由加减法17次，乘除法12次，调用函数5次，优化为加减法11次，乘除法8次，调用函数4次
169 double x3_x2 = x3-x2;
170 double y3_y2 = y3-y2;
171 double aa = x3_x2*x3_x2+y3_y2*y3_y2;
172 double x3_x1 = x3-x1;
173 double y3_y1 = y3-y1;
174 double bb = x3_x1*x3_x1+y3_y1*y3_y1;
175 double b = Math.sqrt(bb);
176 double x1_x2 = x1-x2;
177 double y1_y2 = y1-y2;
178 double cc = x1_x2*x1_x2+y1_y2*y1_y2;
179 double c = Math.sqrt(cc);
180 double angle = Math.acos((bb+cc-aa)/(2*b*c));//余弦定理，值在0-180之间
181 double pos = Geometry.pos(p1.x, p1.y, p0.x, p0.y, x, y);//判断点相对于0度经线的位置
182 if( pos > 0.0 ) {//位于180-360度区域
183 angle = 2.0 * Math.PI - angle;
184 }
185 double lon = Math.toDegrees(angle);
186
187 double deltaLon = Math.toRadians(lon-standard.x);
188 double sgn = standard.y == 90.0 ? 1.0 : -1.0;
189 double k = (offset.y - y) / (-sgn * scale * 0.04149 * 2.0 * Coordinate.RADIUS * Math.cos(deltaLon));
190 double lat = Math.toDegrees(Math.asin((1-sgn*k*k)/(1+k*k)));//由getPosition方法的公式反推得到
191
192 return(new Point2D.Double(lon, lat));
193 }
194
195 /**
196 * 功能：获得相对于标准经度的偏角(逆时针方向)
197 * 参数：
198 * degreeLon - 经度值
199 * degreeLat - 纬度值
200 * 返回值：
201 * 偏角
202 */
203 public double getAngle(double degreeLon, double degreeLat) {
204 double agl = 360.0+standard.x-degreeLon;
205 while(agl>=360.0) {
206 agl -= 360.0;
207 }
208 while(agl<=-360.0) {
209 agl += 360.0;
210 }
211 return(agl);
212 }
213
214 /**
215 * 功能：
216 * 画经线、纬线
217 * 参数：
218 * g - 图形设备
219 * f - 字体
220 * c - 画线颜色
221 * inc_lon - 经线间隔
222 * inc_lat - 纬线间隔
223 * 返回值：
224 * 无
225 */
226 public void drawGridLine(Graphics2D g, Font f, Color c, int inc_lon, int inc_lat) {
227
228 DecimalFormat df = new DecimalFormat("0.##########");
229 String lonString = df.format(inc_lon);
230 String latString = df.format(inc_lat);
231 DecimalFormat dfLon = new DecimalFormat("0."+("##########".substring(0, lonString.indexOf(".")>0?lonString.length()-lonString.indexOf("."):1)));
232 DecimalFormat dfLat = new DecimalFormat("0."+("##########".substring(0, latString.indexOf(".")>0?latString.length()-latString.indexOf("."):1)));
233
234 Color saveColor = g.getColor();
235 Font saveFont = g.getFont();
236 g.setColor(c);
237 g.setFont(null==f?f:new Font("Times New Roman", Font.PLAIN, 12));
238 FontMetrics fm = g.getFontMetrics();
239 String text;
240 byte tmpByte[];
241 int bytesWidth, bytesHeight = fm.getHeight();
242 Point pos0, pos1, pos2;
243 pos0 = getPosition(standard.x, type==Coordinate.BBQ?90.0:-90.0);
244 double rotateAngle;
245 if( inc_lon > 0 ) {//经线间隔为正值
246 for(double lon=0.0;lon<360.0;lon=lon+inc_lon) {//画经线
247 pos1 = getPosition(lon,0.0);//在赤道标经度值
248 g.translate(pos1.x, pos1.y);
249 rotateAngle = this.getAngle(lon, 0.0);
250 g.rotate(Math.toRadians(rotateAngle)); //偏角
251 text = dfLon.format(lon);
252 tmpByte = text.getBytes();
253 bytesWidth = fm.bytesWidth(tmpByte, 0, tmpByte.length);
254 g.drawString(text, -bytesWidth/2, bytesHeight/3);
255 g.rotate(Math.toRadians(-rotateAngle)); //偏角
256 g.translate(-pos1.x, -pos1.y);
257
258 pos1 = getPosition(lon,-90.0+(inc_lat>0.0?inc_lat:10.0));
259 pos2 = getPosition(lon, 90.0-(inc_lat>0.0?inc_lat:10.0));
260 g.drawLine(pos1.x, pos1.y, pos2.x, pos2.y);//画经线
261 }
262 }
263 if( inc_lat > 0 ) {//纬线间隔为正值
264 int r = 0;
265 for(double lat=90.0-inc_lat;lat>=-90.0+inc_lat;lat=lat-inc_lat) {//画纬线(纬圈)
266 pos1 = getPosition(standard.x, lat);
267 r = Math.abs(pos1.y-pos0.y);
268 g.drawArc(pos0.x-r, pos0.y-r,2*r,2*r,0,360);
269 }
270 for(double lon=0.0;lon<360.0;lon=lon+90.0) {//在0、90、180、270经线上标纬度值
271 for(double lat=90.0-inc_lat;lat>=-90.0+inc_lat;lat=lat-inc_lat) {
272 pos1 = getPosition(lon,lat);
273 g.translate(pos1.x, pos1.y);
274 rotateAngle = this.getAngle(lon, lat);
275 g.rotate(Math.toRadians(rotateAngle)); //偏角
276 text = dfLat.format(lat);
277 tmpByte = text.getBytes();
278 bytesWidth = fm.bytesWidth(tmpByte, 0, tmpByte.length);
279 g.drawString(text, -bytesWidth/2, bytesHeight/3);
280 // g.drawString(text, pos1.x-bytesWidth/2, pos1.y+bytesHeight/3);
281 g.rotate(Math.toRadians(-rotateAngle)); //偏角
282 g.translate(-pos1.x, -pos1.y);
283 }
284 }
285 }
286 g.setFont(saveFont);
287 g.setColor(saveColor);
288 }
289 }
扫一扫
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
