DISPLAYING ALTITUDE ON MAP ON MOUSEMOVE OR CLICK

We already know that there are DEM files available online for the whole world. This files have the AMSL values for every x meters depending on its spatial resolution. Retrieving these values on a web page requires little coding. In this blog, I will be telling how to achieve this using GDAL library, ASTER GDEM DEM files and Java.

Following pre-requisites are to be fulfilled before proceeding to the coding part.

• Set up the GDAL environment on you PC/Server depending on its architecture (x64/x86).
• Copy and paste following four/one dll file(s) from C:\gdal\bin\gdal\java to your jdk bin folder.
  • gdalconstjni.dll
  • gdaljni.dll
  • ogrjni.dll
  • osrjni.dll
  OR
  • gdalalljni.dll (in higher versions of gdal v2.3.0 and above)

<VRTDataset rasterXSize="64801" rasterYSize="43201">
  <SRS>GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],
  AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0],
  UNIT["degree",0.0174532925199433],AUTHORITY["EPSG","4326"]]</SRS>
  <GeoTransform> 6.9999861111111116e+001, 2.7777777777777778e-004, 0.0000000000000000e+000, 
  3.9000138888888891e+001, 0.0000000000000000e+000,-2.7777777777777778e-004</GeoTransform>
  <VRTRasterBand dataType="Int16" band="1">
    <ColorInterp>Gray</ColorInterp>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">ASTGTM2_N27E086_dem.tif</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="3601" RasterYSize="3601" DataType="Int16" BlockXSize="3601" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="3601" ySize="3601" />
      <DstRect xOff="57600" yOff="39600" xSize="3601" ySize="3601" />
    </SimpleSource>
    <SimpleSource>
      <SourceFilename relativeToVRT="1">ASTGTM2_N27E087_dem.tif</SourceFilename>
      <SourceBand>1</SourceBand>
      <SourceProperties RasterXSize="3601" RasterYSize="3601" DataType="Int16" BlockXSize="3601" BlockYSize="1" />
      <SrcRect xOff="0" yOff="0" xSize="3601" ySize="3601" />
      <DstRect xOff="61200" yOff="39600" xSize="3601" ySize="3601" />
    </SimpleSource>
	.
	.
	.
  </VRTRasterBand>
</VRTDataset>

Open a new Web application project in your IDE say amslOnMap. Create a servlet say CalcElevation.java in your xyz package. We will be using only the doget method here to handle the incoming request for AMSL at asked location. We will create a Java Class file say GDALEle.java which will read the vrt file created above and retrieve the AMSL value at any given coordinate.

You need to add gdal.jar file in the JAVA BUILD PATH Library location.

• For Netbeans Right-Click project --> Properties --> Libraries --> Compile tab
• For Eclipse Right-Click project --> Properties --> Java Build Path --> Libraries tab.

Copy the following Code in your GDALEle.java file

package xyz;

import org.gdal.gdal.*;
import org.gdal.gdalconst.gdalconstConstants;

public class GDALEle {

	private final static String vrtLocation = "D:\\GDEM\\GDEM.vrt"; //VRT file location
	private static Dataset hDataset = null;		//static dataset allows hDataset to be read only once.
	
	private static void regDataset(){
	        if(hDataset==null){
	            gdal.AllRegister();   	//Registers GDAL
	            hDataset = gdal.Open(GDEM_location, gdalconstConstants.GA_ReadOnly);
	        }
    	}
    
	public static String getElevation(double lat, double lon){
		regDataset(); // This static method opens the DEM location only once.
		String elevationval = null;	//This is return string
		int iBand; 
		Band hBand;
		int iPixel;
		int iLine;
		double[] adfGeoTransform = new double[6]; 
		double[] value = new double[1]
        
		hDataset = gdal.Open(vrtLocation, gdalconstConstants.GA_ReadOnly);          
		// Opens the VRT file as Dataset object
		if (hDataset == null)
		{
			return "No DEM data";        
			//Error if no VRT file available
		}
		int X = hDataset.getRasterXSize();        
		//number of pixels along X axis
		int Y = hDataset.getRasterYSize();        
		// number of pixels along Y axis
		iBand = hDataset.getRasterCount();        
		// Reads the number of bands in vrt file. This value must be 1.
		
		if( iBand == 1 ){
			hBand = hDataset.GetRasterBand(iBand);       
			//Reads the raster band as Band object
			hDataset.GetGeoTransform(adfGeoTransform);       
			// Reads the Geotransform of the vrt. It has six double values
			
			if (adfGeoTransform[2]==0 && adfGeoTransform[4]==0){       
			// 3rd and 5th values of Geotransform must be zero
				iPixel = (int) Math.floor((lon - adfGeoTransform[0]) / adfGeoTransform[1]);      
				// get the pixel value for the longitude in vrt 
				iLine = (int) Math.floor((lat - adfGeoTransform[3]) / adfGeoTransform[5] );      
				// get the pixel value for the latitude in vrt 
				
				if (iPixel<0 || iLine<0 || iPixel>X || iLine>Y){      
				// if pixel values are 0 or outside vrt
					return "Point out of DEM"; 
				}else{    
					try{
						hBand.ReadRaster(iPixel, iLine, 1, 1, value);       
						// else read the pixel value of the vrt
						elevationval = String.valueOf((int)value[0]);
						hBand.delete();
					}catch(Exception e){
						elevationval = "Error: "+e.toString();
					}
				}
		
			}else{elevationval = "Invalid DEM Data";}
		}else{
			elevationval = "DEM is multi-band";
		}
		hDataset.delete();
		return elevationval ;
	}
}

You will notice the hDataset has been declared as static and null. The value of hDataset is initialised using gdal.Open method. This method of reading VRT DEM file takes little time. If this method is called for each coordinate received from map on mousemove, this time would become large. Once hDataset is initialised, then this variable belongs to class (static variable) and need not be initialised again. This technique allows the altitude to be displayed on map without showing processing gif file.

Copy the following codes in doGet method of CalcElevation.java servlet.

	PrintWriter out = response.getWriter();
	if(request.getParameter("lat")!=null && request.getParameter("lon")){
		String lat = request.getParameter("lat");			
		// retrieve the latitude
		String lon = request.getParameter("lon");			
		// retrieve the longitude
		double[] coordinates = new double[2];
		try{
			coordinates[0] = Double.valueOf(lon);			
			// Check longitude as double and save as coordinates[0]
			coordinates[1] = Double.valueOf(lat);			
			// Check latitude as double and save as coordinates[1]
			out.println(GDALEle.getElevation(coordinates));		
			// Call getElevation method of GDALEle class directly since it is static and print the return value
		}catch(Exception e){
			out.println("Invalid coordinates: "+e.toString());	
			//Error if lat or lon are not numbers
		}
	}else{
		out.println("Missing Parameters lat and lon. Try again");	
		// Error if lat or lon are missing
	}

Run this web application in Java supported server like Tomcat, Glassfish, JBoss etc., and pass the values of latitude and longitude as parameter in url. For example if your server runs on port 8080 on http on localhost then the url will be

http://localhost:8080/amslOnMap/CalcElevation?lat=23.5687&lon=45.9862

This will print the value of AMSL at latitude 23.5687 and longitude 45.9862 provided this coordinate is falling in your VRT.

You call pass parameters through Ajax calls using EventListener either click or mousemove on your GIS web clients. Remember to check for cross-origin CORS if web client and the above codes run on different servers. For example:

• AJAX call through jQuery:
  

  
  var xhr = null;
  function calcElevation(lat, lon)
  {
   if( xhr != null ) {
    xhr.abort(); //abort pending request in case a fresh request comes.
    xhr = null;
   }
   xhr = $.ajax({
    type: "GET",
    url: "http://localhost:8080/amslOnMap/CalcElevation?lat="+lat+"&lon="+lon,
    success: function(msg) {
     if (msg!==''){
      $("#elevation").html = "Elevation (m): "+msg; 
     }else{
      $("#elevation").html = "Elevation (m): <img src='images/processing.gif' width='12px;'>"; 
     }
    }
   });
  }
  

• Openlayers v3+: You can use 'pointermove' or 'click' events. Example for pointermove:
  

  
  map.on('pointermove', function(evt) {
  	if (evt.dragging) {
  		return;
  	}
  	var coordinate = map.getEventCoordinate(evt.originalEvent);
  	console.log ("Lat: "+coordinate[1]+" Lon: "+coordinate[0]);
  });
  

• Leaflet: You can use 'mousemove' or 'click' events. Example for click:
  

  
  map.on('click', function(e) {
  	console.log("Lat: "+e.latlng.lat+" Lon: "+e.latlng.lng);
  });
  

Remember that while working with Java code along with GDAL on web applications, you need to restart the server every time you make some changes with Java codes else it would show UnsatisfiedLinkError from GDAL.

While learning different algorithms, I have learnt that this is not the best way to calculate AMSL on web page because opening the VRT file for each request is time consuming.