PyRaster:Python Tools for Processing Geospatial Raster Data

是 否 +2 论坛币

k人 参与回答

经管之家送您一份

应届毕业生专属福利!

求职就业群

赵安豆老师微信：zhaoandou666

经管之家联合CDA

送您一个全额奖学金名额~ !

立即领取

感谢您参与论坛问题回答

经管之家送您两个论坛币！

+2 论坛币

PyRaster is a collection of Python tools for processing geospatial raster data as NumPy arrays.

本帖隐藏的内容

PyRaster-master.zip (79.26 KB)

2016-9-12 08:05:49 上传

扫码加我 拉你入群

请注明：姓名-公司-职位

以便审核进群资格，未注明则拒绝

分享0 收藏1 回帖

关键词：Processing processI Process Spatial python collection

1. #!/usr/bin/env python
   
2. # Import standard modules and rasterIO
   
3. import sys, os, string, rasterIO
   
4. 
   
5. # Get a list of files in current directory
   
6. # Create a loop for all files in current directory
   
7. 
   
8. for file in flist:
   
9.         # Check file type (in this case Geotiff)   
   
10.         if file.endswith('.tif'):
    
11.         
    
12.         # Open a pointer to the file
    
13.         pointer = rasterIO.opengdalraster(file)
    
14.         
    
15.         # Read the raster metadata for the new output file
    
16.         driver, XSize, YSize, proj_wkt, geo_t_params = rasterIO.readrastermeta(pointer)
    
17.         
    
18.         # Read the first band to a matrix called band_1
    
19.         band_1 =rasterIO.readrasterband(pointer, 1)
    
20.         
    
21.         # Read the second band to a matrix called band_2
    
22.         band_2 = rasterIO.readrasterband(pointer, 2)
    
23.         
    
24.         # Perform the NDVI calculation and put the results into a new matrix
    
25.         new_ndvi_band = ((band_2 - band_1) / (band_2 + band_1))
    
26.         # Get the input file filename without extension and create a new file name
    
27.         parts = string.split(file)
    
28.         newname ='./'+parts[0]+'_ndvi.tif' # ./filename_ndvi.tif
    
29.         # Get the EPSG code from well known text projection
    
30.         epsg =rasterIO.wkt2epsg(proj_wkt)
    
31.         
    
32.         # Write the NDVI matrix to a new raster file
    
33.         rasterIO.writerasterband(new_ndvi_band, newname, XSize, YSize, geo_t_params, epsg)
    
34.    
    
35. # loop will now go to next file in input list

复制代码

1. #!/usr/bin/env python
   
2. # Script to plot raster image histogram
   
3. import sys, os, string, rasterIO
   
4. import numpy.ma as ma
   
5. import matplotlib.pyplot as plt
   
6. 
   
7. # Create a histogram function
   
8. def rasterHistogram(raster_matrix):
   
9.     '''Accepts matrix and generates histogram'''
   
10.      
    
11.     # Line above is function docstring
    
12.     # Flatten 2d-matrix
    
13.     flat_raster = ma.compressed(raster_matrix)
    
14.      
    
15.     # Setup the plot (see matplotlib.sourceforge.net)
    
16.     fig = plt.figure(figsize=(8,11))
    
17.     ax = fig.add_subplot(1,1,1)
    
18.      
    
19.     # Plot histogram
    
20.     ax.hist(flat_raster, 10, normed=0, histtype='bar',
    
21.             align=mid)
    
22.     # Show the plot on screen
    
23.     plt.show()
    
24.      
    
25. # Open a raster using rasterIO
    
26. #...
    
27. # Pass matrix representation of raster to function
    
28. rasterHistogram(band)
    
29. # Shows histogram on screen

复制代码

1. #!/usr/bin/env python
   
2. # Batch processing script with multiprocessing enabled
   
3. 
   
4. import sys, os, string, rasterIO
   
5. from multiprocessing import Process, Queue
   
6. 
   
7. # Create a function for the process
   
8. 
   
9. def processfunction(file_list):
   
10.             
    
11.         # Create a loop for all files in current directory
    
12. 
    
13.             
    
14.         for file in file_list:
    
15.             
    
16.                 # Check file type (in this case Geotiff)
    
17.                 if file.endswith('.tif'):
    
18.                             
    
19.                 # Open a pointer to the file
    
20.                 pointer = rasterIO.opengdalraster(file)
    
21.                             
    
22.                 # Read the raster metadata for the new output file
    
23.                 driver, XSize, YSize, proj_wkt, geo_t_params = rasterIO.readrastermeta(pointer)
    
24.                             
    
25.                 # Read the first band to a matrix called band_1
    
26.                 band_1 = rasterIO.readrasterband(pointer, 1)
    
27.                             
    
28.                 # Read the second band to a matrix called band_2
    
29. 
    
30.                 band_2 =rasterIO.readrasterband(pointer, 2)
    
31.                             
    
32.                 # Perform the NDVI calculation and put the results into a new matrix
    
33.                 new_ndvi_band = ((band_2 - band_1) / (band_2 + band_1))
    
34.                             
    
35.                 # Get the input file filename without extension and create a new file name
    
36.                 parts =string.split(file)
    
37.                 newname = './'+parts[0]+'_ndvi.tif' # ./filename_ndvi.tif
    
38. 
    
39.         # Get the EPSG code from well known text projection
    
40.         epsg = rasterIO.wkt2epsg(proj_wkt)
    
41.        
    
42.         # Write the NDVI matrix to a new raster file
    
43.         rasterIO.writerasterband(new_ndvi_band, newname, XSize, YSize, geo_t_params, epsg)
    
44.    
    
45. # loop will now go to next file in input list
    
46. # Create a run function to accept jobs from queue
    
47. def processRun(q, filename):
    
48.         q.put(processfunc(filename))
    
49. 
    
50. # Create the main function
    
51. def main(arg=sys.argv):
    
52.         # Get a list fo files in the current directory (assumed to contain raster data)
    
53.         file_list = os.listdir(os.getcwd())
    
54.             
    
55.         # Get the length of the file list
    
56.         len_flist = len(file_list)
    
57.         half_len = len_flist / 2
    
58. 
    
59.         # Create a queue object
    
60.         q = Queue()
    
61.             
    
62.         # Create two processes (add more depending on processor availability)
    
63.         p1 = Process(target=processRun, args=(q, flist[:half_len]))
    
64.         p2 = Process(target=processRun, args=(q, flist[half_len:]))
    
65. 
    
66.         # Start processes
    
67.         p1.start()
    
68.         p2.start()
    
69. 
    
70. # Standard Python script execution/exit handling
    
71. 
    
72. if __name__=='__main__':
    
73.         sys.exit(main())

复制代码

1. #The numerical Python package Numpy contains several functions to optimise matrix
   
2. #operations. Using the numpy.vectorize() function it is possible to create an
   
3. #optimised matrix function from a standard function.
   
4. 
   
5. # Test each item ('pixel') in a 2d-matrix
   
6. #...
   
7. # first define a standard non-matrix function...
   
8. def value_test(a, b):
   
9.     if a <= b:
   
10.             return 1
    
11.         else:
    
12.             return 0
    
13. # now vectorize the function
    
14. value_test_vect = numpy.vectorize(value_test)
    
15. # now pass the new vectorized function the matrix
    
16. new_matrix = value_test_vect(matrix, some_value)
    
17. #...

复制代码

1. #Read a raster file
   
2. 
   
3. # Import the rasterIO module
   
4. import rasterIO
   
5. 
   
6. # Open a gdal pointer to a file on disk
   
7. file_pointer = rasterIO.opengdalraster(
   
8.     '/path/surface_temperature_kelvin.tif')
   
9. 
   
10. # Get the spatial reference information from the image
    
11. drive, XSize, YSize, proj_wkt, geo_t_params =
    
12.     rasterIO.readrastermeta(pointer)
    
13. 
    
14. # Read a band from the file to a new matrix in memory,
    
15.     #the number indicates which band to read.
    
16. kelvin_band = rasterIO.readrasterband(file_pointer, 1)

复制代码

谢谢楼主，看看是什么！

hhahhaaaaaaaaaaaaaaa

谢谢！

GIS for Statistic