v.surf.rst does not respect mask

[cmbarton]

When the mask argument is set in the v.surf.rst command, it gives an error and does not produce an interpolated surface. The error is:

ERROR: Input window changed while maps are open for read

I have not checked to see if this affects 6.4.3 too

[helena]

I can confirm this ERROR when mask parameter is given in GRASS7. Also, the output of rescaled tension fi which used to be written out with -t flag has disappeared.

I just checked in GRASS6.4.3 and mask there works as intended - user can define it as parameter mask= or by setting MASK for the current mapset.

Helena

[annakrat]

When you have a look at ​the code calling fatal error, there is some additional functionality after G_fatal_error which is nonsense. Maybe there should be a warning instead or the following line should be removed to avoid confusion.

[glynn]

Replying to annakrat:

When you have a look at ​the code calling fatal error, there is some additional functionality after G_fatal_error which is nonsense. Maybe there should be a warning instead or the following line should be removed to avoid confusion.

That was inadvertently left behind in r42876. Fixed in r53920.

The core issue is that v.surf.rst (or possibly the RST library itself) needs to be updated to correctly handle a split window (different windows for read and write). I updated most of the affected modules when the change was introduced, but I didn't understand v.surf.rst well enough to fix it.

GRASS 7 doesn't allow the window to be changed while maps are open.

Most of the modules which did this were resampling modules (r.resamp.*, r.proj) which used one window for reading and one for writing, alternating between the two for each row. This was highly inefficient, as changing the window recalculates the column mapping for each open map.

So it was changed to have one window for reading and one for writing. Both windows must be set before any maps are opened, and cannot be changed thereafter. In order to catch any bugs, any incompatible or ambiguous behaviour currently results in a fatal error.

One consequence of this change is that it isn't possible to read multiple maps at different resolutions. If this turns out to be necessary, lib/raster could be changed e.g. to allow a specific window to be set for a specific map. But I don't want to do that unless it's strictly necessary, as it introduces the potential for confusion. E.g. the number of rows/columns for a given map will no longer match the result of Rast_{input,output}_window_{rows,cols}, which in turn affects functions such as Rast_allocate_*_buf and Rast_zero_*_buf.

[helena]

Glynn, thanks for the explanation - there is no reason for v.surf.rst to read/write multiple maps at different resolutions, I am guessing that it tries to read the mask at its original resolution, but it should read it at the same resolution as is set by g.region for output. I can try to look at the code when time allows,

Helena

[bhlevca]

I am curios what is the status of this. Basically the current code does not allow to interpolate with a mask because of the ERROR: Input window changed while maps are open for read.

I don't understand Glynn's explanation about multiple resolutions. I have one mask and a vector layer and I want to generate a raster based on the values of the vector layer and to be limited by the boundaries provided by the mask, which in my opinion could be a vector Polygon as well.

I would say this is a pretty major feature.

[neteler]

More generally:

Also in my opinion it would be good to have a mask= parameter available for the other interpolation modules.

See also the comment about v.surf.icw here: ​http://gis.stackexchange.com/questions/85387/how-to-interpolate-a-point-vector-layer-inside-a-boundary-polygon-layer-qgis-g

[cmbarton]

Replying to neteler:

More generally:

Also in my opinion it would be good to have a mask= parameter available for the other interpolation modules.

I definitely agree!

See also the comment about v.surf.icw here: ​http://gis.stackexchange.com/questions/85387/how-to-interpolate-a-point-vector-layer-inside-a-boundary-polygon-layer-qgis-g

[bhlevca]

Replying to cmbarton:

Replying to neteler:

Also in my opinion it would be good to have a mask= parameter available for the other interpolation modules.

I definitely agree!

You would be surprised to know that in ArcGis only two interpolation modes (geoprocessing IDW and Kernel) can use a mask layer as a barrier. If GRASS achieve what Markus suggested it will have the upper hand in this domain.

[neteler]

See also trac #1718 for "ERROR: Input window changed while maps are open for read".

See also r47601 for related fixed in i.topo.corr: Essentially replace

G_get_set_window(&window); --> Rast_get_window(&window);

From ML: Rast_set_window() should be called before the map is opened.

[mmetz]

Fixed in trunk r59737 and relbr70 r59738.

[neteler]

Checked using the example of the manual page as well as this command:

v.surf.rst elevrand zcol=value elev=elev_fields_masked mask=fields 

The mask map is now respected. Closing as fixed.

[cmbarton]

This will automatically fix this same problem in r.fill.nulls too, right?

[neteler]

Replying to cmbarton:

This will automatically fix this same problem in r.fill.nulls too, right?

Seems yes, the manual page example of r.fillnulls runs fine now.

[annakrat]

Additional fixes done r62131. It allows to use NULLs (not only zeros) in the mask and it fixes using standard MASK.