1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
|
/* reprojection.c
*
* Convert OSM lattitude / longitude from degrees to mercator
* so that Mapnik does not have to project the data again
*
*/
#include <stdio.h>
#include <string.h>
#include <unistd.h>
#include <proj_api.h>
#include "reprojection.h"
static projPJ pj_ll, pj_merc;
static int Proj;
const struct Projection_Info Projection_Info[] = {
[PROJ_LATLONG] = {
descr: "Latlong",
proj4text: "(none)",
srs:4326,
option: "-l" },
[PROJ_MERC] = {
descr: "WGS84 Mercator",
proj4text: "+proj=merc +datum=WGS84 +k=1.0 +units=m +over +no_defs",
srs:3395,
option: "" },
[PROJ_SPHERE_MERC] = {
descr: "Spherical Mercator",
proj4text: "+proj=merc +a=6378137 +b=6378137 +lat_ts=0.0 +lon_0=0.0 +x_0=0.0 +y_0=0 +k=1.0 +units=m +nadgrids=@null +no_defs +over",
srs:900913,
option: "-m" }
};
static struct Projection_Info custom_projection;
// Positive numbers refer the to the table above, negative numbers are
// assumed to refer to EPSG codes and it uses the proj4 to find those.
void project_init(int proj)
{
char buffer[16];
Proj = proj;
if( proj == PROJ_LATLONG )
return;
pj_ll = pj_init_plus("+proj=latlong +ellps=GRS80 +no_defs");
if( proj >= 0 && proj < PROJ_COUNT )
pj_merc = pj_init_plus( Projection_Info[proj].proj4text );
else if( proj < 0 )
{
sprintf( buffer, "+init=epsg:%d", -proj );
pj_merc = pj_init_plus( buffer );
if( !pj_merc )
{
fprintf( stderr, "Couldn't read EPSG definition (do you have /usr/share/proj/epsg?)\n" );
exit(1);
}
}
if (!pj_ll || !pj_merc) {
fprintf(stderr, "Projection code failed to initialise\n");
exit(1);
}
if( proj >= 0 )
return;
custom_projection.srs = -proj;
custom_projection.proj4text = pj_get_def( pj_merc, 0 );
sprintf( buffer, "EPSG:%d", -proj );
custom_projection.descr = strdup(buffer);
custom_projection.option = "-E";
return;
}
void project_exit(void)
{
if( Proj == PROJ_LATLONG )
return;
pj_free(pj_ll);
pj_ll = NULL;
pj_free(pj_merc);
pj_merc = NULL;
}
struct Projection_Info const *project_getprojinfo(void)
{
if( Proj >= 0 )
return &Projection_Info[Proj];
else
return &custom_projection;
}
void reproject(double *lat, double *lon)
{
double x[1], y[1], z[1];
if( Proj == PROJ_LATLONG )
return;
x[0] = *lon * DEG_TO_RAD;
y[0] = *lat * DEG_TO_RAD;
z[0] = 0;
pj_transform(pj_ll, pj_merc, 1, 1, x, y, z);
//printf("%.4f\t%.4f -> %.4f\t%.4f\n", *lat, *lon, y[0], x[0]);
*lat = y[0];
*lon = x[0];
}
|
