В продолжение темы Векторизация растра.
Задача
В корневой директории s:/data/MOD14A1/Baikal находятся каталоги, содержащие растровые файлы, формат имени которых AYYYY_DDD.tif (например, A2008_230.tif). Необходимо векторизовать весь этот набор данных, сохранив результат в шейп-файл, причем в шейп файл должны попасть только области растров со значениями 8 и 9. Кроме того, в шейп-файл необходимо добавить 2 поля, содержащих год и дату, которые нужно брать из имен файлов исходных растров.
Вариант решения
В ходе решения использовались библиотека GDAL 1.7.0dev, released 2008/11/26 и Python 2.5.4, установленные с помощью OSGeo4W.
Запускает командную строку и выполняем команды:
Код: Выделить всё
PATH=C:\OSGeo4W\apps\gdal-dev\bin\;%PATH%
set PYTHONPATH=C:\OSGeo4W\apps\gdal-dev\pymod
Код: Выделить всё
import csv
import os
root_path='s://data//MOD14A1//Baikal'
merged=root_path+'//'+'merged.shp'
create_merged=1
for dir in os.walk(root_path):
for file_name in dir[2]:
if file_name[-3:]=='tif':
os.chdir(dir[0])
year=file_name[1:5]
day=str(int(file_name[6:9]))
fname=file_name[0:-4]
step1='gdal_polygonize.py'+' '+file_name+' '+ '-f "ESRI Shapefile"'+' '+fname+'_temp.shp "Active_fires" "fire_mask"'
step2='ogr2ogr -f "ESRI Shapefile" -sql "SELECT * FROM'+' '+fname+'_temp'+' '+ 'WHERE ((fire_mask=8) or (fire_mask=9))"'+' '+fname+'_sel.shp'+' '+fname+'_temp.shp'
step4='ogr2ogr -f CSV outputcsv'+' '+fname+'_sel.shp -lco GEOMETRY=AS_WKT'
step5='del'+' '+fname+'_temp.shp'+' '+fname+'_temp.shx'+' '+fname+'_temp.prj'+' '+fname+'_temp.dbf'+' '+fname+'_sel.shp'+' '+fname+'_sel.shx'+' '+fname+'_sel.prj'+' '+fname+'_sel.dbf'
os.system(step1)
os.system(step2)
os.system(step4)
os.system(step5)
os.chdir(dir[0]+'//'+'outputcsv')
csvnamein=fname+'_sel.csv'
csvnameout=fname+'_out.csv'
f = open(csvnamein, 'rt')
if open(csvnamein).read()!='':
reader=csv.reader(f)
tr=1
f_tmp=open(csvnameout,"wt")
for row in reader:
if tr==1:
f_tmp.write(row[0]+','+row[1]+',Year,Day\n');
tr=0
else:
f_tmp.write('"'+row[0]+'"'+','+row[1]+','+year+','+day+'\n');
f.close()
f_tmp.close()
os.system('ogr2ogr'+' '+fname+'.shp'+' '+csvnameout)
os.system('del'+' '+csvnamein+' '+csvnameout)
os.system('move *.* ..//')
os.chdir('..//')
os.system('rmdir outputcsv')
if create_merged==1:
os.system('ogr2ogr'+' '+merged+' '+fname+'.shp')
create_merged=0
else:
step6='ogr2ogr -update -append'+' '+merged+' '+fname+'.shp'+' '+'-nln merged'
os.system(step6)
else:
f.close()
os.chdir('..//')
os.system('rmdir outputcsv /q /s')
P.S.: Для тех, кто хочет поэкспериментировать - один из растров во вложении.