vectoraccess.py 2.6 KB

1234567891011121314151617181920212223242526272829303132333435363738394041424344454647484950515253545556575859606162636465666768697071727374757677787980818283848586878889909192939495969798
  1. #!/usr/bin/python
  2. # run within GRASS Spearfish session
  3. # run this before starting python to append module search path:
  4. # export PYTHONPATH=/usr/src/grass70/swig/python
  5. # check with "import sys; sys.path"
  6. # or:
  7. # sys.path.append("/usr/src/grass70/swig/python")
  8. # FIXME: install the grass bindings in $GISBASE/lib/ ?
  9. import os, sys
  10. from grass.lib import grass
  11. from grass.lib import vector as grassvect
  12. import grass.script as grassscript
  13. if not os.environ.has_key("GISBASE"):
  14. print "You must be in GRASS GIS to run this program."
  15. sys.exit(1)
  16. if len(sys.argv)==2:
  17. input = sys.argv[1]
  18. else:
  19. input = raw_input("Vector Map Name? ")
  20. # initialize
  21. grass.G_gisinit('')
  22. # find map in search path
  23. mapset = grass.G_find_vector2(input,'')
  24. # define map structure
  25. map = grassvect.Map_info()
  26. # define open level (level 2: topology)
  27. grassvect.Vect_set_open_level (2)
  28. # open existing map
  29. grassvect.Vect_open_old(map, input, mapset)
  30. # query
  31. print 'Vect map: ', input
  32. print 'Vect is 3D: ', grassvect.Vect_is_3d (map)
  33. print 'Vect DB links: ', grassvect.Vect_get_num_dblinks(map)
  34. print 'Map Scale: 1:', grassvect.Vect_get_scale(map)
  35. # vector box tests
  36. box = grassvect.bound_box()
  37. c_easting1 = 599505.0
  38. c_northing = 4921010.0
  39. c_easting2 = 4599505.0
  40. grassvect.Vect_get_map_box(map, box)
  41. print 'Position 1 in box? ', grassvect.Vect_point_in_box(c_easting1, c_northing, 0, box)
  42. print 'Position 2 in box? ', grassvect.Vect_point_in_box(c_easting2, c_northing, 0, box)
  43. print 'Vector line 2 in box? ', grassvect.Vect_get_line_box(map, 2, box)
  44. # misleading:
  45. # print 'Number of lines:', grassvect.Vect_get_num_lines(map)
  46. # how to access GV_POINT?
  47. # print 'Number of points: ', grassvect.Vect_get_num_primitives(map,GV_POINT)
  48. # confusing:
  49. #print 'Number of lines: ', grassvect.Vect_get_num_primitives(map,GV_LINE)
  50. #print 'Number of areas:', grassvect.Vect_get_num_primitives(map,GV_AREA)
  51. print 'Number of areas:', grassvect.Vect_get_num_areas(map)
  52. layer = 1
  53. tmp = grassscript.tempfile()
  54. tmpf = file(tmp, 'w')
  55. try:
  56. f = grassscript.vector_db(input)[int(layer)]
  57. except KeyError:
  58. grassscript.fatal("There is no table connected to this map. Run v.db.connect or v.db.addtable first.")
  59. table = f['table']
  60. database = f['database']
  61. driver = f['driver']
  62. # call GRASS command
  63. grassscript.run_command('v.db.select', flags = 'c', map = input,
  64. stdout = tmpf)
  65. tmpf.close()
  66. # check if result is empty
  67. tmpf = file(tmp)
  68. if tmpf.read(1) == '':
  69. grassscript.fatal("Table <%s> contains no data.", table)
  70. # print table to stdout
  71. for line in tmpf:
  72. print line
  73. tmpf.close()
  74. os.remove(tmp)
  75. # close map
  76. grassvect.Vect_close(map)
  77. ## end of the python script