1#!/usr/local/bin/python3.8 2 3""" 4MODULE: m.proj 5 6AUTHOR(S): M. Hamish Bowman, Dept. Marine Science, Otago University, 7 New Zealand 8 Converted to Python by Glynn Clements 9 10PURPOSE: cs2cs reprojection frontend for a list of coordinates. 11 Replacement for m.proj2 from GRASS 5 12 13COPYRIGHT: (c) 2006-2019 Hamish Bowman, and the GRASS Development Team 14 This program is free software under the GNU General Public 15 License (>=v2). Read the file COPYING that comes with GRASS 16 for details. 17""" 18 19# notes: 20# - cs2cs expects "x y" data so be sure to send it "lon lat" not "lat lon" 21# - if you send cs2cs a third data column, beware it might be treated as "z" 22# todo: 23# - `cut` away x,y columns into a temp file, feed to cs2cs, then `paste` 24# back to input file. see method in v.in.garmin.sh. that way additional 25# numeric and string columns would survive the trip, and 3rd column would 26# not be modified as z. 27 28#%module 29#% description: Converts coordinates from one projection to another (cs2cs frontend). 30#% keyword: miscellaneous 31#% keyword: projection 32#% keyword: transformation 33#%end 34#%option G_OPT_M_COORDS 35#% description: Input coordinates to reproject 36#% guisection: Input coordinates 37#%end 38#%option G_OPT_F_INPUT 39#% label: Name of input coordinate file 40#% description: '-' for standard input 41#% required: no 42#% guisection: Input coordinates 43#%end 44#%option G_OPT_F_OUTPUT 45#% description: Name for output coordinate file (omit to send to stdout) 46#% required : no 47#% guisection: Output 48#%end 49#%option G_OPT_F_SEP 50#% label: Field separator (format: input[,output]) 51#% required : no 52#% guisection: Input coordinates 53#%end 54#%option 55#% key: proj_in 56#% type: string 57#% description: Input projection parameters (PROJ.4 style) 58#% required : no 59#% guisection: Projections 60#%end 61#%option 62#% key: proj_out 63#% type: string 64#% description: Output projection parameters (PROJ.4 style) 65#% required : no 66#% guisection: Projections 67#%end 68#%flag 69#% key: i 70#% description: Use LL WGS84 as input and current location as output projection 71#% guisection: Projections 72#%end 73#%flag 74#% key: o 75#% description: Use current location as input and LL WGS84 as output projection 76#% guisection: Projections 77#%end 78#%flag 79#% key: d 80#% description: Output long/lat in decimal degrees, or other projections with many decimal places 81#% guisection: Output 82#%end 83#%flag 84#% key: e 85#% description: Include input coordinates in output file 86#% guisection: Output 87#%end 88#%flag 89#% key: c 90#% description: Include column names in output file 91#% guisection: Output 92#%end 93#%rules 94#% required: coordinates, input 95#% exclusive: coordinates, input 96#% exclusive: proj_in, -i 97#% exclusive: proj_out, -o 98#% exclusive: -i, -o 99#%end 100 101import sys 102import os 103import threading 104from grass.script.utils import separator, parse_key_val, encode, decode 105from grass.script import core as gcore 106 107 108class TrThread(threading.Thread): 109 110 def __init__(self, ifs, inf, outf): 111 threading.Thread.__init__(self) 112 self.ifs = ifs 113 self.inf = inf 114 self.outf = outf 115 116 def run(self): 117 while True: 118 line = self.inf.readline() 119 if not line: 120 break 121 line = line.replace(self.ifs, ' ') 122 line = encode(line) 123 self.outf.write(line) 124 self.outf.flush() 125 126 self.outf.close() 127 128 129def main(): 130 coords = options['coordinates'] 131 input = options['input'] 132 output = options['output'] 133 fs = options['separator'] 134 proj_in = options['proj_in'] 135 proj_out = options['proj_out'] 136 ll_in = flags['i'] 137 ll_out = flags['o'] 138 decimal = flags['d'] 139 copy_input = flags['e'] 140 include_header = flags['c'] 141 142 # check for cs2cs 143 if not gcore.find_program('cs2cs'): 144 gcore.fatal(_( 145 "cs2cs program not found, install PROJ first: \ 146 https://proj.org")) 147 148 # parse field separator 149 # FIXME: input_x,y needs to split on multiple whitespace between them 150 if fs == ',': 151 ifs = ofs = ',' 152 else: 153 try: 154 ifs, ofs = fs.split(',') 155 except ValueError: 156 ifs = ofs = fs 157 158 ifs = separator(ifs) 159 ofs = separator(ofs) 160 161 # set up projection params 162 s = gcore.read_command("g.proj", flags='j') 163 kv = parse_key_val(s) 164 if "XY location" in kv['+proj'] and (ll_in or ll_out): 165 gcore.fatal(_("Unable to project to or from a XY location")) 166 167 in_proj = None 168 169 if ll_in: 170 in_proj = "+proj=longlat +datum=WGS84" 171 gcore.verbose( 172 "Assuming LL WGS84 as input, current projection as output ") 173 174 if ll_out: 175 in_proj = gcore.read_command('g.proj', flags='jf') 176 177 if proj_in: 178 if '+' in proj_in: 179 in_proj = proj_in 180 else: 181 gcore.fatal(_("Invalid PROJ.4 input specification")) 182 183 if not in_proj: 184 gcore.verbose("Assuming current location as input") 185 in_proj = gcore.read_command('g.proj', flags='jf') 186 187 in_proj = in_proj.strip() 188 gcore.verbose("Input parameters: '%s'" % in_proj) 189 190 out_proj = None 191 192 if ll_out: 193 out_proj = "+proj=longlat +datum=WGS84" 194 gcore.verbose( 195 "Assuming current projection as input, LL WGS84 as output ") 196 197 if ll_in: 198 out_proj = gcore.read_command('g.proj', flags='jf') 199 200 if proj_out: 201 if '+' in proj_out: 202 out_proj = proj_out 203 else: 204 gcore.fatal(_("Invalid PROJ.4 output specification")) 205 206 if not out_proj: 207 gcore.fatal(_("Missing output projection parameters ")) 208 out_proj = out_proj.strip() 209 gcore.verbose("Output parameters: '%s'" % out_proj) 210 211 # set up input file 212 if coords: 213 x, y = coords.split(',') 214 tmpfile = gcore.tempfile() 215 fd = open(tmpfile, "w") 216 fd.write("%s%s%s\n" % (x, ifs, y)) 217 fd.close() 218 inf = open(tmpfile) 219 else: 220 if input == '-': 221 infile = None 222 inf = sys.stdin 223 else: 224 infile = input 225 if not os.path.exists(infile): 226 gcore.fatal(_("Unable to read input data")) 227 inf = open(infile) 228 gcore.debug("input file=[%s]" % infile) 229 230 # set up output file 231 if not output: 232 outfile = None 233 outf = sys.stdout 234 else: 235 outfile = output 236 outf = open(outfile, 'w') 237 gcore.debug("output file=[%s]" % outfile) 238 239 # set up output style 240 if not decimal: 241 outfmt = ["-w5"] 242 else: 243 outfmt = ["-f", "%.8f"] 244 if not copy_input: 245 copyinp = [] 246 else: 247 copyinp = ["-E"] 248 249 # do the conversion 250 # Convert cs2cs DMS format to GRASS DMS format: 251 # cs2cs | sed -e 's/d/:/g' -e "s/'/:/g" -e 's/"//g' 252 253 cmd = ['cs2cs'] + copyinp + outfmt + \ 254 in_proj.split() + ['+to'] + out_proj.split() 255 256 p = gcore.Popen(cmd, stdin=gcore.PIPE, stdout=gcore.PIPE) 257 258 tr = TrThread(ifs, inf, p.stdin) 259 tr.start() 260 261 if not copy_input: 262 if include_header: 263 outf.write("x%sy%sz\n" % (ofs, ofs)) 264 for line in p.stdout: 265 try: 266 xy, z = decode(line).split(' ', 1) 267 x, y = xy.split('\t') 268 except ValueError: 269 gcore.fatal(line) 270 271 outf.write('%s%s%s%s%s\n' % 272 (x.strip(), ofs, y.strip(), ofs, z.strip())) 273 else: 274 if include_header: 275 outf.write("input_x%sinput_y%sx%sy%sz\n" % (ofs, ofs, ofs, ofs)) 276 for line in p.stdout: 277 inXYZ, x, rest = decode(line).split('\t') 278 inX, inY = inXYZ.split(' ')[:2] 279 y, z = rest.split(' ', 1) 280 outf.write('%s%s%s%s%s%s%s%s%s\n' % 281 (inX.strip(), ofs, inY.strip(), ofs, x.strip(), 282 ofs, y.strip(), ofs, z.strip())) 283 284 p.wait() 285 286 if p.returncode != 0: 287 gcore.warning(_( 288 "Projection transform probably failed, please investigate")) 289 290if __name__ == "__main__": 291 options, flags = gcore.parser() 292 main() 293