1"""Coordinate reference systems and functions
2
3PROJ.4 is the law of this land: http://proj.osgeo.org/. But whereas PROJ.4
4coordinate reference systems are described by strings of parameters such as
5
6    +proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs
7
8here we use mappings:
9
10    {'proj': 'longlat', 'ellps': 'WGS84', 'datum': 'WGS84', 'no_defs': True}
11"""
12
13from six import string_types
14
15
16def to_string(crs):
17    """Turn a parameter mapping into a more conventional PROJ.4 string.
18
19    Mapping keys are tested against the ``all_proj_keys`` list. Values of
20    ``True`` are omitted, leaving the key bare: {'no_defs': True} -> "+no_defs"
21    and items where the value is otherwise not a str, int, or float are
22    omitted.
23    """
24    items = []
25    for k, v in sorted(filter(
26            lambda x: x[0] in all_proj_keys and x[1] is not False and (
27                isinstance(x[1], (bool, int, float)) or
28                isinstance(x[1], string_types)),
29            crs.items())):
30        items.append(
31            "+" + "=".join(
32                map(str, filter(
33                    lambda y: (y or y == 0) and y is not True, (k, v)))))
34    return " ".join(items)
35
36
37def from_string(prjs):
38    """Turn a PROJ.4 string into a mapping of parameters.
39
40    Bare parameters like "+no_defs" are given a value of ``True``. All keys
41    are checked against the ``all_proj_keys`` list.
42    """
43    parts = [o.lstrip('+') for o in prjs.strip().split()]
44
45    def parse(v):
46        try:
47            return int(v)
48        except ValueError:
49            pass
50        try:
51            return float(v)
52        except ValueError:
53            return v
54    items = map(
55        lambda kv: len(kv) == 2 and (kv[0], parse(kv[1])) or (kv[0], True),
56        (p.split('=') for p in parts))
57    return dict((k, v) for k, v in items if k in all_proj_keys)
58
59
60def from_epsg(code):
61    """Given an integer code, returns an EPSG-like mapping.
62
63    Note: the input code is not validated against an EPSG database.
64    """
65    if int(code) <= 0:
66        raise ValueError("EPSG codes are positive integers")
67    return {'init': "epsg:%s" % code, 'no_defs': True}
68
69
70# Below is the big list of PROJ4 parameters from
71# http://trac.osgeo.org/proj/wiki/GenParms.
72# It is parsed into a list of paramter keys ``all_proj_keys``.
73
74_param_data = """
75+a         Semimajor radius of the ellipsoid axis
76+alpha     ? Used with Oblique Mercator and possibly a few others
77+axis      Axis orientation (new in 4.8.0)
78+b         Semiminor radius of the ellipsoid axis
79+datum     Datum name (see `proj -ld`)
80+ellps     Ellipsoid name (see `proj -le`)
81+init      Initialize from a named CRS
82+k         Scaling factor (old name)
83+k_0       Scaling factor (new name)
84+lat_0     Latitude of origin
85+lat_1     Latitude of first standard parallel
86+lat_2     Latitude of second standard parallel
87+lat_ts    Latitude of true scale
88+lon_0     Central meridian
89+lonc      ? Longitude used with Oblique Mercator and possibly a few others
90+lon_wrap  Center longitude to use for wrapping (see below)
91+nadgrids  Filename of NTv2 grid file to use for datum transforms (see below)
92+no_defs   Don't use the /usr/share/proj/proj_def.dat defaults file
93+over      Allow longitude output outside -180 to 180 range, disables wrapping (see below)
94+pm        Alternate prime meridian (typically a city name, see below)
95+proj      Projection name (see `proj -l`)
96+south     Denotes southern hemisphere UTM zone
97+to_meter  Multiplier to convert map units to 1.0m
98+towgs84   3 or 7 term datum transform parameters (see below)
99+units     meters, US survey feet, etc.
100+vto_meter vertical conversion to meters.
101+vunits    vertical units.
102+x_0       False easting
103+y_0       False northing
104+zone      UTM zone
105+a         Semimajor radius of the ellipsoid axis
106+alpha     ? Used with Oblique Mercator and possibly a few others
107+azi
108+b         Semiminor radius of the ellipsoid axis
109+belgium
110+beta
111+czech
112+e         Eccentricity of the ellipsoid = sqrt(1 - b^2/a^2) = sqrt( f*(2-f) )
113+ellps     Ellipsoid name (see `proj -le`)
114+es        Eccentricity of the ellipsoid squared
115+f         Flattening of the ellipsoid (often presented as an inverse, e.g. 1/298)
116+gamma
117+geoc
118+guam
119+h
120+k         Scaling factor (old name)
121+K
122+k_0       Scaling factor (new name)
123+lat_0     Latitude of origin
124+lat_1     Latitude of first standard parallel
125+lat_2     Latitude of second standard parallel
126+lat_b
127+lat_t
128+lat_ts    Latitude of true scale
129+lon_0     Central meridian
130+lon_1
131+lon_2
132+lonc      ? Longitude used with Oblique Mercator and possibly a few others
133+lsat
134+m
135+M
136+n
137+no_cut
138+no_off
139+no_rot
140+ns
141+o_alpha
142+o_lat_1
143+o_lat_2
144+o_lat_c
145+o_lat_p
146+o_lon_1
147+o_lon_2
148+o_lon_c
149+o_lon_p
150+o_proj
151+over
152+p
153+path
154+proj      Projection name (see `proj -l`)
155+q
156+R
157+R_a
158+R_A       Compute radius such that the area of the sphere is the same as the area of the ellipsoid
159+rf        Reciprocal of the ellipsoid flattening term (e.g. 298)
160+R_g
161+R_h
162+R_lat_a
163+R_lat_g
164+rot
165+R_V
166+s
167+south     Denotes southern hemisphere UTM zone
168+sym
169+t
170+theta
171+tilt
172+to_meter  Multiplier to convert map units to 1.0m
173+units     meters, US survey feet, etc.
174+vopt
175+W
176+westo
177+x_0       False easting
178+y_0       False northing
179+zone      UTM zone
180+wktext    Marker
181"""
182
183_lines = filter(lambda x: len(x) > 1, _param_data.split("\n"))
184all_proj_keys = list(
185    set(line.split()[0].lstrip("+").strip() for line in _lines)) + ['no_mayo']
186