Re: [GRASS-user] Re: Importing Dems with r.in.xyz
2010/6/5, Hamish hamis...@yahoo.com: Hanlie wrote: At this point, g.region reports 1146474 cells in the region, while I have 1146370 lines of coordinates in my file. ... So it looks like there are about 100 coordinates missing from the ASCII ASCII file. 0.01% .. Maybe holes in the data? perhaps this: https://trac.osgeo.org/grass/ticket/123 ?? I don't think it's this bug because this bug discards only one line of data. I don't get any data in because the number of coordinate pairs in the file is less than the number of cells in the defined region. I was thinking perhaps importing the points as vectors, converting them to raster and then doing a nearest neighbour or IDW interpolation to fill the gaps. At least then I'll be able to see where the gaps are and limit the interpolated pixels using a mask? No need to do anything different to find the missing pixels. Inspecting the output of r.univar with r.in.xyz's method=n maps can be very useful for troubleshooting. from the help page: Gridded data If data is known to be on a regular grid r.in.xyz can reconstruct the map perfectly as long as some care is taken to set up the region correctly and that the data's native map projection is used. A typical method would involve determining the grid resolution either by examining the data's associated documentation or by studying the text file. Next scan the data with r.in.xyz's -s (or -g) flag to find the input data's bounds. GRASS uses the cell-center raster convention where data points fall within the center of a cell, as opposed to the grid-node convention. Therefore you will need to grow the region out by half a cell in all directions beyond what the scan found in the file. After the region bounds and resolution are set cor- rectly with g.region, run r.in.xyz using the n method and verify that n=1 at all places. r.univar can help. Once you are confident that the region exactly matches the data proceed to run r.in.xyz using one of the mean, min, max, or median methods. With n=1 throughout, the result should be identical regardless of which of those methods are used. with the n map you might use r.mapcalc to extract the NULL cells as some value, then r.out.xyz or r.to.vect on th extracts to highlight where they are. Or maybe you get lucky with r.colors with nv set to bright magenta on the original data. Thanks, I'll try this to find where the holes in the data are. Hamish ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Re: Importing Dems with r.in.xyz
Hanlie Pretorius wrote: 2010/6/5, Hamish hamis...@yahoo.com: Hanlie wrote: At this point, g.region reports 1146474 cells in the region, while I have 1146370 lines of coordinates in my file. ... So it looks like there are about 100 coordinates missing from the ASCII ASCII file. 0.01% .. Maybe holes in the data? perhaps this: https://trac.osgeo.org/grass/ticket/123 ?? I don't think it's this bug because this bug discards only one line of data. I don't get any data in because the number of coordinate pairs in the file is less than the number of cells in the defined region. Weird. In 6.4, r.in.xyz does import a file where the number of coordinate pairs is far less than the number of cells in the defined region. (I just did a simple test with two input lines and a region with 26.5 million cells, got clean import and correct result) You can interpolate NULL cells and a the same time keep non NULL cell values with r.fillnulls. I would suggest to set the region to north: -49312.5 south: -74587.5 west: -3015862.5 east: -2987512.5 nsres: 25 ewres: 25 have clean 25m resolution for both ns and ew in order to make life easier. The 1cm difference in the last input lines you posted can not possibly make a difference with 25m point spacing and is most probably some floating point rounding error introduced by some preprocessing to generate the ascii xyz input file. HTH, Markus M I was thinking perhaps importing the points as vectors, converting them to raster and then doing a nearest neighbour or IDW interpolation to fill the gaps. At least then I'll be able to see where the gaps are and limit the interpolated pixels using a mask? No need to do anything different to find the missing pixels. Inspecting the output of r.univar with r.in.xyz's method=n maps can be very useful for troubleshooting. from the help page: Gridded data If data is known to be on a regular grid r.in.xyz can reconstruct the map perfectly as long as some care is taken to set up the region correctly and that the data's native map projection is used. A typical method would involve determining the grid resolution either by examining the data's associated documentation or by studying the text file. Next scan the data with r.in.xyz's -s (or -g) flag to find the input data's bounds. GRASS uses the cell-center raster convention where data points fall within the center of a cell, as opposed to the grid-node convention. Therefore you will need to grow the region out by half a cell in all directions beyond what the scan found in the file. After the region bounds and resolution are set cor- rectly with g.region, run r.in.xyz using the n method and verify that n=1 at all places. r.univar can help. Once you are confident that the region exactly matches the data proceed to run r.in.xyz using one of the mean, min, max, or median methods. With n=1 throughout, the result should be identical regardless of which of those methods are used. with the n map you might use r.mapcalc to extract the NULL cells as some value, then r.out.xyz or r.to.vect on th extracts to highlight where they are. Or maybe you get lucky with r.colors with nv set to bright magenta on the original data. Thanks, I'll try this to find where the holes in the data are. Hamish ___ 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] Re: Importing Dems with r.in.xyz
Hanlie Pretorius wrote: I don't think it's this bug because this bug discards only one line of data. I don't get any data in because the number of coordinate pairs in the file is less than the number of cells in the defined region. Markus Metz: Weird. In 6.4, r.in.xyz does import a file where the number of coordinate pairs is far less than the number of cells in the defined region. (I just did a simple test with two input lines and a region with 26.5 million cells, got clean import and correct result) sorry, I'm not really grasping what the problem is. If you suspect there is something weird going on can you make the .xyz.bz2 file available for me to download, and supply the output of r.univar `wc -l` and `r.in.xyz --verbose`? Hamish ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Re: Importing Dems with r.in.xyz
2010/6/7, Hamish hamis...@yahoo.com: Hanlie Pretorius wrote: I don't think it's this bug because this bug discards only one line of data. I don't get any data in because the number of coordinate pairs in the file is less than the number of cells in the defined region. Markus Metz: Weird. In 6.4, r.in.xyz does import a file where the number of coordinate pairs is far less than the number of cells in the defined region. (I just did a simple test with two input lines and a region with 26.5 million cells, got clean import and correct result) sorry, I'm not really grasping what the problem is. If you suspect there is something weird going on can you make the .xyz.bz2 file available for me to download, and supply the output of r.univar `wc -l` and `r.in.xyz --verbose`? Hamish Ok, the text file (6.1MB compressed) that I use as input to r.in.xyz is available at: http://www.nedbib.za.net/dems/ Output of r.in.xyz --verbose: - Reading data ... Writing to map ... r.in.xyz complete. 0 points found in region. - Output of r.univar: - total null and non-null cells: 1146474 total null cells: 1146474 Of the non-null cells: -- n: 0 minimum: -nan maximum: -nan range: -nan mean: -nan mean of absolute values: -nan standard deviation: -nan variance: -nan variation coefficient: -nan % sum: 0 - output of wc -l: 1146370 Output of r.info: - ++ | Layer:dem_2728ab_25m Date: Mon Jun 7 12:09:34 2010 | Mapset: PERMANENT Login of Creator: hanlie | Location: sa_tm_29deg_E | DataBase: /media/0847147784/grassdata | Title:Raw x,y,z data binned into a raster grid by cell mean ( dem_2728 | Timestamp: none | | | Type of Map: raster Number of Categories: 255 | Data Type:FCELL | Rows: 1011 | Columns: 1134 | Total Cells: 1146474 |Projection: Transverse Mercator |N: -49312.5S: -74587.5 Res:25 |E: -2987512.5W: -3015862.5 Res:25 | Range of data:min = -nan max = -nan | | Data Source: |/media/0847147784/data/CD-NGI/DEMS/25m dems/2728/ab/2728AB.ort.grass | | | Data Description: |generated by r.in.xyz | | Comments: |r.in.xyz input=/media/0847147784/data/CD-NGI/DEMS/25m dems/2728/ab/\ |2728AB.ort.grass output=dem_2728ab_25m method=mean type=FCELL\ | x=1 y=2 z=3 zscale=1.0 percent=100 - Thanks Hanlie ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Re: Importing Dems with r.in.xyz
g.region -p north: -49312.5 south: -74587.5 west: -3015862.5 east: -2987512.5 nsres: 25 ewres: 25 rows: 1011 cols: 1134 cells: 1146474 From the region settings you posted previously, I assume that the 1. column is north and the 2. is east, so: r.in.xyz input=2728AB.ort.grass output=2728AB.ort fs=, x=2 y=1 r.in.xyz complete. 1146370 points found in region. r.univar map=2728AB.ort total null and non-null cells: 1146474 total null cells: 104 Of the non-null cells: -- n: 1146370 minimum: 1473.52 maximum: 1618.85 range: 145.33 mean: 1519.76 mean of absolute values: 1519.76 standard deviation: 26.6106 variance: 708.125 variation coefficient: 1.75098 % sum: 1742202717.5849609375 Everything looks ok. Markus M On Mon, Jun 7, 2010 at 12:33 PM, Hanlie Pretorius hanlie.pretor...@gmail.com wrote: 2010/6/7, Hamish hamis...@yahoo.com: Hanlie Pretorius wrote: I don't think it's this bug because this bug discards only one line of data. I don't get any data in because the number of coordinate pairs in the file is less than the number of cells in the defined region. Markus Metz: Weird. In 6.4, r.in.xyz does import a file where the number of coordinate pairs is far less than the number of cells in the defined region. (I just did a simple test with two input lines and a region with 26.5 million cells, got clean import and correct result) sorry, I'm not really grasping what the problem is. If you suspect there is something weird going on can you make the .xyz.bz2 file available for me to download, and supply the output of r.univar `wc -l` and `r.in.xyz --verbose`? Hamish Ok, the text file (6.1MB compressed) that I use as input to r.in.xyz is available at: http://www.nedbib.za.net/dems/ Output of r.in.xyz --verbose: - Reading data ... Writing to map ... r.in.xyz complete. 0 points found in region. - Output of r.univar: - total null and non-null cells: 1146474 total null cells: 1146474 Of the non-null cells: -- n: 0 minimum: -nan maximum: -nan range: -nan mean: -nan mean of absolute values: -nan standard deviation: -nan variance: -nan variation coefficient: -nan % sum: 0 - output of wc -l: 1146370 Output of r.info: - ++ | Layer: dem_2728ab_25m Date: Mon Jun 7 12:09:34 2010 | Mapset: PERMANENT Login of Creator: hanlie | Location: sa_tm_29deg_E | DataBase: /media/0847147784/grassdata | Title: Raw x,y,z data binned into a raster grid by cell mean ( dem_2728 | Timestamp: none | | | Type of Map: raster Number of Categories: 255 | Data Type: FCELL | Rows: 1011 | Columns: 1134 | Total Cells: 1146474 | Projection: Transverse Mercator | N: -49312.5 S: -74587.5 Res: 25 | E: -2987512.5 W: -3015862.5 Res: 25 | Range of data: min = -nan max = -nan | | Data Source: | /media/0847147784/data/CD-NGI/DEMS/25m dems/2728/ab/2728AB.ort.grass | | | Data Description: | generated by r.in.xyz | | Comments: | r.in.xyz input=/media/0847147784/data/CD-NGI/DEMS/25m dems/2728/ab/\ | 2728AB.ort.grass output=dem_2728ab_25m method=mean type=FCELL\ | x=1 y=2 z=3 zscale=1.0 percent=100 - Thanks Hanlie ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
[GRASS-user] Re: Importing Dems with r.in.xyz
2010/6/7, Markus Metz markus.metz.gisw...@googlemail.com: g.region -p north: -49312.5 south: -74587.5 west: -3015862.5 east: -2987512.5 nsres: 25 ewres: 25 rows: 1011 cols: 1134 cells: 1146474 From the region settings you posted previously, I assume that the 1. column is north and the 2. is east, so: r.in.xyz input=2728AB.ort.grass output=2728AB.ort fs=, x=2 y=1 r.in.xyz complete. 1146370 points found in region. r.univar map=2728AB.ort total null and non-null cells: 1146474 total null cells: 104 Of the non-null cells: -- n: 1146370 minimum: 1473.52 maximum: 1618.85 range: 145.33 mean: 1519.76 mean of absolute values: 1519.76 standard deviation: 26.6106 variance: 708.125 variation coefficient: 1.75098 % sum: 1742202717.5849609375 Everything looks ok. Markus M Ahh, you have found the error. In fact, the columns should not be swopped, the region specifications should. I've retried the import and it works now despite the missing data and the DEM sits in the right place compared to the rest of my data. Thanks for everyone's help and my apologies for the confusion. Our national projection swops X and Y and reverses the direction of the positive coordinates, so it can get very confusing. Hanlie ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
[GRASS-user] Re: Importing Dems with r.in.xyz
2010/6/4, Micha Silver mi...@arava.co.il: On 06/04/2010 04:20 PM, Hanlie Pretorius wrote: Hi, I've been following the procedure at http://grass.osgeo.org/wiki/Import_XYZ to import DEMs into GRASS. It has worked for two DEMS, but I have a problem with the third one at the step where one verifies that the number of rows in the ASCII file corresponds to the number of cells in the enlarged region. At this point, g.region reports 1146474 cells in the region, while I have 1146370 lines of coordinates in my file. Assuming you already did the r.in.xyz -s ... step to get and set your region to match the input data. So it looks like there are about 100 coordinates missing from the ASCII file. Maybe holes in the data? One way to work around this might be to import the point data as a 3D vector, then run the usual v.surf.rst interpolation. v.in.xyz -z in=ascii.txt out=elev_pts z=3 fs=, then g.region vect=elev_pts res=choose appropriate resolution v.surf.rst elev_pts elev=dem layer=0 Thanks for the suggestion, Micha. I've tried the interpolation route before and found that the difference between the resulting interpolated surface and the original DEM was up to 7m. I want to use the DEMs for a flood application, so they need to be as accurate as possible. I was thinking perhaps importing the points as vectors, converting them to raster and then doing a nearest neighbour or IDW interpolation to fill the gaps. At least then I'll be able to see where the gaps are and limit the interpolated pixels using a mask? The output of g.region is: - projection: 99 (Transverse Mercator) zone: 0 datum: ** unknown (default: WGS84) ** ellipsoid: wgs84 north: -49312.49 south: -74587.5 west: -3015862.5 east: -2987512.5 nsres: 25.0989 ewres: 25 rows: 1011 cols: 1134 cells: 1146474 - The first three and last three coordinates in my text file are: - -74575.00,-3015850.00,1548.83 -74575.00,-3015825.00,1548.33 -74575.00,-3015800.00,1547.50 . . . -49324.99,-2987575.01,1510.98 -49324.99,-2987550.01,1511.24 -49324.99,-2987525.01,1511.47 - Can someone help me to figure out what's going on? Thanks Hanlie ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user This mail was received via Mail-SeCure System. -- Micha Silver Arava Development Co. +972-52-3665918 http://www.surfaces.co.il ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
Re: [GRASS-user] Re: Importing Dems with r.in.xyz
Hanlie wrote: At this point, g.region reports 1146474 cells in the region, while I have 1146370 lines of coordinates in my file. ... So it looks like there are about 100 coordinates missing from the ASCII ASCII file. 0.01% .. Maybe holes in the data? perhaps this: https://trac.osgeo.org/grass/ticket/123 ?? I was thinking perhaps importing the points as vectors, converting them to raster and then doing a nearest neighbour or IDW interpolation to fill the gaps. At least then I'll be able to see where the gaps are and limit the interpolated pixels using a mask? No need to do anything different to find the missing pixels. Inspecting the output of r.univar with r.in.xyz's method=n maps can be very useful for troubleshooting. from the help page: Gridded data If data is known to be on a regular grid r.in.xyz can reconstruct the map perfectly as long as some care is taken to set up the region correctly and that the data's native map projection is used. A typical method would involve determining the grid resolution either by examining the data's associated documentation or by studying the text file. Next scan the data with r.in.xyz's -s (or -g) flag to find the input data's bounds. GRASS uses the cell-center raster convention where data points fall within the center of a cell, as opposed to the grid-node convention. Therefore you will need to grow the region out by half a cell in all directions beyond what the scan found in the file. After the region bounds and resolution are set cor- rectly with g.region, run r.in.xyz using the n method and verify that n=1 at all places. r.univar can help. Once you are confident that the region exactly matches the data proceed to run r.in.xyz using one of the mean, min, max, or median methods. With n=1 throughout, the result should be identical regardless of which of those methods are used. with the n map you might use r.mapcalc to extract the NULL cells as some value, then r.out.xyz or r.to.vect on th extracts to highlight where they are. Or maybe you get lucky with r.colors with nv set to bright magenta on the original data. Hamish ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user
[GRASS-user] Re: Importing Dems with r.in.xyz
On 06/05/2010 09:45 AM, Hanlie Pretorius wrote: 2010/6/4, Micha Silvermi...@arava.co.il: On 06/04/2010 04:20 PM, Hanlie Pretorius wrote: Hi, I've been following the procedure at http://grass.osgeo.org/wiki/Import_XYZ to import DEMs into GRASS. It has worked for two DEMS, but I have a problem with the third one at the step where one verifies that the number of rows in the ASCII file corresponds to the number of cells in the enlarged region. At this point, g.region reports 1146474 cells in the region, while I have 1146370 lines of coordinates in my file. Assuming you already did the r.in.xyz -s ... step to get and set your region to match the input data. So it looks like there are about 100 coordinates missing from the ASCII file. Maybe holes in the data? One way to work around this might be to import the point data as a 3D vector, then run the usual v.surf.rst interpolation. v.in.xyz -z in=ascii.txt out=elev_pts z=3 fs=, then g.region vect=elev_pts res=choose appropriate resolution v.surf.rst elev_pts elev=dem layer=0 Thanks for the suggestion, Micha. I've tried the interpolation route before and found that the difference between the resulting interpolated surface and the original DEM was up to 7m. I want to use the DEMs for a flood application, so they need to be as accurate as possible. I wonder if playiing with the tension and smooth parameters would help? Higher tension (40) means that the interpolation acts more like a rubber sheet, so each point influences a maller area, and smooth=0 means that the final dem must go exactly thru each point. I was thinking perhaps importing the points as vectors, converting them to raster and then doing a nearest neighbour or IDW interpolation to fill the gaps. At least then I'll be able to see where the gaps are and limit the interpolated pixels using a mask? Worth a try. But as you probably know, nearest neighbor won't take trends in the surface into account. IDW the same. And IDW usually gives worse terrain results than RST. The output of g.region is: - projection: 99 (Transverse Mercator) zone: 0 datum: ** unknown (default: WGS84) ** ellipsoid: wgs84 north: -49312.49 south: -74587.5 west: -3015862.5 east: -2987512.5 nsres: 25.0989 ewres: 25 rows: 1011 cols: 1134 cells: 1146474 - The first three and last three coordinates in my text file are: - -74575.00,-3015850.00,1548.83 -74575.00,-3015825.00,1548.33 -74575.00,-3015800.00,1547.50 . . . -49324.99,-2987575.01,1510.98 -49324.99,-2987550.01,1511.24 -49324.99,-2987525.01,1511.47 - Can someone help me to figure out what's going on? Thanks Hanlie ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user This mail was received via Mail-SeCure System. -- Micha Silver Arava Development Co. +972-52-3665918 http://www.surfaces.co.il This mail was received via Mail-SeCure System. -- Micha Silver Arava Development Co. +972-52-3665918 http://surfaces.co.il ___ grass-user mailing list grass-user@lists.osgeo.org http://lists.osgeo.org/mailman/listinfo/grass-user