Changes between Version 1 and Version 2 of UserDocs/RasterProcTutorial

Timestamp:

Sep 13, 2007, 12:23:55 PM (17 years ago)

Author:

warmerdam

Comment:

--

UserDocs/RasterProcTutorial

@@ -25 +25 @@
 ''Note to self, need to add geogratis and/or geobase license info on download site for cded and CanMatrix data.''
 
-== Exploring Your Data ==
+== Tutorial ==
+
+''format nicely for wiki later - keep plain ascii for now for easy printout at conference.''
 
 {{{
-- Use gdalinfo to review several datasets. 
-- Highlight:
-   image size
-   geotransform
-   gcps
-   coordinate system
-   metadata
-   band type
-   block size
-   nodata value
-   color table
-
-- Use OpenEV to view datasets. 
-- Highlight: 
-   finding exact location of a pixel (and subpixel position issues)
-   Seeing pixel greyscale values (and for multiple bands)
-   "see through" due nodata and alpha bands. 
-   Overlay flipping
-
-Perhaps also demonstrate something similar with QGIS? 
-}}}
-
-== Simple Format Translation ==
-
-=== Know your Drivers ===
-
-The --formats commandline switch of gdal_translate can be used to see a list of available format drivers.  Each format reports if it is read only (ro), read/write (rw) or read/write/update (rw+).
-
-{{{
+GDAL Overview
+-----------------
+
+5 minute introduction to GDAL. 
+
+
+1 Exploring Your Image Data
+---------------------------
+
+1a) Use gdalinfo to get information about the dataset.
+
+  gdalinfo world.png
+
+  Note: 
+   - Driver is "PNG"
+   - Size is 2048x1024
+   - 3 Bands of type Byte. 
+   - No coordinate system and bounds are just pixel/line based.
+
+1b) More with gdalinfo.
+ 
+   gdalinfo -mm east.dem
+
+   Note: 
+    - Coordinate system is NAD83 (expressed as WKT). 
+    - Origin, Pixel Size and bounds are all provided.
+    - Band is now of type "Int16"
+    - Nodata Value is -32767
+    - Unit type is "m" 
+    - Min/max elevation values shown (computed due to -mm switch)
+
+1c) Even more with gdalinfo: 
+
+   gdalinfo 092b05.tif
+
+   Note: 
+    - Extra bits of TIFF metadata. 
+    - Notes the compression 
+    - UTM corner coordinates also reported as lat/long. 
+    - Color table
+
+1d) Look at your data!
+
+    openev world.png
+
+  Note: No georeferenced coordinates 
+   
+    openev east.dem
+
+  Note: You can point at pixels and see the elevations, and everything is
+  georeferenced.
+
+    openev 092b05.tif
+
+  Right click on the layer in layer dialog and look at the "Coordinate
+  System" and "Image Info" tabs for some info similar to gdalinfo.
+
+1e) If you have a minute, try looking at some of the files with QGIS,
+    uDig, etc 
+
+
+2 Simple Format Translation
+---------------------------
+
+-- Know your Drivers --
+
+2a) The --formats commandline switch of gdal_translate can be used to see a
+    list of available format drivers.  Each format reports if it is read 
+    only (ro), read/write (rw) or read/write/update (rw+).
+
  gdal_translate --formats
-}}}
-
 
 2b) The --format commandline switch can be used to query details about a 
@@ -100 +141 @@
     gdal_translate -of JPEG -scale -100 678 0 255 east.dem east.jpg
 
+
+3 Correcting Metadata
+---------------------
+
+-- Set Coordinate System and Bounds --
+
+3a) Use gdal_translate to attach bounds, and coordinate system:
+
+  gdal_translate -a_srs WGS84 -a_ullr -180 90 180 -90 world.png geoworld.tif
+
+3b) Confirm reporting of coordinate system and georeferencing with gdalinfo:
+
+  gdalinfo geoworld.tif
+
+3c) Confirm reporting of georeferenced locations with OpenEV (optional)
+
+  openev geoworld.tif
+
+-- Reset NODATA values --
+
+3d) Use gdal_translate to attach/change the nodata value of a file.
+
+  gdal_translate -a_nodata 0 west.dem zeronodata.tif
+
+3e) Use gdalinfo to verify nodata value differs between the two files.
+
+  gdalinfo west.dem
+  gdalinfo zeronodata.tif
+
+
+Reprojecting
+------------
+
+For this process we assume that geoworld.tif has been properly created with
+bounds and coordinate system info in the previous steps.  
+
+a) The gdalwarp command can be used to reproject images.  Here we reproject
+   the WGS84 geographic image to the Mercator projection:
+
+   gdalwarp -t_srs '+proj=merc +datum=WGS84' geoworld.tif mercator.tif
+
+   Compare the images with OpenEV
+
+   openev mercator.tif
+   openev geoworld.tif
+
+b) Here we reproject to the Ortho projection.  
+
+   gdalwarp -t_srs '+proj=ortho +datum=WGS84' geoworld.tif ortho.tif
+  
+   openev ortho.tif
+
+c) Note how the poles are clipped?  This is because the edges at the pole
+   can't be reprojected gdalwarp does not read all the data.  We can force
+   gdalwarp to read a bunch of surplus data around chunks as one way to 
+   resolve this. 
+
+   gdalwarp -wo SOURCE_EXTRA=125 -t_srs '+proj=ortho +datum=WGS84' 
+             geoworld.tif ortho.tif
+
+d) gdalwarp can also be cause to treat particular values as nodata and to 
+   produce alpha values in the output.  In this example we cause oceans to
+   be treated as transparent, and generate alpha in the output. 
+
+   gdalwarp -wo SOURCE_EXTRA=125 -srcnodata "11 10 50" -dstalpha
+             -t_srs '+proj=ortho +datum=WGS84' 
+             geoworld.tif ortho.tif
+
+   Hmm, it seems the water isn't quite as uniformly blue as I hoped.  It
+   must have lived as a lossily compressed image at one point!
+
+
+Mosaicing
+---------
+
+a) gdal_merge.py is a python script that can be used for simple mosaicing
+   tasks.  Mosaic the east.dem and west.dem into a single file:
+
+   gdal_merge  east.dem west.dem -o mergeddem.tif
+
+   For extra credit, open gdal_merge.py in notepad and review.  It's scary,
+   but not as scary as you might think!
+
+   notepad "C:\Program Files\FWTools-1.3.6\pymod\gdal_merge.py"
+
+b) The same task can be accomplished with gdalwarp.  gdalwarp has a variety 
+   of advantages over gdal_merge, but can be slow to merge many files:
+
+   gdalwarp east.dem west.dem warpmerged.tif
+
+c) Now lets mosaic two CanMatrix image for our area:
+
+   gdalwarp 092b05.tif 092b06.tif mosaic1.tif
+
+   Seems successful, but if we look at the image with OpenEV ... a mess!
+
+   openev mosaic1.tif
+
+d) The problem is that the images have different pseudocolor tables.  So 
+   lets convert from pseudocolor to RGB and then mosaic.  To speed up 
+   later steps we will also reduce the resolution by half.
+
+   gdal_translate -outsize 50% 50% 092b05.tif small5.tif
+   gdal_translate -outsize 50% 50% 092b06.tif small6.tif
+
+   pct2rgb small5.tif rgb5.tif
+   pct2rgb small6.tif rgb6.tif
+
+   gdalwarp rgb5.tif rgb6.tif mosaic2.tif
+
+   Hmm, better but clearly we have a problem with collar overlaps. 
+
+e) So, lets try stripping the collars off.  We can do this manually by
+   digitizing a polygon representing the "valid" area of the map.  This
+   can be done in OpenEV, QGIS, etc, but for our purposes we provide a 
+   pre-created cutline.  Review it with:
+
+   openev rgb06.tif cutline06.shp
+
+   "Burn away" the inverse of the polygon using the command:
+   
+gdal_rasterize -i -burn 17 -b 1 -b 2 -b 3 cutline06.shp -l cutline06 rgb6.tif 
+
+   Review the result, we should have "17,17,17" for all outside pixels.
+
+f) Now mosaic again, but request that "17" be treated as a nodata value.
+
+   gdalwarp -srcnodata "17 17 17" rgb5.tif rgb6.tif mosaic3.tif
+
+   A pretty good result, especially if we were to also strip the collar 
+   off the first image.  
+
+g) For extra points, rerun the warp with the -dstalpha option to add an
+   alpha (transparency) band for nodata areas. 
+
+   gdalwarp -dstalpha -srcnodata "17 17 17" rgb5.tif rgb6.tif mosaic4.tif 
+
+
+Optimizing  
+----------
+
+... produce an optimized version of the previous mosaic. 
+
+ideally demonstrate effect using a mapfile and shp2img? 
+
+
+Virtual Files
+-------------
+
+}}}
+
+