[GRASS-user] sample a strds at specific locations (areas)

Margherita Di Leo diregola at gmail.com
Wed Sep 4 03:01:37 PDT 2019


Hi,

I have an update on this. It seems to work if the vector file resides in
the same mapset! I g.copied big_areas2 at PERMANENT to current mapset and the
excluded polygons are much less. It surely needs further investigation...

Kind regards,

On Fri, Aug 9, 2019 at 4:50 PM Margherita Di Leo <diregola at gmail.com> wrote:

> Hi,
>
> excuse me if I return on this. I have again the same problem, of v.strds.stats
> skipping lots of polygons. Now I'm using a brand new dataset - strds of
> Sentinel 1 and a brand now vector of polygons - and the skipped polygons
> are not narrow, I'm sure that there are cells centroids in it. See for
> example screenshot attached, depicting a polygon that was skipped. I also
> tried to run v.strds.stats on that polygon alone, like:
>
> v.strds.stats in=big_areas2 at PERMANENT where="cat == '1'"
> strds=db_cross_pol at scaled out=test method=average
>
> v.info test
>
>  +----------------------------------------------------------------------------+
>  | Name:            test
>    |
>  | Mapset:          stats
>     |
>  | Location:        S1
>    |
>  | Database:        /media/madi/TOSHIBA EXT/S1/grassdata
>    |
>  | Title:           Output from v.patch
>     |
>  | Map scale:       1:1
>     |
>  | Name of creator: madi
>    |
>  | Organization:
>    |
>  | Source date:     Fri Aug  9 12:12:49 2019
>    |
>  | Timestamp (first layer): none
>    |
>
>  |----------------------------------------------------------------------------|
>  | Map format:      native
>    |
>
>  |----------------------------------------------------------------------------|
>  |   Type of map: vector (level: 2)
>     |
>  |
>    |
>  |   Number of points:       0               Number of centroids:  0
>    |
>  |   Number of lines:        0               Number of boundaries: 0
>    |
>  |   Number of areas:        0               Number of islands:    0
>    |
>  |
>    |
>  |   Map is 3D:              No
>     |
>  |   Number of dblinks:      0
>    |
>  |
>    |
>  |   Projection: UTM (zone 34)
>    |
>  |
>    |
>  |               N:                 0    S:                 0
>     |
>  |               E:                 0    W:                 0
>     |
>  |
>    |
>  |   Digitization threshold: 0
>    |
>  |   Comment:
>     |
>  |
>    |
>
>  +----------------------------------------------------------------------------+
>
> I confess I hadn't look much further into it because I thought it was a
> problem with nodata, but this is not the case and I think it's  worth of
> further investigation. Thanks for any pointers.
>
> Regards,
>
>
> On Thu, Feb 7, 2019 at 4:10 PM Margherita Di Leo <diregola at gmail.com>
> wrote:
>
>> Thanks for testing, Vero. I assume it's due to a local problem then.
>>
>> Cheers,
>>
>> On Thu, Feb 7, 2019 at 3:29 PM Veronica Andreo <veroandreo at gmail.com>
>> wrote:
>>
>>> Hi Stefan and Madi,
>>>
>>> Thanks Stefan for the explanations :) Indeed I agree that a flag to
>>> avoid the alignment to input rasters (and just use region settings) sounds
>>> good.
>>> I tested what Madi said, but cannot reproduce in the climate NC location
>>> [0]. This is the command I used:
>>>
>>> v.strds.stats in=boundary_county where="cat == '261'" strds=tempmean
>>> t_where="start_time >= '2012-01-01'" out=test
>>>
>>> to make it faster (it feels indeed kinda slow for the whole vector and
>>> full time series), I selected only one polygon and a range of dates. I get
>>> the table as expected while leaving methods by default.
>>>
>>> best,
>>> Vero
>>>
>>> [0]
>>> http://courses.ncsu.edu/mea592/common/media/02/nc_climate_spm_2000_2012.zip
>>>
>>>
>>>
>>>
>>> El jue., 7 feb. 2019 a las 14:39, Margherita Di Leo (<diregola at gmail.com>)
>>> escribió:
>>>
>>>> Thank you Stefan for your help! I figured what happens. v.strds.stats
>>>> with default method produces in my case a corrupted output, topology is
>>>> there but there's no table associated to it. If I specify method=average, I
>>>> do obtain the table, and the mystery is solved: some of the polygons fall
>>>> into a nodata (due to cloud mask). If anyone else can reproduce the
>>>> corrupted table issue, I can file a ticket for that.
>>>>
>>>> Thanks for help!
>>>> Regards,
>>>>
>>>> On Thu, Feb 7, 2019 at 12:33 PM Stefan Blumentrath <
>>>> Stefan.Blumentrath at nina.no> wrote:
>>>>
>>>>> Hi Madi,
>>>>>
>>>>>
>>>>>
>>>>> With this combination (polygon size vs. raster resolution), the shape
>>>>> of the polygons can be an issue (narrow areas that do not cover the center
>>>>> of any pixel).
>>>>>
>>>>>
>>>>>
>>>>> Debugging should be simple with v.db.select or v.extract.
>>>>>
>>>>>
>>>>>
>>>>> Areas that did not get rasterized should be NULL in the column with
>>>>> statistics computed with v.rast.stats.
>>>>>
>>>>>
>>>>>
>>>>> In verbose mode v.rast.stats (or probably even v.to.rast) should
>>>>> probably give a more informative Warning message (e.g. listing categories
>>>>> not rasterized). It also would help if you can rasterize the areas yourself
>>>>> and provide a raster with categories as (optional) input to v.rast.stats…
>>>>>
>>>>>
>>>>>
>>>>> For high resolution data like yours, the speed improvement of multiple
>>>>> raster input might help quite a bit esp. with many maps in the time series.
>>>>> Will see if I can come up with a patch rather soon…
>>>>>
>>>>>
>>>>>
>>>>> Cheers
>>>>>
>>>>> Stefan
>>>>>
>>>>>
>>>>>
>>>>> *From:* Margherita Di Leo <diregola at gmail.com>
>>>>> *Sent:* torsdag 7. februar 2019 10:37
>>>>> *To:* Stefan Blumentrath <Stefan.Blumentrath at nina.no>
>>>>> *Cc:* Veronica Andreo <veroandreo at gmail.com>; grass-user <
>>>>> grass-user at lists.osgeo.org>
>>>>> *Subject:* Re: [GRASS-user] sample a strds at specific locations
>>>>> (areas)
>>>>>
>>>>>
>>>>>
>>>>> Hi,
>>>>>
>>>>>
>>>>>
>>>>> thank you for your replies. To give a little more context: I selected
>>>>> my polygon areas to be > 0.5 ha each (this would be 50000 mq if I'm not
>>>>> mistaken) and I'm sampling NDVI maps at 10m resolution (the region being
>>>>> the same as NDVI maps). So I think I need an idea on how to debug the areas
>>>>> that were excluded to check them individually to see what could be the
>>>>> problem...
>>>>>
>>>>> Regarding the alignment problem, if I understand it correctly: if the
>>>>> polygon doesn't include the *center* of the raster beneath it, can't
>>>>> retrieve the value and the polygon is discarded? But a value exists, so it
>>>>> would be correct that it returned a value in any case. But I admit I don't
>>>>> have a full grasp of the problem.
>>>>>
>>>>>
>>>>>
>>>>> On Wed, Feb 6, 2019 at 10:52 PM Stefan Blumentrath <
>>>>> Stefan.Blumentrath at nina.no> wrote:
>>>>>
>>>>> Hi Vero,
>>>>>
>>>>>
>>>>>
>>>>> I think there is a little misunderstanding.
>>>>>
>>>>> v.rast.stats did not change it behaviour with regards to the
>>>>> computational region (at least not if only one raster map is used). The
>>>>> alignment to the input raster (resolution) has been around since the module
>>>>> got ported to Python (like 10 years ago):
>>>>>
>>>>>
>>>>> https://trac.osgeo.org/grass/browser/grass/trunk/scripts/v.rast.stats/v.rast.stats.py?rev=33522#L148
>>>>>
>>>>>
>>>>>
>>>>> So, adding a flag for skipping the alignment was more an idea for an
>>>>> enhancement that allows the behaviour you seem to prefer (too).
>>>>>
>>>>>
>>>>>
>>>>> Cheers
>>>>>
>>>>> Stefan
>>>>>
>>>>>
>>>>>
>>>>> *From:* Veronica Andreo <veroandreo at gmail.com>
>>>>> *Sent:* onsdag 6. februar 2019 21:38
>>>>> *To:* Stefan Blumentrath <Stefan.Blumentrath at nina.no>
>>>>> *Cc:* Margherita Di Leo <diregola at gmail.com>; grass-user <
>>>>> grass-user at lists.osgeo.org>
>>>>> *Subject:* Re: [GRASS-user] sample a strds at specific locations
>>>>> (areas)
>>>>>
>>>>>
>>>>>
>>>>> I had a similar problem some time ago, just that it was not raster
>>>>> resolution, but region resolution that I changed to solve my problem (see
>>>>> this thread and MM's answer:
>>>>> http://osgeo-org.1560.x6.nabble.com/v-to-rast-for-polygons-not-overlapping-center-of-raster-cell-td5355686.html#a5355729
>>>>> )
>>>>>
>>>>>
>>>>>
>>>>> IIUC, MM's proposed solution to my case then does not work anymore
>>>>> because v.to.rast call inside v.rast.stats is affected by the region
>>>>> alignment to the raster to be queried. So, the solution is indeed now, to
>>>>> change raster resolution... ? Then the region would be aligned to it (them)?
>>>>>
>>>>>
>>>>>
>>>>> If one has large areas or long time series and has to resample all
>>>>> rasters to get smallish polygons rasterized, I do not see the advantage of
>>>>> this new behavior... but maybe I'm missing something
>>>>>
>>>>>
>>>>>
>>>>> Cheers,
>>>>>
>>>>> Vero
>>>>>
>>>>>
>>>>>
>>>>> El mié., 6 feb. 2019 16:54, Stefan Blumentrath <
>>>>> Stefan.Blumentrath at nina.no> escribió:
>>>>>
>>>>> Ciao Madi, Vero,
>>>>>
>>>>>
>>>>>
>>>>> Starting with GRASS 7.6, also centroids are used to get the raster
>>>>> representation of your area vector map. That increases the likelihood of
>>>>> smaller areas to be rasterized.
>>>>>
>>>>> Increasing the resolution of the current region alone does not help,
>>>>> because v.rast.stats temporarily changes the computational region to align
>>>>> with the input raster map(s) (see also:
>>>>> https://trac.osgeo.org/grass/ticket/3523 and
>>>>> https://trac.osgeo.org/grass/ticket/3598 for discussion) Even if the
>>>>> first ticket is closed, comments are welcome.
>>>>>
>>>>> It might make sense to add a flag to v.rast.stats like in
>>>>> r.slope.aspect to not align the computational region.
>>>>>
>>>>>
>>>>>
>>>>> Furthermore, with regards to efficiency, v.strds.stats could gain some
>>>>> speed if multi-raster support in v.rast.stats - added in G 7.6 - would be
>>>>> handed down to the addon. Might almost double the speed for larger STRDS…
>>>>>
>>>>>
>>>>>
>>>>> Cheers
>>>>>
>>>>> Stefan
>>>>>
>>>>>
>>>>>
>>>>> *From:* grass-user <grass-user-bounces at lists.osgeo.org> *On Behalf Of
>>>>> *Veronica Andreo
>>>>> *Sent:* onsdag 6. februar 2019 17:20
>>>>> *To:* Margherita Di Leo <diregola at gmail.com>
>>>>> *Cc:* GRASS user list <grass-user at lists.osgeo.org>
>>>>> *Subject:* Re: [GRASS-user] sample a strds at specific locations
>>>>> (areas)
>>>>>
>>>>>
>>>>>
>>>>> Hi Madi
>>>>>
>>>>>
>>>>>
>>>>> El mié., 6 feb. 2019 a las 16:31, Margherita Di Leo (<
>>>>> diregola at gmail.com>) escribió:
>>>>>
>>>>> I have a question regarding v.strds.stats. I get the following warning
>>>>> message:
>>>>>
>>>>>
>>>>>
>>>>> WARNING: Not all vector categories converted to raster. Converted 120
>>>>> of 265.
>>>>>
>>>>>
>>>>>
>>>>> What could be the reason for that?
>>>>>
>>>>>
>>>>>
>>>>> Some vector areas might not be converted because they are too small
>>>>> with respect to the pixel size that you try to query. Others will tell
>>>>> better but I think the polygon must overlap the center of the pixel in
>>>>> order to be converted into raster. One solution could be to resample your
>>>>> rasters to a higher resolution.
>>>>>
>>>>> HTH,
>>>>>
>>>>> Vero
>>>>>
>>>>>
>>>>>
>>>>> --
>>>>>
>>>>> Margherita Di Leo
>>>>>
>>>>
>>>>
>>>> --
>>>> Margherita Di Leo
>>>>
>>>
>>
>> --
>> Margherita Di Leo
>>
>
>
> --
> Margherita Di Leo
>


-- 
Margherita Di Leo
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/grass-user/attachments/20190904/62a64783/attachment-0001.html>


More information about the grass-user mailing list