1

I have a global raster of railways networks in Mollweide CRS (resolution 1000x1000) and I want to reproject it to Europe (EPSG:3035), crop it to Europe, and change the resolution to 30x30 km (0.25°). I downloaded the raster map from here https://datadryad.org/stash/dataset/doi:10.5061/dryad.3tx95x6d9 (folder Human_footprint_pressure_layers, size 438.71 MB, file railways.tif).
What I would usually do would be to reproject it (externally from GRASS GIS) using gdalwarp, first to EPSG:4326 (as I still have this problem and, apparently, first transforming from EPSG:54009 to EPSG:4326 and then do the rest works) with:

gdalwarp -overwrite -t_srs EPSG:4326 -tr 0.25 0.25 /mypath/railways.tif /mypath/railways_WGS84.tif

Then I access GRASS, enter my EPSG:4326 mapset and load the file:

grass /mypath/grass/4326/PERMANENT 
g.region -p res=0.25
r.in.gdal --o input=/mypath/railways_WGS84.tif output=railways_WGS84

Then I'd switch mapset to my EPSG:3035 one, create a raster mask of Europe, and reproject from EPSG:4326 to EPSG:3035, having therefore the railways_WGS84 reprojected and clipped to Europe:

g.mapset mapset=PERMANENT location=3035
g.region -p raster=EU_3035 res=30000
r.mask --o raster=EU_3035
r.proj --o location=4326 mapset=PERMANENT input=railways_WGS84 output=railways_3035
r.out.gdal --o input=railways_3035 output=/mypath/railways_3035.tif

This creates a raster map with desired resolution and CRS but 255 colour categories, and the whole raster looks the same (it doesn't have different values depending on the part of Europe), like this:

enter image description here

I think what is failing here is the masking part, because until the reprojection from railways_WGS84 to railways_3035 and before the masking part, the files look fine, like this:

enter image description here

I tried also r.mask.rast in QGIS but produced the same result:

processing.run("grass7:r.mask.rast", {'raster':'mypath/EU_3035.tif','input':'mypath/railways_3035.tif','maskcats':'1','-i':False,'output':'TEMPORARY_OUTPUT','GRASS_REGION_PARAMETER':None,'GRASS_REGION_CELLSIZE_PARAMETER':0,'GRASS_RASTER_FORMAT_OPT':'','GRASS_RASTER_FORMAT_META':''})

I am using GRASS GIS 7.8.7 from a Linux server and QGIS 3.28.2 from a Windows 10 local machine.

4
  • Why not just reproject from Mollewide directly to EPSG:3035 with gdalwarp? The gdal utility also offers a -te flag to set "target extent", that you can use to get just the region of Europe. Or you can get a polygon surrounding Europe, and use the gdalwarp parameter -cutline to clip out just the countries in Europe. ??
    – Micha
    Commented Jun 22, 2023 at 13:33
  • Hi Micha, apparently I cannot reproject directly from Mollweide to EPSG:3035 (see gis.stackexchange.com/questions/437018/…) and that is the only workaround I found that worked for me. I will definitely try gdalwarp -cutline tho. Will it be enough to obtain the lines of my Europe vector?
    – LT17
    Commented Jun 22, 2023 at 13:50
  • I ws able to run: gdalwarp -t_srs EPSG:3035 -te 2400000 1380000 7650000 5600000 railways.tif railways_EU.tif with no problem. Got the railways raster, projected to 3035, and clipped to EU.
    – Micha
    Commented Jun 22, 2023 at 14:08
  • I can run it without problems as well, the thing is that with those boundaries I would get also the Northern part of Africa and other areas I don't need. Finally I got the desired result with gdalwarp -overwrite -t_srs EPSG:3035 -tr 10000 10000 -of GTiff -cutline EU_3035.shp -crop_to_cutline railways_3035.tif railways_3035_crop.tif Thanks :)
    – LT17
    Commented Jun 23, 2023 at 9:41

0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.