Monday, September 27, 2010

Google Gadget for monitoring earthquakes (USGS source)

This gadget will allow you to subscribe to the US Geological Survey (USGS) real-time, worldwide earthquake list. It shows the latest magnitude 1 and above earthquakes for the past day as color coded (by depth and age) icons with tool tips on a Google Maps backdrop - the older the earthquake, the more translucent the icon will look. 

The gadget has a list tab which displays the earthquake information as a text list with hyperlinks for more details - the list can be sorted by the earthquake magnitude, depth as well as the date-time. Using the list, you can easily find and locate the largest, the latest or the deepest earthquakes for the past day. 

You can leave this gadget on the screen and it will update itself as it receives the latest feed. Upon receiving the data, it will automatically zoom to the latest earthquake epicenter. Clicking on an earthquake icon will bring up more information about the earthquake from the USGS site.

Monday, September 20, 2010

Isolating hot spots by reclassification and zone area calculations in SAGA GIS

Another way to summarize and visualize a continuous grid surface e.g. a smooth surface that modeled the density of occurrences of point data incidents, is to reclassify the continuous data into discrete groups. An example of a continuous data surface is shown below.

In SAGA GIS, the Reclassify Grid Values module can do the zoning. Once the continuous surface has been classified into groups, meaningful numerical information such as the area of the hot spots can be derived.

The following describes how I use SAGA GIS to reclassify and calculate the zone area values.

Reclassifying into discrete zones
  1. Start up SAGA GIS. Load in a continuous grid layer e.g. hotspot.sgrd.
  2. In the Workspace pane on the left, select the grid layer. In the Options pane on the right, click the Description pane.
  3. Examine the Value Range and note down the values, e.g. 0 - 10. Note: this will be used to determine how to group the data.
  4. Select Modules | Grid | Tools | Values | Reclassify Grid Values.

    The Reclassify Grid Values dialog box appears.
  5. In the Grid system drop down list, choose the grid system e.g. "30;354x420y...".
  6. In the Grid drop down list, choose the grid layer for reclassification, e.g. "01.hotspot".
  7. In the Method drop down list, choose simple table.

  8. Click the Method simple table's Lookup Table row. Click [...].

    The Lookup Table dialog box appears.
  9. Use the minimum, maximum, and new columns to define the input ranges to reclassify, e.g. 0~3 to 1, 3~6 to 2 and so on.

  10. Click Okay twice.

    The grid layer is reclassified.
Visualizing the discrete zones
  1. In the Workspace pane, double click the discrete zone grid layer e.g. "03.Reclassified Grid" and display in a new map window.

    Note: the default colors may not be optimal in showing off the zones as shown below.
  2. On the Options pane on the right, click the Settings tab. Change the Display:Color Classification Type to Lookup Table.

  3. Click the Lookup Table's Table row. Click [...].

    The Table dialog box appears.
  4. Enter the NAME, DESCRIPTION, MINIMUM, MAXIMUM values appropriate for the discrete zones grid layer and define unique colors for each range, as shown below.

  5. Click Okay.
  6. Click Apply.

    The map window is updated with the new colors which show off the zones better.
Calculating the zone areas
  1. Select Modules | Geostatistics | Grids | Zonal Grid Statistics.

    The Zonal Grid Statistics dialog box appears.
  2. Click the Grid system drop down list and choose the discrete grid layer's system.
  3. In the Zone Grid drop down list, choose the discrete zone grid layer, e.g. "03.Reclassified Grid".

  4. Click Okay.

    The zonal statistics table is created.
  5. In the Workspace pane, double click the newly created zonal statistics table e.g. "01.Zonal Statistics" to display it.

    Note: the Count column counts the number of grid cells within a zone.
  6. Select Modules | Table | Calculus | Table Calculator.

    The Table calculator dialog box appears.
  7. In the Table drop down list, choose the newly create table, e.g. "01.Zonal Statistics".
  8. In the Formula field, type in the formula to calculate the zone area, e.g. if the area of a grid cell is 2 square meter, then the zone area formula is b*2, where b refers to the second column.
  9. In the Field Name, type in the new field name, e.g. area.

  10. Click Okay.

    The new table is generated.
  11. In the Workspace pane, double click on the newly created table.

    The new table is displayed. Note the new area field showing the calculated zone area.

Monday, September 13, 2010

Isolating hot spots using Isolines in Saga GIS

In my previous post, I generated a surface grid layer of total incident occurrences where a graduated color sequence indicate areas with the highest density of occurrences i.e. the hot spots in red, as shown in the figure below.

We can further demarcate the hot spots with SAGA GIS's Contour Lines from Grid module.


  1. Start up SAGA GIS and load in the grid layer of the total incident occurrences.
  2. Select Modules | Shapes | Grid | Contour Lines from Grid.

    The Contour Lines from Grid dialog box appears.
  3. In the Grid system drop down list, choose the grid system of the grid layer, e.g. "30; 354x 420y...".
  4. Click the Grid row drop down list and choose the grid layer, e.g. "03.ifelse(lt(a,0),0,a)".
  5. If necessary, change the Minimum Contour Value, Maximum Contour Value and Equidistance values.

  6. Click Okay.

    The contour lines are generated.

    Note: in the Workspace pane, double click the contour line shape line layer to display it overlaid on the incident occurrence grid layer, as shown below.


    Labeling the contour lines
  7. In the Workspace pane, select the Line shape layer of the contour lines e.g. "01.ifelse(lt(a,0),0,a)".
  8. In the Options pane on the right, click the Display:Label Attribute drop down list and choose the attribute field to use for labeling, e.g. "ifelse(lt(a,0),0,a)".

  9. Click Apply.

    The contour lines are labeled.

Wednesday, September 8, 2010

Intergraph introduces GeoMedia 3D

I haven't got a chance  to try out GeoMedia 3D but it certainly looks promising. It's an add-on to the core GeoMedia product, which adds in an additional 3D Map window seamlessly into the framework and you could manipulate and query the features just like in the basic 2D map window.

Monday, September 6, 2010

Detecting hot spot patterns in point incidents data with SAGA GIS

SAGA GIS has a number of grid analysis commands that can be used to detect patterns (i.e. find hot spots) in point data e.g. accident incidents. An example incident distribution on a map is shown in the figure below.

The point vector data may be storing additional attributes about the incident such as the number of occurrences.

Just by looking at the vector point data, it is impossible to detect a pattern or hot spot incidents in the area. By forming a graphical occurrence distribution surface from the vector point data with a grid analysis tool like SAGA GIS, it would be so much easier to visualize the areas of high incident density. An example of this surface is shown below, overlaid with other map features for reference.



Below are the steps I took to generate the grid surface for visualizing the incident patterns.

Import vector point data

  1. Start up SAGA GIS
  2. Select Modules | File | GDAL/OGR | OGR: Import Vector Data.

    The OGR: Import Vector Data dialog box appears.
  3. Click the Files row. Click the [...] button on the right. Browse and select the point vector data file e.g. Incidents.shp. Click Open.

    The selected file is displayed in the OGR: Import Vector Data dialog box.

  4. Click Okay.

    The point data is imported into SAGA GIS.

Convert the vector point data to a grid layer
  1. Select Modules | Grid | Gridding | Shapes to Grid.

    The Shapes to Grid dialog box appears.
  2. Click the Shapes drop down list and select the previously imported vector points e.g. 01.Incidents.
  3. Click the Attribute drop down list and select the field to get the new grid cell values, e.g. OCCURRENCE.
  4. Click the Preferred Target Grid Type drop down list and select a suitable data type e.g. Integer (2 byte).


  5. Click Okay.

    The User Defined Grid dialog box appears.
  6. In the Cellsize field, type in the new cell size e.g. 30.


  7. Click Okay.

    The vector points are rasterized into a new grid layer.
Workaround the Aggregate module's Sum No Data cell value bug
SAGA GIS 2.0.5's Aggregate module has a problem summing grid cells with negative values, i.e. cells with no data and usually represented with a large negative number e.g. -99999. In order for the SUM function to work, those large negative values have to be set to 0. This can be done by using the Calculus module as described below.
  1. In the Workspace pane, select the newly created point grid layer e.g. 01.Incidents [OCCURRENCE].


  2. In the Object Properties pane on the right, click the Settings tab. Then change the No Data field values to "0; 0".


  3. Click Apply.
  4. Select Modules | Grid | Calculus | Grid Calculator.

    The Grid Calculator dialog box appears.

  5. In the Grid system drop down list, select the newly created grid layer's system, e.g. "30;354x 420y;...".
  6. Select the Grids row. Click the [...] button

    The Grids dialog box appears.

  7. On the list on the left side, select the newly created grid layer, e.g. 01.Incidents[OCCURRENCE]. Click the >> button.

    The selected grid layer is moved to the list on the right.
  8. Click Okay.

    The Grids dialog box is closed.
  9. In the Formula field, type in ifelse(lt(a,0),0,a).

    Note: this formula tells the Calculator to replace all grid cell values below 0 with the value 0.

  10. Click Okay

    A new grid layer is created with no negative values.

Using Aggregate to sum adjacent grid cells
  1. Select Modules | Grid | Construction | Aggregate.

    The Aggregate dialog box appears.
  2. In the Grid system drop down list, choose the previously created grid layer's system, e.g. "30;354x 420y;...".
  3. In the Grid drop down list, select the newly created grid layer without negative values, e.g. ifelse(lt(a,0),0,a).
  4. In the Aggregation Size field, type in an appropriate size (in cells), e.g. 20. In the Method drop down list, choose Sum.


  5. Click Okay.

    The aggregated grid layer is created.

    This is how the aggregated layer looks like in the map view.

Resample the aggregated grid layer to match the original point grid layer.
  1. Select Modules | Grid | Construction | Resampling.

    The Resampling dialog box appears.
  2. In the Grid system drop down list, select the aggregated grid layer's system, e.g. "600;17x21y...".
  3. In the Grid drop down list, choose the aggregated grid layer, e.g. "ifelse(lt(a,0),0,a)".
  4. In the Target Grid drop down list, choose grid.


  5. Click Okay.

    The Choose Grid dialog box appears.
  6. In the Grid system drop down list, choose the original point grid layer's system, e.g. "30;354x420y...".


  7. Click Okay.

    The Down-Scaling dialog box appears.

  8. Click Okay.

    The aggregated grid layer is resampled to match the original point grid layer's dimensions.