Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

2014-01-04 Thread Blumentrath, Stefan
Actually, a second vrt does the trick:

gdalbuild vrt dem_wcs_bb.vrt dem_wcs.vrt -te 381464.564874 7282848.63791 
422460.634623 7310018.58218

This worked for me.
I shall add the solution to the Wiki.

Cheers
Stefan

-Original Message-
From: martin.zbin...@gmail.com [mailto:martin.zbin...@gmail.com] On Behalf Of 
Martin Zbinden
Sent: 4. januar 2014 23:32
To: Blumentrath, Stefan; nete...@osgeo.org
Cc: grass-user@lists.osgeo.org
Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

Hi Stefan, hi Markus,
I just started to think over everything again after Stefans first maill... 
thanks for the fast reaction. I just tried the proposed HEIGHT and WIDTH 
parameter with gdal_translate without success.

I started to observe the HTTP-transfer with mitmproxy [1]. I can see, that BBOX 
parameter is just ignored. Also i can see, that r.in.gdal by itself seems to 
know about tile-size-limitations and starts to download small 2MB-tiles, while 
gdal_translate just tries to get it all in one file. If you are interested in 
more details, I can share these observations on the new WCS-wikipage?

Concerning the "limit import to current region flag" I agree absolutely with 
you two, I would just take it with a handkiss. Would this -r flag be something 
easy and quick to implement? I see that r.in.gdal is a compiled binary, so 
nothing I believe I can do myself.

Thanks to for the python-based-approach, which hopefully will respect tile-size 
limits. I'll give it a try, until the -r flag will be reality.

Cheers,
Martin



[1] http://mitmproxy.org/
--
Martin Zbinden
Riedacker 523
3154 Rüschegg Heubach
+41 78 628 28 82


2014/1/4 Blumentrath, Stefan :
> Hi Martin,
> Since you mentioned that you are familiar with python, here is a python-based 
> approach:
> http://gis.stackexchange.com/questions/41123/how-to-extract-a-wcs-laye
> r-with-gdal-translate
>
> Cheers
> Stefan
>
> -Original Message-
> From: grass-user-boun...@lists.osgeo.org 
> [mailto:grass-user-boun...@lists.osgeo.org] On Behalf Of Blumentrath, 
> Stefan
> Sent: 4. januar 2014 22:17
> To: Martin Zbinden; grass-user@lists.osgeo.org
> Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region 
> extent only
>
> Hi Martin,
>
> Just add a line like this (adjusted to your case) after the CoverageName line 
> in the .vrt.
> &SERVICE=WCS&VERSION=1.0.0&REQUEST=GetCoverage&FORMA
> T=GeoTIFF&COVERAGE=layername&BBOX=381464.564874,7282848.63791,7310018.
> 58218,422460.634623&CRS=EPSG:32633&WIDTH=800&HEIGHT=600 tra>
>
> Likely you do not need all these GetCoverage specifications.
>
> However, with the BBOX set to your region you can import the .vrt directly to 
> GRASS with r.in.gdal.
>
> (See http://www.gdal.org/frmt_wcs.html or 
> http://wiki.esipfed.org/images/2/2c/BestPractices_WCS.pdf page 24)
>
> Cheers
> Stefan
>
>
> -Original Message-
> From: martin.zbin...@gmail.com [mailto:martin.zbin...@gmail.com] On 
> Behalf Of Martin Zbinden
> Sent: 4. januar 2014 21:01
> To: Blumentrath, Stefan
> Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region 
> extent only
>
> Hey Stefan,
> Thanks, you helped me already a lot. I'm looking forward to build some 
> python script which does these steps automatically. It seems to me 
> that these .vrt-files can do kind of magic :-)
>
> However, one concern is still left. When I choose a too big projwin, then 
> Mapserver gives an error.
>
> gdal_translate -co "COMPRESS=LZW" -projwin 2592000 1182500 2597500
> 1179000  dem_wcs.vrt dem_wcs.tif
> Computed -srcwin 53295 56717 2750 1750 from projected window.
>
> 0ERROR 1: msWCSGetCoverage(): WCS server error. Raster size out of range, 
> width and height of resulting coverage must be no more than MAXSIZE=2048.
>
> It seems that r.in.gdal ist conscious about this and downloads data tile by 
> tile. Maybe you can give me a hint me how to do the same with my own script. 
> No way to tell the projwin through the GDAL-VRT directly, so I could load it 
> directly with r.in.gdal?
>
> Cheers,
> Martin
> --
> Martin Zbinden
> Riedacker 523
> 3154 Rüschegg Heubach
> +41 78 628 28 82
>
>
> 2014/1/3 Blumentrath, Stefan :
>> Hei Martin,
>>
>> You could try making a (temporary) virtual GDAL-VRT dataset (it’s a simple 
>> textfile (XML) containing information like this, in the example below named 
>> "dem_wcs.vrt"):
>>
>> 
>>  http://your_wcs_server/wcs.dtm?
>>  layername
>> 
>>
>> For more info on vrt-format see:
>> http://www.gdal.org/gdal_vrttut.html
>>
>> After that you could use gdal_translate 
>> (http://www.gdal.org/gdal_translate.html), which has a projwin parameter 
>> where you can feed you GRASS region boundaries from g.region -up:
>> gdal_translate -co "COMPRESS=LZW" -projwin ulx uly lrx lry 
>> dem_wcs.vrt dem_wcs.tif
>>
>> Hope that helps,
>>
>> Cheers
>> Stefan
>>
>> -Original Message-
>> From: grass-user-boun...@lists.osgeo.org
>> [mailto:grass-user-boun...@lists.osgeo.org] On Behalf Of Martin 
>> Zbinden
>> Sent: 3. 

Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

2014-01-04 Thread Martin Zbinden
Hi Stefan, hi Markus,
I just started to think over everything again after Stefans first
maill... thanks for the fast reaction. I just tried the proposed
HEIGHT and WIDTH parameter with gdal_translate without success.

I started to observe the HTTP-transfer with mitmproxy [1]. I can see,
that BBOX parameter is just ignored. Also i can see, that r.in.gdal by
itself seems to know about tile-size-limitations and starts to
download small 2MB-tiles, while gdal_translate just tries to get it
all in one file. If you are interested in more details, I can share
these observations on the new WCS-wikipage?

Concerning the "limit import to current region flag" I agree
absolutely with you two, I would just take it with a handkiss. Would
this -r flag be something easy and quick to implement? I see that
r.in.gdal is a compiled binary, so nothing I believe I can do myself.

Thanks to for the python-based-approach, which hopefully will respect
tile-size limits. I'll give it a try, until the -r flag will be
reality.

Cheers,
Martin



[1] http://mitmproxy.org/
--
Martin Zbinden
Riedacker 523
3154 Rüschegg Heubach
+41 78 628 28 82


2014/1/4 Blumentrath, Stefan :
> Hi Martin,
> Since you mentioned that you are familiar with python, here is a python-based 
> approach:
> http://gis.stackexchange.com/questions/41123/how-to-extract-a-wcs-layer-with-gdal-translate
>
> Cheers
> Stefan
>
> -Original Message-
> From: grass-user-boun...@lists.osgeo.org 
> [mailto:grass-user-boun...@lists.osgeo.org] On Behalf Of Blumentrath, Stefan
> Sent: 4. januar 2014 22:17
> To: Martin Zbinden; grass-user@lists.osgeo.org
> Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region extent 
> only
>
> Hi Martin,
>
> Just add a line like this (adjusted to your case) after the CoverageName line 
> in the .vrt.
> &SERVICE=WCS&VERSION=1.0.0&REQUEST=GetCoverage&FORMAT=GeoTIFF&COVERAGE=layername&BBOX=381464.564874,7282848.63791,7310018.58218,422460.634623&CRS=EPSG:32633&WIDTH=800&HEIGHT=600
>
> Likely you do not need all these GetCoverage specifications.
>
> However, with the BBOX set to your region you can import the .vrt directly to 
> GRASS with r.in.gdal.
>
> (See http://www.gdal.org/frmt_wcs.html or 
> http://wiki.esipfed.org/images/2/2c/BestPractices_WCS.pdf page 24)
>
> Cheers
> Stefan
>
>
> -Original Message-
> From: martin.zbin...@gmail.com [mailto:martin.zbin...@gmail.com] On Behalf Of 
> Martin Zbinden
> Sent: 4. januar 2014 21:01
> To: Blumentrath, Stefan
> Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region extent 
> only
>
> Hey Stefan,
> Thanks, you helped me already a lot. I'm looking forward to build some python 
> script which does these steps automatically. It seems to me that these 
> .vrt-files can do kind of magic :-)
>
> However, one concern is still left. When I choose a too big projwin, then 
> Mapserver gives an error.
>
> gdal_translate -co "COMPRESS=LZW" -projwin 2592000 1182500 2597500
> 1179000  dem_wcs.vrt dem_wcs.tif
> Computed -srcwin 53295 56717 2750 1750 from projected window.
>
> 0ERROR 1: msWCSGetCoverage(): WCS server error. Raster size out of range, 
> width and height of resulting coverage must be no more than MAXSIZE=2048.
>
> It seems that r.in.gdal ist conscious about this and downloads data tile by 
> tile. Maybe you can give me a hint me how to do the same with my own script. 
> No way to tell the projwin through the GDAL-VRT directly, so I could load it 
> directly with r.in.gdal?
>
> Cheers,
> Martin
> --
> Martin Zbinden
> Riedacker 523
> 3154 Rüschegg Heubach
> +41 78 628 28 82
>
>
> 2014/1/3 Blumentrath, Stefan :
>> Hei Martin,
>>
>> You could try making a (temporary) virtual GDAL-VRT dataset (it’s a simple 
>> textfile (XML) containing information like this, in the example below named 
>> "dem_wcs.vrt"):
>>
>> 
>>  http://your_wcs_server/wcs.dtm?
>>  layername
>> 
>>
>> For more info on vrt-format see:
>> http://www.gdal.org/gdal_vrttut.html
>>
>> After that you could use gdal_translate 
>> (http://www.gdal.org/gdal_translate.html), which has a projwin parameter 
>> where you can feed you GRASS region boundaries from g.region -up:
>> gdal_translate -co "COMPRESS=LZW" -projwin ulx uly lrx lry dem_wcs.vrt
>> dem_wcs.tif
>>
>> Hope that helps,
>>
>> Cheers
>> Stefan
>>
>> -Original Message-
>> From: grass-user-boun...@lists.osgeo.org
>> [mailto:grass-user-boun...@lists.osgeo.org] On Behalf Of Martin
>> Zbinden
>> Sent: 3. januar 2014 21:46
>> To: grass-user@lists.osgeo.org
>> Subject: [GRASS-user] WCS import into GRASS GIS - select region extent
>> only
>>
>> Hi,
>>
>> For my Bachelor-Thesis I want to optimize the erosion risk map of 
>> Switzerland for use in agricultural extension. Although I concentrate on 3 
>> cantons only (BE, FR, SO),  I have to manage huge amounts of data. The fact 
>> that the raster data I got from the federal office are overlapping each 
>> other (following the defined fieldblocks) doesn't make things easier, but 
>> I've at least fo

Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

2014-01-04 Thread Blumentrath, Stefan
Hi Martin,
Since you mentioned that you are familiar with python, here is a python-based 
approach:
http://gis.stackexchange.com/questions/41123/how-to-extract-a-wcs-layer-with-gdal-translate

Cheers
Stefan

-Original Message-
From: grass-user-boun...@lists.osgeo.org 
[mailto:grass-user-boun...@lists.osgeo.org] On Behalf Of Blumentrath, Stefan
Sent: 4. januar 2014 22:17
To: Martin Zbinden; grass-user@lists.osgeo.org
Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

Hi Martin,

Just add a line like this (adjusted to your case) after the CoverageName line 
in the .vrt.
&SERVICE=WCS&VERSION=1.0.0&REQUEST=GetCoverage&FORMAT=GeoTIFF&COVERAGE=layername&BBOX=381464.564874,7282848.63791,7310018.58218,422460.634623&CRS=EPSG:32633&WIDTH=800&HEIGHT=600

Likely you do not need all these GetCoverage specifications.

However, with the BBOX set to your region you can import the .vrt directly to 
GRASS with r.in.gdal.

(See http://www.gdal.org/frmt_wcs.html or 
http://wiki.esipfed.org/images/2/2c/BestPractices_WCS.pdf page 24)

Cheers
Stefan


-Original Message-
From: martin.zbin...@gmail.com [mailto:martin.zbin...@gmail.com] On Behalf Of 
Martin Zbinden
Sent: 4. januar 2014 21:01
To: Blumentrath, Stefan
Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

Hey Stefan,
Thanks, you helped me already a lot. I'm looking forward to build some python 
script which does these steps automatically. It seems to me that these 
.vrt-files can do kind of magic :-)

However, one concern is still left. When I choose a too big projwin, then 
Mapserver gives an error.

gdal_translate -co "COMPRESS=LZW" -projwin 2592000 1182500 2597500
1179000  dem_wcs.vrt dem_wcs.tif
Computed -srcwin 53295 56717 2750 1750 from projected window.

0ERROR 1: msWCSGetCoverage(): WCS server error. Raster size out of range, width 
and height of resulting coverage must be no more than MAXSIZE=2048.

It seems that r.in.gdal ist conscious about this and downloads data tile by 
tile. Maybe you can give me a hint me how to do the same with my own script. No 
way to tell the projwin through the GDAL-VRT directly, so I could load it 
directly with r.in.gdal?

Cheers,
Martin
--
Martin Zbinden
Riedacker 523
3154 Rüschegg Heubach
+41 78 628 28 82


2014/1/3 Blumentrath, Stefan :
> Hei Martin,
>
> You could try making a (temporary) virtual GDAL-VRT dataset (it’s a simple 
> textfile (XML) containing information like this, in the example below named 
> "dem_wcs.vrt"):
>
> 
>  http://your_wcs_server/wcs.dtm?
>  layername
> 
>
> For more info on vrt-format see:
> http://www.gdal.org/gdal_vrttut.html
>
> After that you could use gdal_translate 
> (http://www.gdal.org/gdal_translate.html), which has a projwin parameter 
> where you can feed you GRASS region boundaries from g.region -up:
> gdal_translate -co "COMPRESS=LZW" -projwin ulx uly lrx lry dem_wcs.vrt 
> dem_wcs.tif
>
> Hope that helps,
>
> Cheers
> Stefan
>
> -Original Message-
> From: grass-user-boun...@lists.osgeo.org
> [mailto:grass-user-boun...@lists.osgeo.org] On Behalf Of Martin 
> Zbinden
> Sent: 3. januar 2014 21:46
> To: grass-user@lists.osgeo.org
> Subject: [GRASS-user] WCS import into GRASS GIS - select region extent 
> only
>
> Hi,
>
> For my Bachelor-Thesis I want to optimize the erosion risk map of Switzerland 
> for use in agricultural extension. Although I concentrate on 3 cantons only 
> (BE, FR, SO),  I have to manage huge amounts of data. The fact that the 
> raster data I got from the federal office are overlapping each other 
> (following the defined fieldblocks) doesn't make things easier, but I've at 
> least found r.patch which can combine them one by one.
>
> I thought, with Mapserver  (V 6.4.0) I'd found the perfect solution to all my 
> problems. Now I'm stuck after 5+ hours of trying to import the data from my 
> local WCS-Service. I have seen the hints on GRASS Wiki [1] and r.in.gdal 
> starts to download geotiff-tiles from Mapserver. The problem is, that I don't 
> want all the 12 GB of data, but only the grass region extent.
>
> The following http-query does what I expect, it makes me a GeoTIFF-File from 
> my region defined by BBOX=... .
> http://localhost/mywcs?SERVICE=WCS&VERSION=1.0.0&REQUEST=GetCoverage&c
> overage=blw.eros_z&FORMAT=image/tiff&CRS=EPSG:2056&BBOX=2596000,118100
> 0,2597500,1182500&RESX=2&RESY=2
>
>
> I tried to restrict the download to the  region by defining a BBOX in this 
> WCS_GDAL xml file. However, this didn't work and isn't supposed to work [2]. 
> Is there a way to say r.in.gdal to only get the data from the region extent?
>
> I know there is this new, very clever r.in.wms2 script, which does exactly 
> what I want -- for WMS. Is there a way to get WMS-Mapserver to serve RAW-data 
> (erosion in t per annum, elevation data)?  Or maybe better create a python 
> script to download map by a custom http-request?
>
> Thanks in advance for any advice! I hope I can get it to work someho

Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

2014-01-04 Thread Markus Neteler
Hi Stefan,

On Sat, Jan 4, 2014 at 10:40 PM, Blumentrath, Stefan
 wrote:
> Hi Markus,
>
> I should not write emails offline ;-), so sorry for reposting parts of your 
> answer...
> And thanks for the hints. I was not aware that the BBOX had no effect in the 
> XML.

... this is what I understood from FrankW's email.

> Hopefully, the HEIGHT and WIDTH parameter work, as Martin seemed to have 
> problems to fetch the data for his region due to tile-size limits...

Please add your findings in this new Wiki page:

http://grasswiki.osgeo.org/wiki/WCS

> If you consider enhancing r.external (or r.in.gdal) in the direction of a 
> projwin option, maybe a "limit import to current region flag" would be more 
> easy to use and maybe also more consistent (like e.g in v.in.ogr).

Indeed, this is a long-term wish for me :-)
The code in
http://svn.osgeo.org/gdal/trunk/gdal/apps/gdal_translate.cpp
should give some inspiration...

Best
Markus
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] GRASS Addons manual pages online

2014-01-04 Thread Markus Neteler
On Tue, Dec 10, 2013 at 1:10 PM, Martin Landa  wrote:
> Hi all,
>
> manual pages of GRASS Addons modules are now online [1,2,3]. There are
> few remaining issues. Footer links (index, full index, ...) are
> broken. No index.html page containing list of modules.

At least for GRASS 7 Addons this has been now solved.

I have made it a "News" item on the Web site:
http://grass.osgeo.org/news/29/15/GRASS-GIS-7-Addons-Manuals-online/

Markus


> Other todo is related to wxGUI Extension manager which could display
> these online manual pages to give the user better idea what selected
> addon module is doing.
>
> Any other ideas welcome of course. Martin
>
> [1] http://grass.osgeo.org/grass70/manuals/addons/
> [2] http://grass.osgeo.org/grass65/manuals/addons/
> [3] http://grass.osgeo.org/grass64/manuals/addons/
>
> --
> Martin Landa  * http://geo.fsv.cvut.cz/~landa
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

2014-01-04 Thread Blumentrath, Stefan
Hi Markus,

I should not write emails offline ;-), so sorry for reposting parts of your 
answer...
And thanks for the hints. I was not aware that the BBOX had no effect in the 
XML.

Hopefully, the HEIGHT and WIDTH parameter work, as Martin seemed to have 
problems to fetch the data for his region due to tile-size limits...

If you consider enhancing r.external (or r.in.gdal) in the direction of a 
projwin option, maybe a "limit import to current region flag" would be more 
easy to use and maybe also more consistent (like e.g in v.in.ogr).

Cheers
Stefan

-Original Message-
From: neteler.os...@gmail.com [mailto:neteler.os...@gmail.com] On Behalf Of 
Markus Neteler
Sent: 4. januar 2014 21:49
To: Blumentrath, Stefan
Cc: grass-user@lists.osgeo.org
Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

Hi,

just for fun I have made a test with r.external (of GRASS GIS 7, it potentially 
works with G6 as well):

# create ASCII file with this content, named e.g. "wcs_geoserver_LL.wcs":

  http://demo.opengeo.org/geoserver/wcs?
  Img_Sample


or download this file from
http://svn.osgeo.org/gdal/trunk/autotest/gdrivers/data/geoserver.wcs

Start GRASS in latlong location, then

GRASS 7.0.svn (latlong):~ > r.external input=wcs_geoserver_LL.wcs 
output=wcs_data --o Projection of input dataset and current location appear to 
match Reading band 1 of 3...
r.external complete. Link to raster map  created.
Reading band 2 of 3...
r.external complete. Link to raster map  created.
Reading band 3 of 3...
r.external complete. Link to raster map  created.
Imagery group  created

Interestingly the content of wcs_geoserver_LL.wcs is expanded on the fly, i.e. 
the file is rewritten.
BTW I had to run r.external twice on the same file to avoid some GDAL errors 
(perhaps the demo server is simply stressed).

GRASS 7.0.svn (latlong):~ > g.region rast=wcs_data.1 -p
projection: 3 (Latitude-Longitude)
zone:   0
datum:  wgs84
ellipsoid:  wgs84
north:  54:06:50.76N
south:  20:42:18.72N
west:   130:51:06.048W
east:   62:00:19.44W
nsres:  0:03:21.123813
ewres:  0:04:12.132867
rows:   598
cols:   983
cells:  587834

Voilà, the WCS layer is loaded (then I used the RGB viewer in the wxGUI).

Eventually I then played around with the bounding box
(BBOX=minx,miny,maxx,maxy) and "injected" it in above small file 
&bbox=-103.93,33.74-93.4940.61

so that it becomes:


  http://demo.opengeo.org/geoserver/wcs?
  Img_Sample
  &bbox=-103.93,33.74-93.4940.61


BUT, as stated in
http://lists.osgeo.org/pipermail/gdal-dev/2008-March/016552.html
it does not have any effect.
I don't know if it would make sense to implement a "projwin" parameter in 
r.external.

BTW: Another file/WCS server to play with is 
http://svn.osgeo.org/gdal/trunk/autotest/gdrivers/data/srtmplus.wcs

Not really helpful but at least we know that r.external works with WCS servers.

Markus
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

2014-01-04 Thread Blumentrath, Stefan
Hi Martin,

Just add a line like this (adjusted to your case) after the CoverageName line 
in the .vrt.
&SERVICE=WCS&VERSION=1.0.0&REQUEST=GetCoverage&FORMAT=GeoTIFF&COVERAGE=layername&BBOX=381464.564874,7282848.63791,7310018.58218,422460.634623&CRS=EPSG:32633&WIDTH=800&HEIGHT=600

Likely you do not need all these GetCoverage specifications.

However, with the BBOX set to your region you can import the .vrt directly to 
GRASS with r.in.gdal.

(See http://www.gdal.org/frmt_wcs.html or 
http://wiki.esipfed.org/images/2/2c/BestPractices_WCS.pdf page 24)

Cheers
Stefan


-Original Message-
From: martin.zbin...@gmail.com [mailto:martin.zbin...@gmail.com] On Behalf Of 
Martin Zbinden
Sent: 4. januar 2014 21:01
To: Blumentrath, Stefan
Subject: Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

Hey Stefan,
Thanks, you helped me already a lot. I'm looking forward to build some python 
script which does these steps automatically. It seems to me that these 
.vrt-files can do kind of magic :-)

However, one concern is still left. When I choose a too big projwin, then 
Mapserver gives an error.

gdal_translate -co "COMPRESS=LZW" -projwin 2592000 1182500 2597500
1179000  dem_wcs.vrt dem_wcs.tif
Computed -srcwin 53295 56717 2750 1750 from projected window.

0ERROR 1: msWCSGetCoverage(): WCS server error. Raster size out of range, width 
and height of resulting coverage must be no more than MAXSIZE=2048.

It seems that r.in.gdal ist conscious about this and downloads data tile by 
tile. Maybe you can give me a hint me how to do the same with my own script. No 
way to tell the projwin through the GDAL-VRT directly, so I could load it 
directly with r.in.gdal?

Cheers,
Martin
--
Martin Zbinden
Riedacker 523
3154 Rüschegg Heubach
+41 78 628 28 82


2014/1/3 Blumentrath, Stefan :
> Hei Martin,
>
> You could try making a (temporary) virtual GDAL-VRT dataset (it’s a simple 
> textfile (XML) containing information like this, in the example below named 
> "dem_wcs.vrt"):
>
> 
>  http://your_wcs_server/wcs.dtm?
>  layername
> 
>
> For more info on vrt-format see:
> http://www.gdal.org/gdal_vrttut.html
>
> After that you could use gdal_translate 
> (http://www.gdal.org/gdal_translate.html), which has a projwin parameter 
> where you can feed you GRASS region boundaries from g.region -up:
> gdal_translate -co "COMPRESS=LZW" -projwin ulx uly lrx lry dem_wcs.vrt 
> dem_wcs.tif
>
> Hope that helps,
>
> Cheers
> Stefan
>
> -Original Message-
> From: grass-user-boun...@lists.osgeo.org 
> [mailto:grass-user-boun...@lists.osgeo.org] On Behalf Of Martin 
> Zbinden
> Sent: 3. januar 2014 21:46
> To: grass-user@lists.osgeo.org
> Subject: [GRASS-user] WCS import into GRASS GIS - select region extent 
> only
>
> Hi,
>
> For my Bachelor-Thesis I want to optimize the erosion risk map of Switzerland 
> for use in agricultural extension. Although I concentrate on 3 cantons only 
> (BE, FR, SO),  I have to manage huge amounts of data. The fact that the 
> raster data I got from the federal office are overlapping each other 
> (following the defined fieldblocks) doesn't make things easier, but I've at 
> least found r.patch which can combine them one by one.
>
> I thought, with Mapserver  (V 6.4.0) I'd found the perfect solution to all my 
> problems. Now I'm stuck after 5+ hours of trying to import the data from my 
> local WCS-Service. I have seen the hints on GRASS Wiki [1] and r.in.gdal 
> starts to download geotiff-tiles from Mapserver. The problem is, that I don't 
> want all the 12 GB of data, but only the grass region extent.
>
> The following http-query does what I expect, it makes me a GeoTIFF-File from 
> my region defined by BBOX=... .
> http://localhost/mywcs?SERVICE=WCS&VERSION=1.0.0&REQUEST=GetCoverage&c
> overage=blw.eros_z&FORMAT=image/tiff&CRS=EPSG:2056&BBOX=2596000,118100
> 0,2597500,1182500&RESX=2&RESY=2
>
>
> I tried to restrict the download to the  region by defining a BBOX in this 
> WCS_GDAL xml file. However, this didn't work and isn't supposed to work [2]. 
> Is there a way to say r.in.gdal to only get the data from the region extent?
>
> I know there is this new, very clever r.in.wms2 script, which does exactly 
> what I want -- for WMS. Is there a way to get WMS-Mapserver to serve RAW-data 
> (erosion in t per annum, elevation data)?  Or maybe better create a python 
> script to download map by a custom http-request?
>
> Thanks in advance for any advice! I hope I can get it to work somehow.
> Liebe Grüsse
> Martin Zbinden
>
> PS: the WCS server given in GRASS wiki seems to be down this time.
>
> [1] 
> http://grasswiki.osgeo.org/wiki/Global_datasets#OGC_WCS_-_Albedo_examp
> le [2] 
> http://lists.osgeo.org/pipermail/gdal-dev/2008-March/016552.html
>
>
> --
> Martin Zbinden
> Riedacker 523
> 3154 Rüschegg Heubach
> +41 78 628 28 82
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] WCS import into GRASS GIS - select region extent only

2014-01-04 Thread Markus Neteler
Hi,

just for fun I have made a test with r.external (of GRASS GIS 7, it
potentially works with G6 as well):

# create ASCII file with this content, named e.g. "wcs_geoserver_LL.wcs":

  http://demo.opengeo.org/geoserver/wcs?
  Img_Sample


or download this file from
http://svn.osgeo.org/gdal/trunk/autotest/gdrivers/data/geoserver.wcs

Start GRASS in latlong location, then

GRASS 7.0.svn (latlong):~ > r.external input=wcs_geoserver_LL.wcs
output=wcs_data --o
Projection of input dataset and current location appear to match
Reading band 1 of 3...
r.external complete. Link to raster map  created.
Reading band 2 of 3...
r.external complete. Link to raster map  created.
Reading band 3 of 3...
r.external complete. Link to raster map  created.
Imagery group  created

Interestingly the content of wcs_geoserver_LL.wcs is expanded on the
fly, i.e. the file is rewritten.
BTW I had to run r.external twice on the same file to avoid some GDAL
errors (perhaps the demo server is simply stressed).

GRASS 7.0.svn (latlong):~ > g.region rast=wcs_data.1 -p
projection: 3 (Latitude-Longitude)
zone:   0
datum:  wgs84
ellipsoid:  wgs84
north:  54:06:50.76N
south:  20:42:18.72N
west:   130:51:06.048W
east:   62:00:19.44W
nsres:  0:03:21.123813
ewres:  0:04:12.132867
rows:   598
cols:   983
cells:  587834

Voilà, the WCS layer is loaded (then I used the RGB viewer in the wxGUI).

Eventually I then played around with the bounding box
(BBOX=minx,miny,maxx,maxy) and "injected" it in above small file
&bbox=-103.93,33.74-93.4940.61

so that it becomes:


  http://demo.opengeo.org/geoserver/wcs?
  Img_Sample
  &bbox=-103.93,33.74-93.4940.61


BUT, as stated in
http://lists.osgeo.org/pipermail/gdal-dev/2008-March/016552.html
it does not have any effect.
I don't know if it would make sense to implement a "projwin" parameter
in r.external.

BTW: Another file/WCS server to play with is
http://svn.osgeo.org/gdal/trunk/autotest/gdrivers/data/srtmplus.wcs

Not really helpful but at least we know that r.external works with WCS servers.

Markus
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] r.contour fails to close a contour at the region's border

2014-01-04 Thread Markus Neteler
On Fri, Jan 3, 2014 at 12:48 AM, John Ciolek  wrote:
> This is a re-post.
>
> I have automatically generated raster data from which I create a contour
> using r.contour.  Sometimes the "feature" to be contoured is not completely
> contained within the specified region.  When this happens, r.contour
> generates an open contour - basically the beginning line segment point does
> not equal the ending line segment point.

Could you please send a small screenshot (to me) or put it somewhere?

> Is there a way to get r.contour to close the contour at the region's border?
>
> If not, is there an easy way to close the contour?

Would it be feasible to work with a slightly larger region?
Note that g.region supports syntax like (say, relative coordinates):

g.region w=w-1000 e=e+1000 n=n+2000 s=s-2000

Best
Markus
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] r.sunmask cast shadow question

2014-01-04 Thread Markus Neteler
On Sat, Jan 4, 2014 at 12:17 PM, Tim Michelsen
 wrote:
...
> Perhaps r.sun.hourly should get
...
> * average sunlight duration per day during a month, etc. (h/d)

This needs some modifications in the r.sun code. I was discussing this
issue some time ago with Jaro Hofierka who gave the following hints:

On Thu, Apr 16, 2009 at 18:34, Jaro Hofierka  wrote:
> Markus,
>
> sunrise, sunset times are calculated for each cell, however, not including
> shadows by surrounding terrain features. This would need another test which
> in fact is done during radiation calculations (there is a time step to sum
> up radiation), but the time is not stored.
> To minimize a number of output maps we have stored the information of
> sunrise/sunset times using min-max values in a history file (or a single
> value when a single value of latitude is used for the map). So it would need
> a modification of the code to put there extra files storing sunrise-sunset
> times. If you need a quicker solution, you would need to make tests just for
> a small part of the day (e.g. 1-2 hours around reported sunrise/sunset
> times).
>
> Jaro

Perhaps useful as a pointer.

Markus
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] r.sunmask cast shadow question

2014-01-04 Thread Tim Michelsen
Am 04.01.2014 15:12, schrieb Vaclav Petras:
> 
> 
> 
> On Sat, Jan 4, 2014 at 6:17 AM, Tim Michelsen
> mailto:timmichel...@gmx-topmail.de>> wrote:
> 
> Am 02.01.2014 23:13, schrieb Markus Neteler:
> > You may consider r.sun as well in order to compute cast shadow. I
> > played around with it recently, see
> > http://courses.neteler.org/will-the-sun-shine-on-us/
> Thanks for sharing this.
> 
> >Perhaps r.sun.hourly should get an additional flag to generate the
> binary “sun-yes” – >“sun-no” maps directly.
> Additionally, it would be nice to have a agregrate flaf
> * sum up the daily/monthly/annual radiation (kWh/m^2)
> * average sunlight duration per day during a month, etc. (h/d)
> 
> Would r.sun.hourly and r.sun.daily help to satisfy these requests?
>  
> http://grass.osgeo.org/grass70/manuals/addons/r.sun.hourly.html
> http://grass.osgeo.org/grass70/manuals/addons/r.sun.daily.html
>  
> 
> Where can such ideas be collected for GRASS addons?
> In the normal bugtracker?
> 
> 
> Never thought about that. It is an feature request, so I would put it to
> bugtracker. If the reporter thinks that it should go to addons (and not
> to main code)
I added these isses but currently have not time to work on this:
add flag for binary “sun-yes” – >“sun-no” output
http://trac.osgeo.org/grass/ticket/2154

r.sun.hourly: add flag for aggregating a series of output files
http://trac.osgeo.org/grass/ticket/2155

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user


Re: [GRASS-user] r.sunmask cast shadow question

2014-01-04 Thread Vaclav Petras
On Sat, Jan 4, 2014 at 6:17 AM, Tim Michelsen
wrote:

> Am 02.01.2014 23:13, schrieb Markus Neteler:
> > You may consider r.sun as well in order to compute cast shadow. I
> > played around with it recently, see
> > http://courses.neteler.org/will-the-sun-shine-on-us/
> Thanks for sharing this.
>
> >Perhaps r.sun.hourly should get an additional flag to generate the
> binary “sun-yes” – >“sun-no” maps directly.
> Additionally, it would be nice to have a agregrate flaf
> * sum up the daily/monthly/annual radiation (kWh/m^2)
> * average sunlight duration per day during a month, etc. (h/d)
>
> Would r.sun.hourly and r.sun.daily help to satisfy these requests?

http://grass.osgeo.org/grass70/manuals/addons/r.sun.hourly.html
http://grass.osgeo.org/grass70/manuals/addons/r.sun.daily.html


> Where can such ideas be collected for GRASS addons?
> In the normal bugtracker?
>

Never thought about that. It is an feature request, so I would put it to
bugtracker. If the reporter thinks that it should go to addons (and not to
main code), he or she can write this to ticket description (trac component
'Addons' is appropriate). However, the reporter is not obliged to know or
decide if main code or addons are appropriate.

Vaclav

>
> Kind regards,
> Timmie
>
> ___
> grass-user mailing list
> grass-user@lists.osgeo.org
> http://lists.osgeo.org/mailman/listinfo/grass-user
>
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user

Re: [GRASS-user] r.sunmask cast shadow question

2014-01-04 Thread Tim Michelsen
Am 02.01.2014 23:13, schrieb Markus Neteler:
> You may consider r.sun as well in order to compute cast shadow. I
> played around with it recently, see
> http://courses.neteler.org/will-the-sun-shine-on-us/
Thanks for sharing this.

>Perhaps r.sun.hourly should get an additional flag to generate the
binary “sun-yes” – >“sun-no” maps directly.
Additionally, it would be nice to have a agregrate flaf
* sum up the daily/monthly/annual radiation (kWh/m^2)
* average sunlight duration per day during a month, etc. (h/d)

Where can such ideas be collected for GRASS addons?
In the normal bugtracker?

Kind regards,
Timmie

___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user