d.correlate.py 2.9 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106
  1. #!/usr/bin/env python
  2. #
  3. ############################################################################
  4. #
  5. # MODULE: d.correlate for GRASS 6; based on dcorrelate.sh for GRASS 4,5
  6. # AUTHOR(S): CERL - Michael Shapiro; updated to GRASS 6 by Markus Neteler 5/2005
  7. # Converted to Python by Glynn Clements
  8. # PURPOSE: prints a graph of the correlation between data layers (in pairs)
  9. # derived from <grass5>/src.local/d.correlate.sh
  10. # COPYRIGHT: (C) 2005, 2008 by the GRASS Development Team
  11. #
  12. # This program is free software under the GNU General Public
  13. # License (>=v2). Read the file COPYING that comes with GRASS
  14. # for details.
  15. #
  16. # TODO GRASS 7: rename parameters to map1, map2 or better multiple map=map1[,map2[...]]
  17. #############################################################################
  18. #%Module
  19. #% description: Prints a graph of the correlation between data layers (in pairs).
  20. #% keywords: display
  21. #% keywords: diagram
  22. #%End
  23. #%option
  24. #% key: layers
  25. #% type: string
  26. #% gisprompt: old,cell,raster
  27. #% description: raster input map
  28. #% required: yes
  29. #% multiple: yes
  30. #%end
  31. import sys
  32. import os
  33. from grass.script import core as grass
  34. def main():
  35. layers = options['layers'].split(',')
  36. if len(layers) < 2:
  37. grass.error("At least 2 layers are required")
  38. tmpfile = grass.tempfile()
  39. for map in layers:
  40. if not grass.find_file(map, element = 'cell')['file']:
  41. grass.fatal(_("Input map <%s> not found") % map)
  42. grass.write_command('d.text', color = 'black', size = 4, line = 1, stdin = "CORRELATION")
  43. colors = "red black blue green gray violet".split()
  44. line = 2
  45. iloop = 0
  46. jloop = 0
  47. for iloop, i in enumerate(layers):
  48. for jloop, j in enumerate(layers):
  49. if i != j and iloop <= jloop:
  50. color = colors[0]
  51. colors = colors[1:]
  52. colors.append(color)
  53. grass.write_command('d.text', color = color, size = 4, line = line, stdin = "%s %s" % (i, j))
  54. line += 1
  55. ofile = file(tmpfile, 'w')
  56. grass.run_command('r.stats', flags = 'cnA', input = (i, j), stdout = ofile)
  57. ofile.close()
  58. ifile = file(tmpfile, 'r')
  59. first = True
  60. for l in ifile:
  61. f = l.rstrip('\r\n').split(' ')
  62. x = float(f[0])
  63. y = float(f[1])
  64. if first:
  65. minx = maxx = x
  66. miny = maxy = y
  67. first = False
  68. if minx > x: minx = x
  69. if maxx < x: maxx = x
  70. if miny > y: miny = y
  71. if maxy < y: maxy = y
  72. ifile.close()
  73. kx = 100.0/(maxx-minx+1)
  74. ky = 100.0/(maxy-miny+1)
  75. p = grass.feed_command('d.graph', color = color)
  76. ofile = p.stdin
  77. ifile = file(tmpfile, 'r')
  78. for l in ifile:
  79. f = l.rstrip('\r\n').split(' ')
  80. x = float(f[0])
  81. y = float(f[1])
  82. ofile.write("icon + 0.1 %f %f\n" % ((x-minx+1) * kx, (y-miny+1) * ky))
  83. ifile.close()
  84. ofile.close()
  85. p.wait()
  86. grass.try_remove(tmpfile)
  87. if __name__ == "__main__":
  88. options, flags = grass.parser()
  89. main()