代码如下:
#include
#include
#include
#include
#define Num 512
#define PI 3.1415926535898
typedef struct POINT
{
double x, y;
}point;
void read(char * InFILEname, point **a)
{
FILE * fp = fopen(InFILEname, "r");
if (fp == NULL)
printf("Can not open it.");
char s[100];
int i = 1;
fgets(s, 100, fp); //s不是线编号就是文件读取结束的标志“END”
while (s[0] != 'E')
{
int n = 0; //记录该线的点数
int j = 1;
fgets(s, 100, fp); //s不是点对就是线读取结束的标志“END”
while (s[0] != 'E')
{
fseek(fp, -strlen(s), SEEK_CUR);
fscanf(fp, "%lf,%lf", &a[i][j].x, &a[i][j].y);
j++;
fgets(s, 100, fp);
}
a[i][0].x = j - 1; //第i条线有j-1个点
i++;
fgets(s, 100, fp);
}
a[0][0].x = i - 1; //有i-1条线
}
void LambertPro(point ** p, char * OutFILEname)
{
FILE *fp = fopen(OutFILEname, "w");
if (fp == NULL)
printf("Can not open it.");
for (int i = 1; i <= p[0][0].x; i++) //第i条线
{
fprintf(fp, "%d\n", i);
for (int j = 1; j < p[i][0].x; j++) //第j个点
{
double lon = p[i][j].x;
double lat = p[i][j].y; //被投影的点的经纬度
double a = 6378245; //长半轴
double b = 6356863.0188; //短半轴
double e = sqrt(a*a - b*b) / a; //第一偏心率
double originLon = 105 * PI / 180;
double originLat = 0; //原始经纬度
double SP1 = 20 * PI / 180;
double SP2 = 40 * PI / 180; //第一、二标准纬线
double Ef = 609600;
double Nf = 0; //假东和假北
double t = tan(PI / 4 - lat / 2) / pow(((1 - e*sin(lat)) / (1 + e*sin(lat))), e / 2);
double t1 = tan(PI / 4 - SP1 / 2) / pow(((1 - e*sin(SP1)) / (1 + e*sin(SP1))), e / 2);
double t2 = tan(PI / 4 - SP2 / 2) / pow(((1 - e*sin(SP2)) / (1 + e*sin(SP2))), e / 2);
double tF = tan(PI / 4 - originLat / 2) / pow(((1 - e*sin(originLat)) / (1 + e*sin(originLat))), e / 2);
double m1 = cos(SP1) / pow((1 - e*e*sin(SP1)*sin(SP1)), 0.5);
double m2 = cos(SP2) / pow((1 - e*e*sin(SP2)*sin(SP2)), 0.5);
double n = (log(m1) - log(m2)) / (log(t1) - log(t2));
double F = m1 / (n*pow(t1, n));
double A = n*(lon - originLon);
double r = a*F*pow(t, n);
double rF = a*F*pow(tF, n);
double E = Ef + r*sin(A);
double N = Nf + rF - r*cos(A);
fprintf(fp, "%lf,%lf\n", E, N);
}
fprintf(fp, "END\n");
}
fprintf(fp, "END\n");
}
void MercatoPro(point ** p, char * OutFILEname)
{
FILE *fp = fopen(OutFILEname, "w");
if (fp == NULL)
printf("Can not open it.");
for (int i = 1; i <= p[0][0].x; i++) //第i条线
{
fprintf(fp, "%d\n", i);
for (int j = 1; j < p[i][0].x; j++) //第j个点
{
double lon = p[i][j].x;
double lat = p[i][j].y; //被投影的点的经纬度
double a = 6378245; //长半轴
double b = 6356863.0188; //短半轴
double e = sqrt(a*a - b*b) / a; //第一偏心率
double e2 = sqrt(a*a - b*b) / b; //第二偏心率
double L0 = 0; //原始经度
double B0 = 42 * PI / 180; //标准纬度
double K = (a*a / b) / sqrt(1 + e2*e2*cos(B0)*cos(B0)) * cos(B0);
double N = K*log(tan(PI / 4 + lat / 2))*(pow((1 - e*sin(lat)) / (1 + e*sin(lat)), e / 2));
double E = K*(lon - L0);
fprintf(fp, "%lf,%lf\n", E, N);
}
fprintf(fp, "END\n");
}
fprintf(fp, "END\n");
}
void main()
{
point ** a;
a = (point **)malloc(Num*sizeof(point *));
for (int i = 0; i < Num; i++)
{
for (int j = 0; j < 300; j++)
{
a[j] = (point *)malloc(300 * sizeof(point)); //a[0]放数组内元素个数
a[i][j].x = 0;
a[i][j].y = 0;
}
}
read("CHINA_Arc.gen", a);
char outL[15] = "LambertPro";
char outM[15] = "MercatoPro";
LambertPro(a, outL);
MercatoPro(a, outM);
return;
}
读入文件的格式大致如下:(举个例子,随便输入的数字哈)
1
1.1,1.2
2.1,2.1
END
2
3.2,3.2
4.2,4.5
4.3,2.3
END
3
1.3,1.2
END
END