A script that worked fine one week ago, has errors and warnings now. I must have changed something.
import grass.script as grass
import numpy as np
def main():
grass.read_command('r.import',input='D:\\dhm75_lcc.img',output='dhm75_lcc', overwrite=True)
f = open("D:\\Glaciers_RGI_60_selected_classified.csv", 'r')
fl = f.readlines()
for x in range(1, len(fl)):
nf = fl[x].split(';')[1] # Get name of file to import
grass.run_command('g.proj', flags='c', georef='D:\\Farinotti_composites\\' + nf + '.tif', location="test")
grass.run_command('g.mapset', location="test", mapset="PERMANENT", overwrite=True)
grass.read_command('r.import', input='D:\\Farinotti_composites\\' + nf + '.tif',output='RGI60_11',overwrite=True)
grass.read_command('r.reclass',input='dhm75_lcc',output='dhm75_lcc_rc',rules='D:\\Elev_Reclass_rules.txt',overwrite=True)
grass.read_command('r.reclass',input='RGI60_11',output= 'RGI60_11_rc',rules='D:\\Thickness_Reclass_rules.txt',overwrite=True)
rs = grass.read_command('r.stats', flags='c', input = ['RGI60_11_rc', 'dhm75_lcc_rc'], separator = ' ')
np.save('D:\\Cell_count_composite_'+nf, rs)
if __name__ == "__main__":
main()
This is the output I get:
WARNING: Projection of dataset does not appear to match current location.
Location PROJ_INFO is:
name: WGS 84 / UTM zone 33N
datum: wgs84
ellps: wgs84
proj: utm
zone: 33
no_defs: defined
Dataset PROJ_INFO is:
name: Lambert Conformal Conic
datum: wgs84
ellps: wgs84
proj: lcc
lat_1: 46
lat_2: 49
lat_0: 47.5
lon_0: 13.3333333333
x_0: 400000
y_0: 400000
towgs84: 0,0,0,-0,-0,-0,0
no_defs: defined
ERROR: proj
Importing raster map <dhm75_lcc>...
Estimated target resolution for input band <dhm75_lcc>: 73.3401507889
Reprojecting <dhm75_lcc>...
ERROR: Unable to create location <test>: File exists
As the error seems to be that I cannot overwrite the location, I tried to make a new location for each file:
So I changed two lines of the script: instead of using location="test", I used location="test"+str(x). As x differs for each file I get different location names for each file, so there is no need for overwriting.
...
for x in range(1, len(fl)):
nf = fl[x].split(';')[1] # Get name of file to import
grass.run_command('g.proj', flags='c', georef='D:\\Farinotti_composites\\' + nf + '.tif', location="test"+str(x))
grass.run_command('g.mapset', location="test"+str(x), mapset="PERMANENT", overwrite=True)
...
Now I get this error:
Importing raster map <dhm75_lcc>...
Estimated target resolution for input band <dhm75_lcc>: 73.3401507889
Reprojecting <dhm75_lcc>...
Location <test1> created
You can switch to the new location by
`g.mapset mapset=PERMANENT location=tets1`
WARNING: Concurrent mapset locking is not supported on Windows
Mapset switched. Your shell continues to use the history for the old mapset
Importing raster map <RGI60_11>...
ERROR: Raster map <dhm75_lcc> not found
What I did now was to move the line where I read the DEM (dhm75_lcc) below the line where I switch to the new location. (Note: as this file is large and takes some time to import, I would prefer to read it just once at the beginning of the script and not in each loop):
...
def main():
f = open("D:\\Glaciers_RGI_60_selected_classified.csv", 'r')
fl = f.readlines()
for x in range(1, len(fl)):
nf = fl[x].split(';')[1] # Get name of file to import
grass.run_command('g.proj', flags='c', georef='D:\\Farinotti_composites\\' + nf + '.tif', location="test"+str(x))
grass.run_command('g.mapset', location="test"+str(x), mapset="PERMANENT", overwrite=True)
grass.read_command('r.import',input='D:\\dhm75_lcc.img',output='dhm75_lcc')
...
The script runs now and produces a result, but (i) there are two warnings, (ii) one strange message that I was not getting before: "0.600000 rounded up to 0" and (iii) the results differ (although not much) from the results I got last week, when there were no warnings.
ERROR: proj
Importing raster map <dhm75_lcc>...
Estimated target resolution for input band <dhm75_lcc>: 73.3401507889
Reprojecting <dhm75_lcc>...
Importing raster map <RGI60_11>...
0.600000 rounded up to 0
Location <tett3> created
You can switch to the new location by
`g.mapset mapset=PERMANENT location=tett3`
WARNING: Concurrent mapset locking is not supported on Windows
Mapset switched. Your shell continues to use the history for the old mapset
WARNING: Projection of dataset does not appear to match current location.
So, why does the projection of dataset does not appear to match current location? Do i need to do something about "Concurrent mapset locking is not supported on Windows"?