Provided by: gmt-common_5.4.3+dfsg-1_all bug

NAME

       gmtspatial - Do geospatial operations on lines and polygons

SYNOPSIS

       gmtspatial  [ table ] [  -A[amin_dist][unit]] [  -C ] [  -D[+ffile][+aamax][+ddmax][+c|Ccmax][+sfact] ] [
       -E+|-    ]    [     -F[l]    ]    [     -I[e|i]     ]     [      -Npfile[+a][+pstart][+r][+z]     ]     [
       -Q[[-|+]unit][+cmin[/max]][+h][+l][+p][+s[a|d]]  ]  [   -Rregion  ] [  -Si|u|s|j ] [  -T[clippolygon] ] [
       -V[level] ] [ -bbinary ] [ -dnodata ] [ -eregexp ] [ -fflags ] [ -ggaps ] [ -hheaders ]  [  -iflags  ]  [
       -oflags ] [ -:[i|o] ]

       Note: No space is allowed between the option flag and the associated arguments.

DESCRIPTION

       gmtspatial  reads  one or more data files (which may be multisegment files) that contains closed polygons
       and operates of these polygons in the specified way.  Operations  include  area  calculation,  handedness
       reversals, and polygon intersections.

REQUIRED ARGUMENTS

       None.

OPTIONAL ARGUMENTS

       table  One  or  more  ASCII (or binary, see -bi[ncols][type]) data table file(s) holding a number of data
              columns. If no tables are given then we read from standard input.

       -A[amin_dist][unit]
              Perform spatial nearest neighbor (NN) analysis: Determine the nearest neighbor of each  point  and
              report  the NN distances and the point IDs involved in each pair (IDs are the input record numbers
              starting at 0).  Use -Aa to decimate a data set so that no NN distance is lower than the threshold
              min_dist.  In this case we write out the  (possibly  averaged)  coordinates  and  the  updated  NN
              distances  and  point  IDs.   A  negative  point number means the original point was replaced by a
              weighted average (the absolute ID value gives the ID of the first original point ID to be included
              in the average.).  Note: The input data are assumed to contain (lon, lat) or  (x,  y),  optionally
              followed  by a z and a weight [1] column.  We compute a weighted average of the location and z (if
              present).

       -C     Clips polygons to the map region, including map boundary to the polygon as needed. The result is a
              closed polygon (see -T for truncation instead). Requires -R.

       -D[+ffile][+aamax][+ddmax][+c|Ccmax][+sfact]
              Check for duplicates among the input lines or polygons, or, if file is given via +f, check if  the
              input  features  already  exist  among  the features in file. We consider the cases of exact (same
              number and coordinates) and approximate matches (average distance between nearest  points  of  two
              features  is  less  than a threshold). We also consider that some features may have been reversed.
              Features are considered approximate matches if their minimum distance is less than dmax  [0]  (see
              UNITS) and their closeness (defined as the ratio between the average distance between the features
              divided  by  their  average length) is less than cmax [0.01]. For each duplicate found, the output
              record begins with the single letter Y (exact match) or ~ (approximate match). If the two matching
              segments differ in length by more than a factor of 2 then we consider the duplicate to be either a
              subset (-) or a superset (+). Finally, we also note if two lines are the  result  of  splitting  a
              continuous  line across the Dateline (|).  For polygons we also consider the fractional difference
              in areas; duplicates must differ by less than amax [0.01]. By default, we compute  the  mean  line
              separation.  Use  +Ccmin  to  instead  compute  the  median line separation and therefore a robust
              closeness value. Also by default we consider all distances between points on one line and another.
              Append +p to limit the comparison to points that project perpendicularly to points  on  the  other
              line (and not its extension).

       -E+|- ]
              Reset  the  handedness  of all polygons to match the given + (counter-clockwise) or - (clockwise).
              Implies -Q+.

       -F[l]  Force input data to become polygons on output, i.e., close them explicitly if not already  closed.
              Optionally, append l to force line geometry.

       -I[e|i]
              Determine  the  intersection  locations  between  all pairs of polygons.  Append i to only compute
              internal (i.e., self-intersecting polygons) crossovers  or  e  to  only  compute  external  (i.e.,
              between paris of polygons) crossovers [Default is both].

       -Npfile[+a][+pstart][+r][+z]
              Determine  if one (or all, with +a) points of each feature in the input data are inside any of the
              polygons given in the pfile. If inside, then report which polygon it is; the polygon ID is  either
              taken  from  the aspatial value assigned to Z, the segment header (first -Z, then -L are scanned),
              or it is assigned the running number that is initialized  to  start  [0].  By  default  the  input
              segment that are found to be inside a polygon are written to stdout with the polygon ID encoded in
              the  segment  header  as  -ZID.  Alternatively,  append +r to just report which polygon contains a
              feature or +z to have the IDs added as an extra data column on output. Segments that  fail  to  be
              inside  a  polygon are not written out. If more than one polygon contains the same segment we skip
              the second (and further) scenario.

       -Q[[-|+]unit][+cmin[/max]][+h][+l][+p][+s[a|d]]
              Measure the area of all polygons or length of line segments. Use -Q+h to append the area  to  each
              polygons  segment  header [Default simply writes the area to stdout]. For polygons we also compute
              the centroid location while for line data we compute the mid-point (half-length) position.  Append
              a distance unit to select the unit used (see UNITS). Note that the area will depend on the current
              setting of PROJ_ELLIPSOID; this should be a recent ellipsoid to get accurate results. The centroid
              is  computed using the mean of the 3-D Cartesian vectors making up the polygon vertices, while the
              area is obtained via an equal-area projection.  For line lengths you may prepend -|+ to  the  unit
              and  the  calculation  will  use Flat Earth or Geodesic algorithms, respectively [Default is great
              circle distances].  Normally, all input segments will  be  be  reflected  on  output.   Use  c  to
              restrict  processing  to those whose length (or area for polygons) fall inside the specified range
              set by min and max.  If max is not set it defaults to infinity.  To sort  the  segments  based  on
              their  lengths  or  area, use s and append a for ascending and d for descending order [ascending].
              By default, we consider open polygons as lines.   Append  +p  to  close  open  polygons  and  thus
              consider all input as polygons, or append +l to consider all input as lines, even if closed.

       -Rwest/east/south/north[/zmin/zmax][+r][+uunit]
              west,  east,  south, and north specify the region of interest, and you may specify them in decimal
              degrees or in [±]dd:mm[:ss.xxx][W|E|S|N] format Append +r  if  lower  left  and  upper  right  map
              coordinates  are  given instead of w/e/s/n. The two shorthands -Rg and -Rd stand for global domain
              (0/360 and -180/+180 in longitude respectively, with -90/+90 in latitude).  Alternatively for grid
              creation, give Rcodelon/lat/nx/ny, where code is a 2-character combination of L, C, R  (for  left,
              center, or right) and T, M, B for top, middle, or bottom. e.g., BL for lower left.  This indicates
              which  point  on a rectangular region the lon/lat coordinate refers to, and the grid dimensions nx
              and ny with grid spacings via -I is used  to  create  the  corresponding  region.   Alternatively,
              specify  the  name  of an existing grid file and the -R settings (and grid spacing, if applicable)
              are copied from the grid. Appending +uunit expects projected  (Cartesian)  coordinates  compatible
              with  chosen  -J  and we inversely project to determine actual rectangular geographic region.  For
              perspective view (-p), optionally append /zmin/zmax.  In case of perspective view (-p), a  z-range
              (zmin,  zmax)  can  be  appended  to indicate the third dimension. This needs to be done only when
              using the -Jz option, not when using only the -p option. In the latter case a perspective view  of
              the  plane  is  plotted,  with no third dimension. Clips polygons to the map region, including map
              boundary to the polygon as needed. The result is a closed polygon.

       -Si|j|s|u
              Spatial processing of polygons. Choose  from  -Si  which  returns  the  intersection  of  polygons
              (closed),  -Su  which  returns  the union of polygons (closed), -Ss which will split polygons that
              straddle the Dateline, and -Sj which will join polygons that were split by  the  Dateline.   Note:
              Only -Ss has been implemented.

       -T[clippolygon]
              Truncate  polygons against the specified polygon given, possibly resulting in open polygons. If no
              argument is given to -T we create a clipping polygon from -R which then  is  required.  Note  that
              when  the  -R  clipping  is in effect we will also look for polygons of length 4 or 5 that exactly
              match the -R clipping polygon.

       -V[level] (more …)
              Select verbosity level [c].

       -bi[ncols][t] (more …)
              Select native binary input. [Default is 2 input columns].

       -bo[ncols][type] (more …)
              Select native binary output. [Default is same as input].

       -d[i|o]nodata (more …)
              Replace input columns that equal nodata with NaN and do the reverse on output.

       -e[~]”pattern” | -e[~]/regexp/[i] (more …)
              Only accept data records that match the given pattern.

       -f[i|o]colinfo (more …)
              Specify data types of input and/or output columns.

       -g[a]x|y|d|X|Y|D|[col]z[+|-]gap[u] (more …)
              Determine data gaps and line breaks.

       -h[i|o][n][+c][+d][+rremark][+rtitle] (more …)
              Skip or produce header record(s).

       -icols[+l][+sscale][+ooffset][,] (more …)
              Select input columns and transformations (0 is first column).

       -ocols[,…] (more …)
              Select output columns (0 is first column).

       -:[i|o] (more …)
              Swap 1st and 2nd column on input and/or output.

       -^ or just -
              Print a short message about the syntax of the command, then exits (NOTE: on Windows just use -).

       -+ or just +
              Print an extensive usage (help) message, including the explanation of any  module-specific  option
              (but not the GMT common options), then exits.

       -? or no arguments
              Print a complete usage (help) message, including the explanation of all options, then exits.

UNITS

       For  map  distance  unit,  append unit d for arc degree, m for arc minute, and s for arc second, or e for
       meter [Default], f for foot, k for km, M for statute mile, n for nautical mile, and u for US survey foot.
       By default we compute such distances using a spherical approximation with great circles. Prepend -  to  a
       distance  (or  the  unit  is no distance is given) to perform “Flat Earth” calculations (quicker but less
       accurate) or prepend + to perform exact geodesic calculations (slower but more accurate).

ASCII FORMAT PRECISION

       The ASCII output formats of numerical data are controlled by parameters in your gmt.conf file.  Longitude
       and  latitude  are  formatted  according  to  FORMAT_GEO_OUT,  absolute  time  is  under  the  control of
       FORMAT_DATE_OUT and FORMAT_CLOCK_OUT, whereas general floating point values are  formatted  according  to
       FORMAT_FLOAT_OUT. Be aware that the format in effect can lead to loss of precision in ASCII output, which
       can  lead  to  various  problems downstream. If you find the output is not written with enough precision,
       consider  switching  to  binary  output  (-bo  if  available)  or  specify  more   decimals   using   the
       FORMAT_FLOAT_OUT setting.

EXAMPLE

       To turn all lines in the multisegment file lines.txt into closed polygons, run

              gmt spatial lines.txt -F > polygons.txt

       To compute the area of all geographic polygons in the multisegment file polygons.txt, run

              gmt spatial polygons.txt -Q > areas.txt

       Same  data,  but  now  orient  all  polygons to go counter-clockwise and write their areas to the segment
       headers, run

              gmt spatial polygons.txt -Q+h -E+ > areas.txt

       To determine the areas of  all  the  polygon  segments  in  the  file  janmayen_land_full.txt,  add  this
       information  to  the  segment  headers,  sort the segments from largest to smallest in area but only keep
       polygons with area larger than 1000 sq. meters, run

              gmt spatial -Qe+h+p+c1000+sd -V janmayen_land_full.txt > largest_pols.txt

       To determine the intersections between the polygons A.txt and B.txt, run

              gmt spatial A.txt B.txt -Ie > crossovers.txt

       To truncate polygons A.txt against polygon B.txt, resulting in an open line segment, run

              gmt gmtspatial A.txt -TB.txt > line.txt

NOTES

       OGR/GMT files are considered complete datasets and thus you cannot specify more than one at a given time.
       This causes problems if you want to examine the intersections of two OGR/GMT files.  The solution  is  to
       convert them to regular datasets via gmtconvert and then run gmtspatial on the converted files.

SEE ALSO

       gmt, gmtconvert, gmtselect, gmtsimplify

COPYRIGHT

       2018, P. Wessel, W. H. F. Smith, R. Scharroo, J. Luis, and F. Wobbe

5.4.3                                             Jan 03, 2018                                  GMTSPATIAL(1gmt)