Conversion between geographic data formats
convert a shapefile from NZGD2000 (found on LINZ) to WGS84 (different source spacial reference system and output srs)
ogr2ogr -f "ESRI Shapefile" wgs84.shp linz.shp -s_srs EPSG:2193 -t_srs EPSG:4326 -lco ENCODING=UTF-8
convert a shapefile from NZGD2000 to GMT (generic mapping tools) format WGS84 (for use with tools like gmt psxy) (changing output format)
ogr2ogr -f "OGR_GMT" wgs84.gmt linz.shp -s_srs EPSG:2193 -t_srs EPSG:4326
NOTE: GMT (5.3) will not plot segments with missing headers gmt.soest.hawaii.edu/issues/1072 which results in missing segments. GRASS will output all headers when exporting OGR_GMT (v.out.ogr). Alternatively run a command like the following to remove all headers
grep -v '^#' wgs84.gmt > wgs84_nohead.gmt
NOTE: when converting the GMT format to a shapefile, always specify -s_src and -t_src otherwise the .prj file won't be created. This is because the source has no projection information and therefore the destination projection will be unknown.
EPSG:2193 is NZGD2000 http://spatialreference.org/ref/epsg/nzgd2000-new-zealand-transverse-mercator-2000/
EPSG:4326 is WGS84 http://spatialreference.org/ref/epsg/4326/
- ogr2ogr is part of GDAL which is packaged with a whole bunch of useful tools (have a look at included binaries)
Cropping coordinates (station locations) using GMT
Option 1: using gmt select
gmt select -Fpolygon.ll hh400.ll > hh400_land.ll
gmt select -F/nesi/projects/nesi00213/PlottingData/Paths/lds-nz-coastlines-and-islands/150k.gmt non-uniform-100m.ll > non-uniform-land.ll
Option 2: using gmt spatial
Convert points to virtual paths by duplicating position and naming segment with 3rd value (station name).
awk '{print "> " $3 "\n" $1, $2 "\n" $1, $2}' hh400.ll > hh400_paths.llUse GMT to truncate paths (points) by polygon.
gmt spatial -Tpolygon.ll hh400_paths.ll > hh400_cropped.ll
Convert paths back to points.
awk '{if ($1 == ">") n=$2; if (NR % 3 == 2) print $1, $2, n}' hh400_cropped.ll > hh400_land.ll