-
-
Save fronteer-kr/14d7f779d52a21ac2f16 to your computer and use it in GitHub Desktop.
| // 소스출처 : http://www.kma.go.kr/weather/forecast/digital_forecast.jsp 내부에 있음 | |
| // 기상청에서 이걸 왜 공식적으로 공개하지 않을까? | |
| // | |
| // (사용 예) | |
| // var rs = dfs_xy_conv("toLL","60","127"); | |
| // console.log(rs.lat, rs.lng); | |
| // | |
| <script language="javascript"> | |
| //<!-- | |
| // | |
| // LCC DFS 좌표변환을 위한 기초 자료 | |
| // | |
| var RE = 6371.00877; // 지구 반경(km) | |
| var GRID = 5.0; // 격자 간격(km) | |
| var SLAT1 = 30.0; // 투영 위도1(degree) | |
| var SLAT2 = 60.0; // 투영 위도2(degree) | |
| var OLON = 126.0; // 기준점 경도(degree) | |
| var OLAT = 38.0; // 기준점 위도(degree) | |
| var XO = 43; // 기준점 X좌표(GRID) | |
| var YO = 136; // 기1준점 Y좌표(GRID) | |
| // | |
| // LCC DFS 좌표변환 ( code : "toXY"(위경도->좌표, v1:위도, v2:경도), "toLL"(좌표->위경도,v1:x, v2:y) ) | |
| // | |
| function dfs_xy_conv(code, v1, v2) { | |
| var DEGRAD = Math.PI / 180.0; | |
| var RADDEG = 180.0 / Math.PI; | |
| var re = RE / GRID; | |
| var slat1 = SLAT1 * DEGRAD; | |
| var slat2 = SLAT2 * DEGRAD; | |
| var olon = OLON * DEGRAD; | |
| var olat = OLAT * DEGRAD; | |
| var sn = Math.tan(Math.PI * 0.25 + slat2 * 0.5) / Math.tan(Math.PI * 0.25 + slat1 * 0.5); | |
| sn = Math.log(Math.cos(slat1) / Math.cos(slat2)) / Math.log(sn); | |
| var sf = Math.tan(Math.PI * 0.25 + slat1 * 0.5); | |
| sf = Math.pow(sf, sn) * Math.cos(slat1) / sn; | |
| var ro = Math.tan(Math.PI * 0.25 + olat * 0.5); | |
| ro = re * sf / Math.pow(ro, sn); | |
| var rs = {}; | |
| if (code == "toXY") { | |
| rs['lat'] = v1; | |
| rs['lng'] = v2; | |
| var ra = Math.tan(Math.PI * 0.25 + (v1) * DEGRAD * 0.5); | |
| ra = re * sf / Math.pow(ra, sn); | |
| var theta = v2 * DEGRAD - olon; | |
| if (theta > Math.PI) theta -= 2.0 * Math.PI; | |
| if (theta < -Math.PI) theta += 2.0 * Math.PI; | |
| theta *= sn; | |
| rs['x'] = Math.floor(ra * Math.sin(theta) + XO + 0.5); | |
| rs['y'] = Math.floor(ro - ra * Math.cos(theta) + YO + 0.5); | |
| } | |
| else { | |
| rs['x'] = v1; | |
| rs['y'] = v2; | |
| var xn = v1 - XO; | |
| var yn = ro - v2 + YO; | |
| ra = Math.sqrt(xn * xn + yn * yn); | |
| if (sn < 0.0) - ra; | |
| var alat = Math.pow((re * sf / ra), (1.0 / sn)); | |
| alat = 2.0 * Math.atan(alat) - Math.PI * 0.5; | |
| if (Math.abs(xn) <= 0.0) { | |
| theta = 0.0; | |
| } | |
| else { | |
| if (Math.abs(yn) <= 0.0) { | |
| theta = Math.PI * 0.5; | |
| if (xn < 0.0) - theta; | |
| } | |
| else theta = Math.atan2(xn, yn); | |
| } | |
| var alon = theta / sn + olon; | |
| rs['lat'] = alat * RADDEG; | |
| rs['lng'] = alon * RADDEG; | |
| } | |
| return rs; | |
| } | |
| //--> | |
| </script> |
기상청에서 제공한 위 경도 -> 자체 격자 변환 소스코드입니다. 혹시나 해서 같이 올립니다.
/********************************************************************************/
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <signal.h>
#include <sys/types.h>
#include <sys/stat.h>
#include <dirent.h>
#include <time.h>
#include <math.h>
#define NX 149 /* X축 격자점 수 /
#define NY 253 / Y축 격자점 수 */
struct lamc_parameter {
float Re; /* 사용할 지구반경 [ km ] /
float grid; / 격자간격 [ km ] /
float slat1; / 표준위도 [degree] /
float slat2; / 표준위도 [degree] /
float olon; / 기준점의 경도 [degree] /
float olat; / 기준점의 위도 [degree] /
float xo; / 기준점의 X좌표 [격자거리] /
float yo; / 기준점의 Y좌표 [격자거리] /
int first; / 시작여부 (0 = 시작) */
};
/******************************************************************************
*
- MAIN
******************************************************************************/
int main (int argc, char *argv[]) {
float lon, lat, x, y;
struct lamc_parameter map;
//
// 인수 확인
//
if (argc != 4) {
printf("[Usage] %s 1 <X-grid><Y-grid>\n", argv[0]);
printf(" %s 0 <longitude><latitude>\n", argv[0]);
exit(0);
}
if (atoi(argv[1]) == 1) {
x = atof(argv[2]);
y = atof(argv[3]);
if (x < 1 || x > NX || y < 1 || y > NY) {
printf("X-grid range [1,%d] / Y-grid range [1,%d]\n", NX, NY);
exit(0);
}
} else if (atoi(argv[1]) == 0) {
lon = atof(argv[2]);
lat = atof(argv[3]);
}
//
// 단기예보 지도 정보
//
map.Re = 6371.00877; // 지도반경
map.grid = 5.0; // 격자간격 (km)
map.slat1 = 30.0; // 표준위도 1
map.slat2 = 60.0; // 표준위도 2
map.olon = 126.0; // 기준점 경도
map.olat = 38.0; // 기준점 위도
map.xo = 210/map.grid; // 기준점 X좌표
map.yo = 675/map.grid; // 기준점 Y좌표
map.first = 0;
//
// 단기예보
//
map_conv(&lon, &lat, &x, &y, atoi(argv[1]), map);
if (atoi(argv[1]))
printf("X = %d, Y = %d --->lon.= %f, lat.= %f\n", (int)x, (int)y, lon, lat);
else
printf("lon.= %f, lat.= %f ---> X = %d, Y = %d\n", lon, lat, (int)x, (int)y);
return 0;
}
/============================================================================
-
좌표변환
============================================================================/
int map_conv
(
float *lon, // 경도(degree)
float *lat, // 위도(degree)
float *x, // X격자 (grid)
float *y, // Y격자 (grid)
int code, // 0 (격자->위경도), 1 (위경도->격자)
struct lamc_parameter map // 지도정보
) {
float lon1, lat1, x1, y1;//
// 위경도 -> (X,Y)
//if (code == 0) {
lon1 = *lon;
lat1 = *lat;
lamcproj(&lon1, &lat1, &x1, &y1, 0, &map);
*x = (int)(x1 + 1.5);
*y = (int)(y1 + 1.5);
}//
// (X,Y) -> 위경도
//if (code == 1) {
x1 = *x - 1;
y1 = *y - 1;
lamcproj(&lon1, &lat1, &x1, &y1, 1, &map);
*lon = lon1;
*lat = lat1;
}
return 0;
}
/***************************************************************************
*
- [ Lambert Conformal Conic Projection ]
- olon, lat : (longitude,latitude) at earth [degree]
- o x, y : (x,y) cordinate in map [grid]
- o code = 0 : (lon,lat) --> (x,y)
- 1 : (x,y) --> (lon,lat)
***************************************************************************/
int lamcproj(lon, lat, x, y, code, map)
float *lon, lat; / Longitude, Latitude [degree] */
float *x, y; / Coordinate in Map [grid] /
int code; / (0) lon,lat ->x,y (1) x,y ->lon,lat */
struct lamc_parameter *map;
{
static double PI, DEGRAD, RADDEG;
static double re, olon, olat, sn, sf, ro;
double slat1, slat2, alon, alat, xn, yn, ra, theta;
if ((*map).first == 0) {
PI = asin(1.0)*2.0;
DEGRAD = PI/180.0;
RADDEG = 180.0/PI;
re = (*map).Re/(*map).grid;
slat1 = (*map).slat1 * DEGRAD;
slat2 = (*map).slat2 * DEGRAD;
olon = (*map).olon * DEGRAD;
olat = (*map).olat * DEGRAD;
sn = tan(PI*0.25 + slat2*0.5)/tan(PI*0.25 + slat1*0.5);
sn = log(cos(slat1)/cos(slat2))/log(sn);
sf = tan(PI*0.25 + slat1*0.5);
sf = pow(sf,sn)*cos(slat1)/sn;
ro = tan(PI*0.25 + olat*0.5);
ro = re*sf/pow(ro,sn);
(*map).first = 1;
}
if (code == 0) {
ra = tan(PI*0.25+(*lat)*DEGRAD*0.5);
ra = re*sf/pow(ra,sn);
theta = (*lon)*DEGRAD - olon;
if (theta > PI) theta -= 2.0*PI;
if (theta < -PI) theta += 2.0*PI;
theta *= sn;
*x = (float)(ra*sin(theta)) + (*map).xo;
*y = (float)(ro - ra*cos(theta)) + (*map).yo;
} else {
xn = *x - (*map).xo;
yn = ro - *y + (*map).yo;
ra = sqrt(xn*xn+yn*yn);
if (sn< 0.0) -ra;
alat = pow((re*sf/ra),(1.0/sn));
alat = 2.0*atan(alat) - PI*0.5;
if (fabs(xn) <= 0.0) {
theta = 0.0;
} else {
if (fabs(yn) <= 0.0) {
theta = PI*0.5;
if(xn< 0.0 ) -theta;
} else
theta = atan2(xn,yn);
}
alon = theta/sn + olon;
*lat = (float)(alat*RADDEG);
*lon = (float)(alon*RADDEG);
}
return 0;
}
@frogio "찾아보면 있어요"라고 해서 찾아봤습니다. 찾기 어렵더만요.
구글에서 아래 쿼리로 하면 나오기는 하네요.
"위경도 값을 직접 좌표 값으로 변환하여 사용하기 원하는 사용자를 위한 예제"