Re: [GRASS-user] trouble with r.in.gdal to georeference a NetCDF file
On Tue, Nov 12, 2013 at 5:37 AM, Lee Eddington wrote: > I'm trying to import a NetCDF file from the WRF weather forecast model with > r.in.gdal, but I can't get it to georeference. Instead I end up with a x,y > coordinate system of rows and columns. The file georeferences fine in a > number of other weather graphics programs. I used some other netCDF files recently. However, I always first convert those files to GeoTIFF using gdalwarp, then I use r.in.gdal. Markus ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Openness Calculation
Hi,
On Mon, Nov 4, 2013 at 2:16 PM, Anna Petrášová wrote:
> Hi,
>
>
> On Mon, Nov 4, 2013 at 1:36 PM, Rebecca Bennett wrote:
>
>> Dear GRASS users,
>>
>> I have recently been made aware of the positive and negative openness
>> vizualisations (Yokoyama et al 2002
>> http://www.asprs.org/a/publications/pers/2002journal/march/2002_mar_257-265.pdf)
>> and would like to try to compute them in GRASS to help identify micro
>> relief.
>>
>> Though I suspect that this has already been done, I can't find a heads up
>> other than this paper (which claims to have creates the web service based
>> on GRASS here
>> http://geobrain.laits.gmu.edu/grassweb/manuals/raster/openness.html)
>>
>> http://www.lpi.usra.edu/science/stepinskiWebPage/pdfFiles/geoinformaticsConf2009.pdf
>>
>> Just a note that although similar in it's aim of highlighting
>> microtopography, this is not quite the same as the r.horizon calculation
>> (image attached for clarification).
>>
>> All the best,
>> Rebecca
>>
>>
>>
>
> I added recently a module r.skyview to grass 7 addons for computing the
> sky-view factor [1]. It uses r.horizon module. I am not sure how much
> different is what r.horizon can do and what you need but it's worth looking
> if r.horizon cannot be extended easily and then used for the openness
> computation.
>
> Anna
>
>
>
I attached a simple patch for r.horizon which enables to output also
negative elevation angles. Addon r.skyview can be then extended to compute
also openness besides the skyview factor (the difference is just the range
of elevation angles if I understood correctly). It looks like the negative
openness can be then computed in the same way, but with inverted elevation.
Regards,
Anna
> [1] https://svn.osgeo.org/grass/grass-addons/grass7/raster/r.skyview/
>
>>
>>
>>
>> ___
>> grass-user mailing list
>> grass-user@lists.osgeo.org
>> http://lists.osgeo.org/mailman/listinfo/grass-user
>>
>
>
Index: main.c
===
--- main.c (revision 58190)
+++ main.c (working copy)
@@ -108,6 +108,7 @@
int ip, jp, ip100, jp100;
int n, m, m100, n100;
int degreeOutput = FALSE;
+int negative_angles = FALSE;
float **z, **z100, **horizon_raster;
double stepx, stepy, stepxhalf, stepyhalf, stepxy, xp, yp, op, dp, xg0, xx0,
yg0, yy0, deltx, delty;
@@ -168,7 +169,7 @@
struct
{
- struct Flag *degreeOutput;
+ struct Flag *degreeOutput, *negative_angles;
}
flag;
@@ -295,6 +296,10 @@
flag.degreeOutput->description =
_("Write output in degrees (default is radians)");
+flag.negative_angles = G_define_flag();
+flag.negative_angles->key = 'n';
+flag.negative_angles->description =
+ _("Allow negative elevation angles");
if (G_parser(argc, argv))
exit(EXIT_FAILURE);
@@ -330,6 +335,7 @@
delty = fabs(cellhd.north - cellhd.south);
degreeOutput = flag.degreeOutput->answer;
+negative_angles = flag.negative_angles->answer;
elevin = parm.elevin->answer;
@@ -712,7 +718,10 @@
{
double height;
-tanh0 = 0.;
+if (negative_angles)
+tanh0 = -BIG;
+else
+tanh0 = 0.;
length = 0;
height = searching();
___
grass-user mailing list
grass-user@lists.osgeo.org
http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Classification of segments
Dear Ferdinando, [cc'ing the grass-user list. I suggest to hold such discussions there, as more people can participate] On 11/11/13 14:04, Ferdinando Urbano wrote: my name is Ferdinando Urbano. I found your email address on some threads on the web. I am writing to you because i am working on a project related to forest cover change assessment in Nepal. I want to test some possible approach that should then implemented by local technicians. I would like to propose an open source software platform (i.e. no e-cognition). I tested visual interpretation and segmentation (with i.segment) + visual classification (and verification). The results are really good. i would like to try to run an automatic (or, if not possible) supervised classification of my segments to compare the results. As far as i can see, this option is not yet available, at least as "simple" tool (e.g. feeding i.maxlik with the result of i.segment). One technique I've thought about, but not tried, yet, is the following: - i.segment - r.to.vect to create vector polygons of segments - fill the attribute table of the vector map with variables of interest - for each variable of interest, create a raster map (v.to.rast) using the values from the attribute table - each pixel of a given polygon would have the same value - group all these raster maps with i.group - run i.cluster on the group I've also posted a script [1] to show the process for classification based on i.segment and v.class.mlpy [1]. MLPY also contains unsupervised classification, so it should be possible to integrate that into v.class.mlpy. If you follow that thread you'll that Pietro Zambelli is working on a more integrated system. I see from various posts that you worked on this kind of approach, so I ask you some suggestion. i also sent an email to Eric Momsen (developer of i.segment) to check if he is planning any new feature. As said, I do not need just a way to do it, but a way that is simple for users, otherwise local technicians will not be able to manage it. i hope you have some recommendation for me, if not, thanks anyway for your attention. I would say that you are just a bit too early to find a complete easy-to-use interface. We're not far from it, and with the elements already in place it should not be hard to cook up a module with graphical interface which combines most of the elements, but it's not there, yet. Other options you might want to explore: - combine the segmentation in GRASS (or other tools) with classification modules in R. This is also something you can quite easily script and so your technicians would just have to run the script. - other free software tools exist for segmentation / classification, such as Orfeo Toolbox, SAGA, Opticks, etc - QGIS has the semi-automatic classification plugin [2] which also looks quite nice (haven't tested it myself, yet), but AFAIK, it does not offer segmentation and object-based classification. In any case, keep us posted about your experiences ! Moritz [1] http://lists.osgeo.org/pipermail/grass-user/2013-October/069189.html [2] http://fromgistors.blogspot.be/p/semi-automatic-classification-plugin.html ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Permission denied to access a mapset in Linux
Dear Glynn, When you said: GRASS won't allow you to select a mapset as the current mapset unless you are the owner of the mapset directory. How can I change the "owner of a mapset directory"? Regards, Antonio On 12-11-2013 14:46, Glynn Clements wrote: António Rocha wrote: I have a Location and Mapset in a folder in Linux Server that was created by me and I runned some scripts and created some maps etc. I have a colleague that is trying enter in the mapset and it gets Permission Denied. How can I allow her to access the mapset without any problems? You can "access" any mapset provided that you have the relevant permissions (read permission for files, read and execute permission for directories). GRASS won't allow you to select a mapset as the current mapset unless you are the owner of the mapset directory. If you just need to access maps from a particular mapset, you can add the mapset to the search path with g.mapset, or reference maps using fully-qualified names (map@mapset format). ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Permission denied to access a mapset in Linux
António Rocha wrote: > I have a Location and Mapset in a folder in Linux Server that was > created by me and I runned some scripts and created some maps etc. I > have a colleague that is trying enter in the mapset and it gets > Permission Denied. How can I allow her to access the mapset without any > problems? You can "access" any mapset provided that you have the relevant permissions (read permission for files, read and execute permission for directories). GRASS won't allow you to select a mapset as the current mapset unless you are the owner of the mapset directory. If you just need to access maps from a particular mapset, you can add the mapset to the search path with g.mapset, or reference maps using fully-qualified names (map@mapset format). -- Glynn Clements ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] trouble with r.in.gdal to georeference a NetCDF file
OK, here it is:
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
Driver: netCDF/Network Common Data Format
Files: none associated
Size is 198, 153
Coordinate System is `'
Metadata:
HGT_M#description=Topography height
HGT_M#FieldType=104
HGT_M#MemoryOrder=XY
HGT_M#sr_x=1
HGT_M#sr_y=1
HGT_M#stagger=M
HGT_M#units=meters MSL
NC_GLOBAL#BOTTOM-TOP_GRID_DIMENSION=0
NC_GLOBAL#CEN_LAT=-8.4095154
NC_GLOBAL#CEN_LON=115.02
NC_GLOBAL#corner_lats={ -10.454437, -6.3537445, -6.3537445, -10.454437,
-10.454437, -6.3537445, -6.3537445, -10.454437, -10.46785, -6.3401871,
-6.3401871, -10.46785, -10.46785, -6.3401871, -6.3401871, -10.46785 }
NC_GLOBAL#corner_lons={ 112.3332, 112.3332, 117.70679, 117.70679,
112.31956, 112.31956, 117.72044, 117.72044, 112.3332, 112.3332, 117.70679,
117.70679, 112.31956, 112.31956, 117.72044, 117.72044 }
NC_GLOBAL#DX=3000
NC_GLOBAL#DY=3000
NC_GLOBAL#DYN_OPT=2
NC_GLOBAL#FLAG_MF_XY=1
NC_GLOBAL#grid_id=3
NC_GLOBAL#GRIDTYPE=C
NC_GLOBAL#i_parent_end=99
NC_GLOBAL#i_parent_start=34
NC_GLOBAL#ISICE=24
NC_GLOBAL#ISLAKE=-1
NC_GLOBAL#ISOILWATER=14
NC_GLOBAL#ISURBAN=1
NC_GLOBAL#ISWATER=16
NC_GLOBAL#j_parent_end=84
NC_GLOBAL#j_parent_start=34
NC_GLOBAL#MAP_PROJ=3
NC_GLOBAL#MMINLU=USGS
NC_GLOBAL#MOAD_CEN_LAT=-8.4095078
NC_GLOBAL#NUM_LAND_CAT=24
NC_GLOBAL#parent_grid_ratio=3
NC_GLOBAL#parent_id=2
NC_GLOBAL#POLE_LAT=90
NC_GLOBAL#POLE_LON=0
NC_GLOBAL#SIMULATION_START_DATE=-00-00_00:00:00
NC_GLOBAL#SOUTH-NORTH_GRID_DIMENSION=154
NC_GLOBAL#SOUTH-NORTH_PATCH_END_STAG=154
NC_GLOBAL#SOUTH-NORTH_PATCH_END_UNSTAG=153
NC_GLOBAL#SOUTH-NORTH_PATCH_START_STAG=1
NC_GLOBAL#SOUTH-NORTH_PATCH_START_UNSTAG=1
NC_GLOBAL#sr_x=1
NC_GLOBAL#sr_y=1
NC_GLOBAL#STAND_LON=115.02
NC_GLOBAL#TITLE=OUTPUT FROM GEOGRID V3.4
NC_GLOBAL#TRUELAT1=-8.4095001
NC_GLOBAL#TRUELAT2=0
NC_GLOBAL#WEST-EAST_GRID_DIMENSION=199
NC_GLOBAL#WEST-EAST_PATCH_END_STAG=199
NC_GLOBAL#WEST-EAST_PATCH_END_UNSTAG=198
NC_GLOBAL#WEST-EAST_PATCH_START_STAG=1
NC_GLOBAL#WEST-EAST_PATCH_START_UNSTAG=1
Corner Coordinates:
Upper Left (0.0,0.0)
Lower Left (0.0, 153.0)
Upper Right ( 198.0,0.0)
Lower Right ( 198.0, 153.0)
Center ( 99.0, 76.5)
Band 1 Block=198x1 Type=Float32, ColorInterp=Undefined
NoData Value=9.96920996838686905e+36
Metadata:
description=Topography height
FieldType=104
MemoryOrder=XY
NETCDF_DIMENSION_Time=1
NETCDF_VARNAME=HGT_M
sr_x=1
sr_y=1
stagger=M
units=meters MSL
On Tue, Nov 12, 2013 at 5:44 AM, Anna Petrášová wrote:
>
>
>
> On Tue, Nov 12, 2013 at 8:38 AM, Lee Eddington
> wrote:
>
>> SUBDATASET_26_NAME=NETCDF:"geo_em.d03.nc":HGT_M
>> SUBDATASET_26_DESC=[1x153x198] HGT_M (32-bit floating-point)
>>
>>
> no, I mean when you run this command (in the directory with the netcdf
> file):
> gdalinfo NETCDF:geo_em.d03.nc:HGT_M
>
>
>
>
>>
>> On Tue, Nov 12, 2013 at 5:28 AM, Anna Petrášová wrote:
>>
>>>
>>>
>>>
>>> On Tue, Nov 12, 2013 at 12:55 AM, Lee Eddington <
>>> lee.w.edding...@gmail.com> wrote:
>>>
Anna,
I tried your suggestion, but now I'm unable to import the data. Here's
some output showing the region and the error messages:
>>>
>>>
>>> What does the gdalinfo says about the subdataset?
>>> gdalinfo NETCDF:geo_em.d03.nc:HGT_M
>>>
>>>
>>>
GRASS 6.4.3 (bali_d03_ll):~ > g.proj -p
-PROJ_INFO-
name : Lat/Lon
proj : ll
datum : wgs84
ellps : wgs84
no_defs: defined
towgs84: 0.000,0.000,0.000
-PROJ_UNITS
unit : degree
units : degrees
meters : 1.0
GRASS 6.4.3 (bali_d03_ll):~ > g.region -p
projection: 3 (Latitude-Longitude)
zone: 0
datum: wgs84
ellipsoid: wgs84
north: 6:19:40.44S
south: 10:27:59.04S
west: 112:19:46.2E
east: 117:43:20.64E
nsres: 0:01:36.12
ewres: 0:01:37.56
rows: 155
cols: 199
cells: 30845
GRASS 6.4.3 (bali_d03_ll):~ > r.in.gdal -o
input="NETCDF:grassdata/bali/geo_em.d03.nc:HGT_M" output="d03_HGT_M"
Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
WARNING: Over-riding projection check
WARNING: G_set_window(): Illegal latitude for North
On Mon, Nov 11, 2013 at 9:21 PM, Anna Petrášová
wrote:
> Hi,
>
> On Mon, Nov 11, 2013 at 11:37 PM, Lee Eddington <
> lee.w.edding...@gmail.com> wrote:
>
>> I'm trying to import a NetCDF file from the WRF weather forecast
>> model with r.in.gdal, but I can't get it to georeference. Instead I end
>> up
>> with a x,y coordinate system of rows and columns.
Re: [GRASS-user] trouble with r.in.gdal to georeference a NetCDF file
On Tue, Nov 12, 2013 at 12:55 AM, Lee Eddington
wrote:
> Anna,
>
> I tried your suggestion, but now I'm unable to import the data. Here's
> some output showing the region and the error messages:
>
What does the gdalinfo says about the subdataset?
gdalinfo NETCDF:geo_em.d03.nc:HGT_M
> GRASS 6.4.3 (bali_d03_ll):~ > g.proj -p
>
> -PROJ_INFO-
>
> name : Lat/Lon
>
> proj : ll
>
> datum : wgs84
>
> ellps : wgs84
>
> no_defs: defined
>
> towgs84: 0.000,0.000,0.000
>
> -PROJ_UNITS
>
> unit : degree
>
> units : degrees
>
> meters : 1.0
>
> GRASS 6.4.3 (bali_d03_ll):~ > g.region -p
>
> projection: 3 (Latitude-Longitude)
>
> zone: 0
>
> datum: wgs84
>
> ellipsoid: wgs84
>
> north: 6:19:40.44S
>
> south: 10:27:59.04S
>
> west: 112:19:46.2E
>
> east: 117:43:20.64E
>
> nsres: 0:01:36.12
>
> ewres: 0:01:37.56
>
> rows: 155
>
> cols: 199
>
> cells: 30845
>
> GRASS 6.4.3 (bali_d03_ll):~ > r.in.gdal -o
> input="NETCDF:grassdata/bali/geo_em.d03.nc:HGT_M" output="d03_HGT_M"
>
> Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
>
> WARNING: Over-riding projection check
>
> WARNING: G_set_window(): Illegal latitude for North
>
>
> On Mon, Nov 11, 2013 at 9:21 PM, Anna Petrášová wrote:
>
>> Hi,
>>
>> On Mon, Nov 11, 2013 at 11:37 PM, Lee Eddington <
>> lee.w.edding...@gmail.com> wrote:
>>
>>> I'm trying to import a NetCDF file from the WRF weather forecast model
>>> with r.in.gdal, but I can't get it to georeference. Instead I end up with
>>> a x,y coordinate system of rows and columns. The file georeferences fine
>>> in a number of other weather graphics programs. I tried to follow the
>>> directions at:
>>>
>>> http://www.gdal.org/frmt_netcdf.html
>>>
>>> but when I run gdalinfo I get the following:
>>>
>>> $ gdalinfo geo_em.d03.nc
>>>
>>> Warning 1: No UNIDATA NC_GLOBAL:Conventions attribute
>>>
>>> Driver: netCDF/Network Common Data Format
>>>
>>> Files: geo_em.d03.nc
>>>
>>> Size is 512, 512
>>>
>>> Coordinate System is `'
>>>
>>> Metadata:
>>>
>>> NC_GLOBAL#BOTTOM-TOP_GRID_DIMENSION=0
>>>
>>> NC_GLOBAL#CEN_LAT=-8.4095154
>>>
>>> NC_GLOBAL#CEN_LON=115.02
>>>
>>> NC_GLOBAL#corner_lats={ -10.454437, -6.3537445, -6.3537445,
>>> -10.454437, -10.454437, -6.3537445, -6.3537445, -10.454437, -10.46785,
>>> -6.3401871, -6.3401871, -10.46785, -10.46785, -6.3401871, -6.3401871,
>>> -10.46785 }
>>>
>>> NC_GLOBAL#corner_lons={ 112.3332, 112.3332, 117.70679, 117.70679,
>>> 112.31956, 112.31956, 117.72044, 117.72044, 112.3332, 112.3332, 117.70679,
>>> 117.70679, 112.31956, 112.31956, 117.72044, 117.72044 }
>>>
>>> NC_GLOBAL#DX=3000
>>>
>>> NC_GLOBAL#DY=3000
>>>
>>> NC_GLOBAL#DYN_OPT=2
>>>
>>> NC_GLOBAL#FLAG_MF_XY=1
>>>
>>> NC_GLOBAL#grid_id=3
>>>
>>> NC_GLOBAL#GRIDTYPE=C
>>>
>>> NC_GLOBAL#i_parent_end=99
>>>
>>> NC_GLOBAL#i_parent_start=34
>>>
>>> NC_GLOBAL#ISICE=24
>>>
>>> NC_GLOBAL#ISLAKE=-1
>>>
>>> NC_GLOBAL#ISOILWATER=14
>>>
>>> NC_GLOBAL#ISURBAN=1
>>>
>>> NC_GLOBAL#ISWATER=16
>>>
>>> NC_GLOBAL#j_parent_end=84
>>>
>>> NC_GLOBAL#j_parent_start=34
>>>
>>> NC_GLOBAL#MAP_PROJ=3
>>>
>>> NC_GLOBAL#MMINLU=USGS
>>>
>>> NC_GLOBAL#MOAD_CEN_LAT=-8.4095078
>>>
>>> NC_GLOBAL#NUM_LAND_CAT=24
>>>
>>> NC_GLOBAL#parent_grid_ratio=3
>>>
>>> NC_GLOBAL#parent_id=2
>>>
>>> NC_GLOBAL#POLE_LAT=90
>>>
>>> NC_GLOBAL#POLE_LON=0
>>>
>>> NC_GLOBAL#SIMULATION_START_DATE=-00-00_00:00:00
>>>
>>> NC_GLOBAL#SOUTH-NORTH_GRID_DIMENSION=154
>>>
>>> NC_GLOBAL#SOUTH-NORTH_PATCH_END_STAG=154
>>>
>>> NC_GLOBAL#SOUTH-NORTH_PATCH_END_UNSTAG=153
>>>
>>> NC_GLOBAL#SOUTH-NORTH_PATCH_START_STAG=1
>>>
>>> NC_GLOBAL#SOUTH-NORTH_PATCH_START_UNSTAG=1
>>>
>>> NC_GLOBAL#sr_x=1
>>>
>>> NC_GLOBAL#sr_y=1
>>>
>>> NC_GLOBAL#STAND_LON=115.02
>>>
>>> NC_GLOBAL#TITLE=OUTPUT FROM GEOGRID V3.4
>>>
>>> NC_GLOBAL#TRUELAT1=-8.4095001
>>>
>>> NC_GLOBAL#TRUELAT2=0
>>>
>>> NC_GLOBAL#WEST-EAST_GRID_DIMENSION=199
>>>
>>> NC_GLOBAL#WEST-EAST_PATCH_END_STAG=199
>>>
>>> NC_GLOBAL#WEST-EAST_PATCH_END_UNSTAG=198
>>>
>>> NC_GLOBAL#WEST-EAST_PATCH_START_STAG=1
>>>
>>> NC_GLOBAL#WEST-EAST_PATCH_START_UNSTAG=1
>>>
>>> Subdatasets:
>>>
>>> SUBDATASET_1_NAME=NETCDF:"geo_em.d03.nc":Times
>>>
>>> SUBDATASET_1_DESC=[1x19] Times (8-bit character)
>>>
>>> SUBDATASET_2_NAME=NETCDF:"geo_em.d03.nc":XLAT_M
>>>
>>> SUBDATASET_2_DESC=[1x153x198] XLAT_M (32-bit floating-point)
>>>
>>> SUBDATASET_3_NAME=NETCDF:"geo_em.d03.nc":XLONG_M
>>>
>>> SUBDATASET_3_DESC=[1x153x198] XLONG_M (32-bit floating-point)
>>>
>>> SUBDATASET_4_NAME=NETCDF:"geo_em.d03.nc":XLAT_V
>>> .
>>> .
>>> .
>>> .
>>>
>>> Corner Coordinates:
>>>
>>> Upper Left (0.0,0.0)
>>>
>>> Lower Left (0.0, 512.0)
>>>
>>> Upper Right ( 512.0,0.0)
>>>
>>> Lower Right ( 512.0, 512.0)
>>>
>>> Center ( 256.0, 256.0)
>>>
>>> The wa
