Package osgeo :: Package utils :: Module gcps2vec |
[hide private]
[frames] | no frames]

# Source Code for Module osgeo.utils.gcps2vec

      1  #!/usr/bin/env python3
      2  # ******************************************************************************
      3  #  $Id$
      4  #
      5  #  Project:  GDAL
      6  #  Purpose:  Convert GCPs to a point layer.
      7  #  Author:   Frank Warmerdam, warmerdam@pobox.com
      8  #
      9  # ******************************************************************************
     10  #  Copyright (c) 2005, Frank Warmerdam <warmerdam@pobox.com>
     11  #
     12  #  Permission is hereby granted, free of charge, to any person obtaining a
     13  #  copy of this software and associated documentation files (the "Software"),
     14  #  to deal in the Software without restriction, including without limitation
     15  #  the rights to use, copy, modify, merge, publish, distribute, sublicense,
     16  #  and/or sell copies of the Software, and to permit persons to whom the
     17  #  Software is furnished to do so, subject to the following conditions:
     18  #
     19  #  The above copyright notice and this permission notice shall be included
     20  #  in all copies or substantial portions of the Software.
     21  #
     22  #  THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
     23  #  OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
     24  #  FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL
     25  #  THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
     26  #  LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING
     27  #  FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER
     28  #  DEALINGS IN THE SOFTWARE.
     29  # ******************************************************************************
     30
     31  import sys
     32
     33  from osgeo import gdal
     34  from osgeo import ogr
     35  from osgeo import osr
     36
     37
     38 -def Usage():
     39      print('Usage: gcps2vec.py [-of <ogr_drivername>] [-p] <raster_file> <vector_file>')
     40      sys.exit(1)
     41
     42
     43 -def main(argv):
     44      out_format = 'GML'
     45      in_file = None
     46      out_file = None
     47      pixel_out = 0
     48
     49      gdal.AllRegister()
     50      argv = gdal.GeneralCmdLineProcessor(argv)
     51      if argv is None:
     52          sys.exit(0)
     53
     54      # Parse command line arguments.
     55      i = 1
     56      while i < len(argv):
     57          arg = argv[i]
     58
     59          if arg == '-of':
     60              i = i + 1
     61              out_format = argv[i]
     62
     63          elif arg == '-p':
     64              pixel_out = 1
     65
     66          elif in_file is None:
     67              in_file = argv[i]
     68
     69          elif out_file is None:
     70              out_file = argv[i]
     71
     72          else:
     73              Usage()
     74
     75          i = i + 1
     76
     77      if out_file is None:
     78          Usage()
     79
     80      # ----------------------------------------------------------------------------
     81      # Open input file, and fetch GCPs.
     82      # ----------------------------------------------------------------------------
     83      ds = gdal.Open(in_file)
     84      if ds is None:
     85          print('Unable to open %s' % in_file)
     86          sys.exit(1)
     87
     88      gcp_srs = ds.GetGCPProjection()
     89      gcps = ds.GetGCPs()
     90
     91      ds = None
     92
     93      if gcps is None or not gcps:
     94          print('No GCPs on file %s!' % in_file)
     95          sys.exit(1)
     96
     97      # ----------------------------------------------------------------------------
     98      # Create output file, and layer.
     99      # ----------------------------------------------------------------------------
    100
    101      drv = ogr.GetDriverByName(out_format)
    102      if drv is None:
    103          print('No driver named %s available.' % out_format)
    104          sys.exit(1)
    105
    106      ds = drv.CreateDataSource(out_file)
    107
    108      if pixel_out == 0 and gcp_srs != "":
    109          srs = osr.SpatialReference()
    110          srs.ImportFromWkt(gcp_srs)
    111      else:
    112          srs = None
    113
    114      if pixel_out == 0:
    115          geom_type = ogr.wkbPoint25D
    116      else:
    117          geom_type = ogr.wkbPoint
    118
    119      layer = ds.CreateLayer('gcps', srs, geom_type=geom_type)
    120
    121      if pixel_out == 0:
    122          fd = ogr.FieldDefn('Pixel', ogr.OFTReal)
    123          layer.CreateField(fd)
    124
    125          fd = ogr.FieldDefn('Line', ogr.OFTReal)
    126          layer.CreateField(fd)
    127      else:
    128          fd = ogr.FieldDefn('X', ogr.OFTReal)
    129          layer.CreateField(fd)
    130
    131          fd = ogr.FieldDefn('Y', ogr.OFTReal)
    132          layer.CreateField(fd)
    133
    134          fd = ogr.FieldDefn('Z', ogr.OFTReal)
    135          layer.CreateField(fd)
    136
    137      fd = ogr.FieldDefn('Id', ogr.OFTString)
    138      layer.CreateField(fd)
    139
    140      fd = ogr.FieldDefn('Info', ogr.OFTString)
    141      layer.CreateField(fd)
    142
    143      # ----------------------------------------------------------------------------
    144      # Write GCPs.
    145      # ----------------------------------------------------------------------------
    146
    147      for gcp in gcps:
    148
    149          feat = ogr.Feature(layer.GetLayerDefn())
    150
    151          if pixel_out == 0:
    152              geom = ogr.Geometry(geom_type)
    153              feat.SetField('Pixel', gcp.GCPPixel)
    154              feat.SetField('Line', gcp.GCPLine)
    155              geom.SetPoint(0, gcp.GCPX, gcp.GCPY, gcp.GCPZ)
    156          else:
    157              geom = ogr.Geometry(geom_type)
    158              feat.SetField('X', gcp.GCPX)
    159              feat.SetField('Y', gcp.GCPY)
    160              feat.SetField('Z', gcp.GCPZ)
    161              geom.SetPoint(0, gcp.GCPPixel, gcp.GCPLine)
    162
    163          feat.SetField('Id', gcp.Id)
    164          feat.SetField('Info', gcp.Info)
    165
    166          feat.SetGeometryDirectly(geom)
    167          layer.CreateFeature(feat)
    168
    169      feat = None
    170      ds.Destroy()
    171
    172
    173  if __name__ == '__main__':
    174      sys.exit(main(sys.argv))
    175