Create a digital elevation model (DEM) from a LAS file with FME

I am quite impressed with the capabilities of the surface generation and raster manipulation transformers in FME. Coupled with the rich set of vector feature transformers, you could almost put together a terrain modelling and image processing package. I tried out generating a terrain model in ESRI ArcGrid ASCII format from a LiDAR LAS file and I am quite happy with the results. To do a simple transformation, just do the following:

1. Start up FME Workbench and create an blank workspace. Choose Source Data | Add Source Dataset.
   
   The Add Source Dataset dialog box appears.
   
   
2. Choose ASPRS Lidar Data Exchange Format (LAS) in the Format drop down list. In the Dataset field, click the Browse button and select a LAS file e.g. 288_4320.LAS. Click OK.
   
   The Select Feature Types dialog box appears.
   
   
3. Toggle off header and variable_length_header (optional). Click OK.
   
   The source dataset is added to the workbench.
   
4. Select Destination Data | Add Dataset.
   
   The Add Destination Dataset dialog box appears.
   
   
5. In the Format drop down list, choose ESRI ASCII Grid. In the Dataset field, click Browse and choose an output folder e.g. C:\data. Click OK to close the Add Destination Dataset dialog box.
   
   A message appears to confirm the addition of a new feature type.
   
6. Click Yes and OK to add a new Feature Type to the workbench.
   
   The destination dataset is added to the workbench.
   
7. Drag and drop the RasterDEMGenerator transformer and connect the source and destination datasets with the transformer as shown below.
   
   
   
8. Open up the RasterDEMGenerator Parameters. In the Output DEM X Cell Spacing and Output DEM Y Cell Spacing fields, type in a cell spacing value e.g. 1 as shown below. Click OK.
   
   
   
9. Run the translation.
   
   The ESRI ASCII ArcGrid DEM file is generated and can be displayed in a terrain model viewer.
   

Create LiDAR intensity images using GeoMedia Grid

I tried to create LiDAR intensity images using GeoMedia Grid. It can be done interactively but I found the workflow to be quite tedious especially if you have a large number of tiles to handle. If you have a large number of intensity images to create, it would be more practical to automate the process either by developing a custom GeoMedia Grid program or use alternative software. To create the intensity images, here are the steps to do the job interactively.

Create a Study Area

1. Start up GeoMedia and open up a geo-workspace and connect to a read/write warehouse with the appropriate coordinate system.
2. Select Grid | Study Area | Define New.
   
   
3. In the Coordinate Precision Key-in tool bar, type in the precise coordinate of one corner of the tile. Press RETURN.
   
   
4. Type in the precise coordinate of the diagonally opposite corner of the tile. Press RETURN.
   
   The Define New Study Area dialog box appears.
   
5. In the Study Area name field, type in the desired name. In the Cell resolution field, type in the output resolution e.g. 1 and set the units e.g. meters. Press OK.
   
   The study area is created.

Import a LAS tile file

1. Follow the steps in my previous post to import a LAS file into the study area created previously.

Triangulate the tile by intensity values

1. Select Grid | Study Area | List. Choose the previously created study area and imported LAS tile layer.
   
   The selected study area and layer are active.
2. Select Grid | Interpolation | Triangulation.
   
   
   The Triangulation dialog box appears.
   
3. Click Add.
   
   The Add Source Feature dialog box appears.
4. In the Feature field, select the imported LAS tile layer in the current study area. Toggle on the Attribute Z-value. In the Attribute drop down list, choose Intensity. In the Units drop down list, choose None.
   
   The Add Source Feature dialog box may look like this at this point.
   
5. Click OK.
   
   The Triangulation dialog box is updated.
   
6. In the Output units field, select None. In the Result layer Name field, type in a meaningful name. Toggle off Place result in map window. Click OK.
   
   The result layer is created and added to the active study area.
   

Define the intensity gray scale range
By default, GeoMedia Grid will assign a grey scale color range to the result layer. Normally, the void or null values will be assigned as white. I prefer to color the voids as black. If you have only one tile, you might as well use the default. But if you have a number of tiles, then in order to ensure that all the tiles look the same, you have to create a grey scale xml file and reuse that xml file for all the tiles.

1. In the Study Area List dialog box, right click on a layer.
   
   A pop up menu appears.
   
   
2. Choose View Legend.
   
   The View Legend dialog box appears.
   
   
3. Double click the box on the left of the VOID label as shown above.
   
   The color picker appears.
   
   
4. Choose black.
   
   The View Legend(s) dialog box is updated.
   
   
5. Click Apply.
   
   The color range is applied to the layer.
   
   Note: Click Save As to save the color range into an xml file. To reuse the xml file, simply click Load and choose the previously saved xml file.
   
6. Click Close.

Create the intensity image as a GeoTiff file

1. In the Study Area List dialog box, click on the newly created layer.
2. Select Grid | Layer | Export to File.
   
   
   
   The Export to File dialog box appears.
   
3. Click OK.
   
   The Save As dialog box appears.
   
4. In the File name field, type in a name. In the Save as type drop down list, choose TIFF. Click Save.
   
   The Export to GeoTiff dialog box appears.
   
5. Choose RGB or Grayscale. Click OK.
   
   The intensity GeoTiff image file is created.
   

Create the coverage shape file of flightline LiDAR data

It is useful to be able to visualize the coverage of the LiDAR data captured during an airborne flight especially for large survey areas. I found this tool lasboundary.exe from the free lastools toolkit to be quite useful to generate the coverage outline from the captured flight LiDAR data. The lasboundary executable generates a linear shape file from the LiDAR data points of the source LAS file by fitting a concave hull of all the points. I would have like it to generate the coverage as a polygon instead of a linear geometry but if I have the time, I might modify the executable to do just that. 

Using the lasboundary executable is easy as described below:

1. In the Windows Command Prompt, to show the help just type in the following:
   
   C:\> lasboundary -h
   
   
   
2. To generate the coverage, type in the following:
   
   C:\> lasboundary -i input.las -o coverage.shp
   
   
   The coverage shape file is created and can be displayed.
   

Performance of GeoMedia Grid's Import LAS File to Feature Class command

I have tried GeoMedia Grid 6.1.1's Import LAS File to Feature Class command and I found it rather painful to use especially when processing hundreds of gigabytes worth of LAS files. Converting binary LiDAR points in LAS format files to read/write database feature class records is an expensive and time consuming process. You might be better off  filtering the LAS files to just what you need before using the Import LAS File to Feature Class command. For example, if I want to generate a bare earth DEM, I might only require LiDAR ground points; there is no point keeping the vegetation points or the building points.

I found the free command line programs LAStools on the Internet. It is quite useful to manipulate LAS files and it comes with the las2las executable. I use it to remove the LiDAR points that I do not want and retain the points that I want.

To see how to use this las2las executable, type in the following in a Windows Command Prompt as shown below:
C:\> las2las -h

To remove all points except the LiDAR ground points from a LAS file, you can type in the following in the Windows Command Prompt as shown below.

C:\> las2las -i input.las -o output.las -keep_class 2


where the value 2 is the standar LiDAR classification for ground points. 

Once I have filtered out the unwanted points, I then run GeoMedia Grid's Import LAS File to Feature Class command, which will run a little faster - because there would be less points to convert than before.

Importing LAS files with GeoMedia Grid 6.1.1

The latest GeoMedia Grid version 6.1.1 comes with a function to import LiDAR data in LAS (*.las) format. With this new function, you will be able to convert LiDAR xyzi points into point geometry features in a read/write warehouse. Once converted into point geometry features, you will then be able to query and manipulate the LiDAR points like any typical GeoMedia feature classes. If you want to generate a digital elevation model (DEM) from the LiDAR points, then you would use the typical GeoMedia Grid workflow to rasterize the points, then interpolate the rasterized points into a surface. I tried out the new command on a small LAS file.

1. Run GeoMedia and create or connect to a read/write warehouse.
   
   
2. Select Grid | Utilities | LAS Files to Feature Classes as shown below.
   
   
   
   The LAS Files to Feature Class dialog box appears.
   
   
   
3. Click Add. Browse and select a LAS file e.g. B3.LAS and click Open.
   
   The LAS Files to Feature Class dialog box is updated with the selected file.
   
   
   
4. In the Coordinate System file name field, click Browse. Select an appropriate coordinate system file (*.csf) for the LAS file e.g. StatePlane83TexasCentral.csf and click Open.
   
   The LAS Files to Feature Class is updated with the selected CSF file.
   
   
5. If the Connection is not the right destination, then choose the correct read/write connection from the drop down list.
   
   
6. If necessary, type in the new feature class name in the Feature Class field.
   
   
7. If you want to expose the LAS file name as an attribute in the new feature class, then toggle on Include the name as attribute.
   
   
8. Click OK.
   
   The status field shows the progress.
   
   
   At the end of the conversion, the summary messages appear.
   
   
9. Click Close.
   
   The LAS file is converted into a GeoMedia feature class.

If you open up a Data Window and display the attribute records of the newly converted LiDAR points, the following fields can be seen (including Filename if the Include Filename as attribute toggle was on during the conversion):

Import ArcGrid ASCII files into GeoMedia Grid

ArcGrid ASCII files are convenient to pass data between application software. Before using terrain data in ArcGrid ASCII format files (*.asc, *.agr), they must be imported into GeoMedia Grid's native format raster files (*.mfm) and placed as layers under what is termed a Study Area. I have no idea what .mfm stands for and I could not find any documentation details on the extension name. 

To import an ArcGrid ASCII file, the following steps can be done. 

1. Run GeoMedia and open up or create a GeoWorkspace with the desired coordinate system and a read/write warehouse connection.
   
   
2. Select Grid | Study Area | Import File(s).
   
   The Import File(s) dialog box appears.
   
   
   
3. Click Browse and select an ArcGrid ASCII file e.g. B3.asc.
   
   
   
   
4. Click Open.
   
   The File Type(s) Found dialog box may appear.
   
   
   
5. Choose ArcGrid ASCII. Click OK.
   
   
6. Click the Coordinate system file Browse button. Choose the appropriate coordinate system file (*.csf) and click Open.
   
   The Import File(s) dialog box is updated with the selected.
   
   
7. Click Next.
   
   The selected file is scanned and then the ASC Input data dialog box appears.
   
   
8. Click the X Coordinate and Y Coordinate units drop down list and select the appropriate horizontal units e.g. feet.
   
   
9. In the Cell resolution text field, type in the correct cell spacing e.g. 2.
   
   
10. In the Cell resolution units drop down list, select the correct vertical units e.g. feet.
    
    At this point, the ASC Input data dialog box may look like this.
    
    
11. Click OK.
    
    The Import File(s) dialog box appears.
    
    
    Note: By default the Study Area name will be Study Area 1. Right click on the name to rename it to something else if necessary.
    
    
12. Click Next.
    
    The Import File(s) dialog box is updated with statistical information about the file.
    
    
13. Click Finish.
    
    The ArcGrid ASCII file is imported into GeoMedia Grid and is ready to be used.

GeoMedia Export to Shape File Annoyances

GeoMedia has been able to export features and queries out into ESRI Shape Files for the longest time. Users of ArcGis can be a little miffed if they find out they have to read in shape files created from GeoMedia because the shape files lack .prj files, which contain the coordinate projection information associated with the exported files.

The shape files created by GeoMedia will only have the following files created: .dbf, .shp, and .shx, as shown in the screen shot below. GeoMedia will also create a file with the extension .csf that contains the coordinate system of the exported shape file. However, this is a proprietary file and can only be read using Intergraph software.

When these GeoMedia created shape files are read in ArcMap with the Layer | Add Data command as shown below, ArcGis will not be able to determine the correct coordinate system and the software will pop up a warning message: 

The shape files can be displayed in ArcMap but the coordinates are plain XY coordinates; the units and projection are unknown as highlighted in red below. The shape layers cannot be reprojected to another projection system. 

In order to resolve this problem, ArcGis users will have to create the .prj file by copying over some existing .prj file from perhaps the ArcGis product folder or through some other means. This task can get very tedious if there are a large number of shape files. 

If the shape files have associated .prj files, then ArcMap can read the shape files including the coordinate system defined in the .prj file successfully, as shown in the figure below. Note the projection information highlighted in red. 

I have heard from Intergraph that the next version of GeoMedia will be able to generate the .prj file when exporting out the shape files. I hope this will come to pass.