Sunday, August 16, 2015

Setting Up GeoDjango II

In my previous post on this topic, I ended up stuck with a database not existing error, which my good friend Kevin Cole very politely (only implying that I'm an idiot, while refraining from directly saying so ;-) pointed out that the documentation I'm using contains the instructions I need to create the database:
(env)$ createdb -T template_postgis geodjango

Unfortunately, running that gave me the following error:

createdb: database creation failed: ERROR:  template database "template_postgis" does not exist
Let me try a modified version (modified because my user doesn't have the privileges needed to create extensions) of the steps laid out here:
(env)$ createdb geodjango
(env)$ sudo su - postgres
$ psql geodjango
psql (9.3.9)
Type "help" for help.

geodjango=# CREATE EXTENSION postgis;
geodjango=# \q
$ exit
Now let me resume where I left off before the error:
(env)$ python sqlmigrate world 0001
CREATE TABLE "world_worldborder" ("id" serial NOT NULL PRIMARY KEY, "name" varchar(50) NOT NULL, "area" integer NOT NULL, "pop2005" integer NOT NULL, "fips" varchar(2) NOT NULL, "iso2" varchar(2) NOT NULL, "iso3" varchar(3) NOT NULL, "un" integer NOT NULL, "region" integer NOT NULL, "subregion" integer NOT NULL, "lon" double precision NOT NULL, "lat" double precision NOT NULL, "mpoly" geometry(MULTIPOLYGON,4326) NOT NULL);
CREATE INDEX "world_worldborder_mpoly_id" ON "world_worldborder" USING GIST ("mpoly" );

(env)$ python migrate
Operations to perform:
  Synchronize unmigrated apps: gis, messages, staticfiles
  Apply all migrations: world, admin, contenttypes, auth, sessions
Synchronizing apps without migrations:
  Creating tables...
    Running deferred SQL...
  Installing custom SQL...
Running migrations:
  Rendering model states... DONE
  Applying contenttypes.0001_initial... OK
  Applying auth.0001_initial... OK
  Applying admin.0001_initial... OK
  Applying contenttypes.0002_remove_content_type_name... OK
  Applying auth.0002_alter_permission_name_max_length... OK
  Applying auth.0003_alter_user_email_max_length... OK
  Applying auth.0004_alter_user_username_opts... OK
  Applying auth.0005_alter_user_last_login_null... OK
  Applying auth.0006_require_contenttypes_0002... OK
  Applying sessions.0001_initial... OK
  Applying world.0001_initial... OK
Progress!  Let me keep going and see if my good fortune holds:
(env)$ python shell
Python 3.4.0 (default, Jun 19 2015, 14:20:21)
[GCC 4.8.2] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import os
>>> import world
>>> world_shp = os.path.abspath(os.path.join(os.path.dirname(world.__file__),
...                             'data/TM_WORLD_BORDERS-0.3.shp'))
>>> from django.contrib.gis.gdal import DataSource
>>> ds = DataSource(world_shp)
>>> print(ds)
/home/[user/geodjango/geodjango/world/data/TM_WORLD_BORDERS-0.3.shp (ESRI Shapefile)
>>> print(len(ds))
>>> lyr = ds[0]
>>> print(lyr)
>>> print(lyr.geom_type)
>>> print(len(lyr))
The tutorial continues with several other interactive examples showing how to use GeoDjango's pythonic interface to the GDAL library.  I began experimenting with python's GDAL wrapper back in April as part of the Introduction to GIS Programming and Algorithms course I took at George Mason University.  I documented the installation of these tools in a post at that time.  The ability to "play" with data at run time is one of the many things I love about Python, and this tutorial, like most python tutorials, is making good use of that powerful pedagogical feature of the language.  There is no need for me to recount the other examples here, however, so I'll skip over them.

The next step in the tutorial is to create a file in the world app named that contains the following:
import os
from django.contrib.gis.utils import LayerMapping
from .models import WorldBorder

world_mapping = {
    'fips' : 'FIPS',
    'iso2' : 'ISO2',
    'iso3' : 'ISO3',
    'un' : 'UN',
    'name' : 'NAME',
    'area' : 'AREA',
    'pop2005' : 'POP2005',
    'region' : 'REGION',
    'subregion' : 'SUBREGION',
    'lon' : 'LON',
    'lat' : 'LAT',
    'mpoly' : 'MULTIPOLYGON',

world_shp = os.path.abspath(os.path.join(os.path.dirname(__file__), 'data/TM_WORLD_BORDERS-0.3.shp'))

def run(verbose=True):
    lm = LayerMapping(WorldBorder, world_shp, world_mapping,
                      transform=False, encoding='iso-8859-1'), verbose=verbose)
Note: The tutorial lists: from models import WorldBorder, which will cause an import error.  models needs to be .models for this to work.

After making that change, I was able to:
(env) python shell
>>> from world import load
and watch as the countries of the world were loaded into the database.

Creating the github repo

Now would a good time to create a github repo.  First, I'll create a .gitignore file inside the top level geodjango directory (where is located) that lists the things I don't want in the repository:
This will tell git not to include the virtual environment, the python byte code files, and any vim swap files.  Next install git and initialize the repository:
(env)$ sudo aptititude install git
(env)$ git init
Now I'll check to see what a git add . would add:
(env)$ git add -n .
Since it looked good, I'll add do it, after configuring my git email and user:
(env)$ git config --global [github email address]
(env)$ git config --global [github name]
(env)$ git add .
(env)$ git commit -a
Initial commit.
Now push it to github (after adding an ssh key and creating a learn_geodjango project on github):
git remote add origin[user]/learn_geodjango.git
git push -u origin master
With the github repo now created, I'll continue the tutorial in a future post.

Thursday, August 6, 2015

Setting Up GeoDjango I

According to the documentation website,
"GeoDjango intends to be a world-class geographic Web framework. Its goal is to make it as easy as possible to build GIS Web applications and harness the power of spatially enabled data."
Since my goal is to become a free software GIS application developer, and since Python and Django are two of the core technologies I hope to utilize, GeoDjango seems like a no-brainer as something I should learn.

In this post I'm going to document the beginning process of setting up a GeoDjango server on a local VirtualBox VM.  In a later post, I'll look into installing it on the WebFaction application hosting service. Throughout this process, I'm going to modify the steps described in the tutorial to facility creating a github repo and migrating the app to WebFaction.

Since GeoDjango comes with Django, I'm going to start by creating a virtualenv on the same virtual machine I used to setup a django CMS virtualenv, the process for which I described in a previous post, Installing Django CMS on Ubuntu 14.04.  I'll install this new virtualenv right alongside the other one.  I'll be following the installation instructions from the GeoDjango Installation page, together with blog posts I've made previously documenting specific steps to install in my Ubuntu 14.04 VM. Installing GeoDjango requires installation of:
  1. Django
  2. Spatial database
I'll tackle each in turn.

Installing Django

$ mkdir geodjango
$ cd geodjango
$ virtualenv  env
$ source env/bin/activate
(env)$ pip install django
(env)$ django-admin --version
(env)$ deactivate

Installing a Spatial Database

According to the documentation, "PostGIS is recommended, because it is the most mature and feature-rich open source spatial database." It also recommends that Ubuntu installation use packages.  Since I'm new to this, I'll follow the recommendation.

$ sudo aptitude install postgresql-9.3 postgresql-9.3-postgis-2.1 postgresql-server-dev-9.3 python3-psycopg2
Next setup a database user who can create databases:
$ sudo su - postgres
$ createuser --createdb [user]
$ exit
Now follow the tutorial to see if it works (note: the command in red has been modified from tutorial with the aim of keeping the virtual environment directory (env) inside the project directory):
$ cd geodjango
$ source env/bin/activate
(env)$ django-admin startproject geodjango .
(env)$ python startapp world
(env)$ vi geodjango/
Change the DATABASES section to match the following:
    'default': {
        'ENGINE': 'django.contrib.gis.db.backends.postgis',
        'NAME': 'geodjango',
        'USER': '[user]',
Also add the last two items to INSTALLED_APPS so that it looks like this:
For the next step, we will need some gdal packages (and I'll grab unzip while I'm at it:
(env)$ sudo aptitude install gdal-bin python3-gdal unzip
(env)$ mkdir world/data
(env)$ cd world/data
(env)$ wget
(env)$ unzip
(env)$ cd ../..
(env)$ ogrinfo world/data/TM_WORLD_BORDERS-0.3.shp
INFO: Open of `world/data/TM_WORLD_BORDERS-0.3.shp'
      using driver `ESRI Shapefile' successful.
1: TM_WORLD_BORDERS-0.3 (Polygon)

The tutorial includes another ogrinfo example command that provides more detail about the world borders shape file, which I'll skip in the interest of space.  I will include the excellent short description of each of the component files in the shapefile:
Holds the vector data for the world borders geometries.
Spatial index file for geometries stored in the .shp.
Database file for holding non-geometric attribute data (e.g., integer and character fields).
Contains the spatial reference information for the geographic data stored in the shapefile.
Next edit the file so that it looks like this:

from django.contrib.gis.db import models

class WorldBorder(models.Model):
    name = models.CharField(max_length=50)
    area = models.IntegerField()
    pop2005 = models.IntegerField('Population 2005')
    fips = models.CharField('FIPS Code', max_length=2)
    iso2 = models.CharField('2 Digit ISO', max_length=2)
    iso3 = models.CharField('3 Digit ISO', max_length=3)
    un = models.IntegerField('United Nations Code')
    region = models.IntegerField('Region Code')
    subregion = models.IntegerField('Sub-Region Code')
    lon = models.FloatField()
    lat = models.FloatField()
    mpoly = models.MultiPolygonField()
    objects = models.GeoManager()

    def __str__(self):              # __unicode__ on Python 2
The default spatial reference system is WGS84.  New to me from the documentation is the existence of Open Geospatial Consortium (OGC) Spatial Reference System Identifier (SRID), which in the case of WGS84 is 4326.

The next step in the instructions were to run:
(env)$ python makemigrations
This didn't work for me, giving me the following error:
ImportError: No module named 'psycopg2'
I fixed this with:
(env)$ pip install psycopg2
(env)$ python makemigrations
Migrations for 'world':
    - Create model WorldBorder
But when I ran the next command:
(env)$ python sqlmigrate world 0001
I got a database error:
psycopg2.OperationalError: FATAL:  database "geodjango" does not exist
Anticipating that I might need knowledge of PostgreSQL databased administration, I began looking into that in an earlier post. I better look into that further before returning to the present task.

Monday, August 3, 2015

Tools for Visualizing Lidar Data I


A friend of mine recently sent me a link to an article on titled   Manipulating data in 3D with LidarViewer.  His timing couldn't have been better, since LidarViewer is a free software tool for visualizing Lidar data, and is thus just the kind of thing I need for my Summer project.  Even better, the downloads page starts with an Ubuntu PPA, so installation should not be a nightmare.

I'll do my installation on a VirtualBox VM running Lubuntu 14.04.  I like to use VMs whenever I am trying new software that is not part of the standard Ubuntu software repositories, since this keeps my host machine stable, while letting me experiment without fear.

Here is what I did:
$ sudo apt-add-repository ppa:keckcaves/ppa
$ sudo aptitude update
The installation page doesn't list the packages included in the PPA, but the Launchpad page for the repo does.  There is a package named lidarviewer, so I'll install that:
$ sudo aptitude install lidarviewer
Running lidarviewer at the command prompt after installation completed gave me a "command not found", so I checked to see what the package had installed:
$ dpkg -L lidarviewer

I tried:
$ LidarViewer
Caught exception LidarViewer::LidarViewer: No octree file name provided
A quick search on the exception led me to the Lidar Viewer Manual.  Since I already installed the application from the Ubuntu PPA (Ubuntu rocks!), I can skip most of the installation instruction section. In the MacOS instructions, however, I found sample data for testing the application.  Downloading and unzipping the sample data, I changed into the LidarViewerExamples directory and ran the following command and got the following error:
$ LidarViewer PtArena.lidar
Cache sizes: 4672 memory nodes, 1170 GPU nodes
libGL error: pci id for fd 12: 80ee:beef, driver (null)
OpenGL Warning: Failed to connect to host. Make sure 3D acceleration is enabled for this VM.
libGL error: core dri or dri2 extension not found
libGL error: failed to load driver: vboxvideo
I'm using a VirtualBox VM, and this message is telling me to enable 3D acceleration.  After enabling 3D acceleration on the VM, the application reported a long list of OpenGL errors.  I found this bug ticket showing I'm not the only one with the issue.

I added the same PPA to a laptop running Ubuntu 14.04 (thus loosing the safety of the virtual machine) and installed both the lidarviewer and  vrui-examples packages, after which I could run the examples on the laptop without error.


Next I'm going to install another set of tools for visualizing and processing Lidar data, lidar2dems.  Installation instructions are found here, and contain a number of utilities such as LAStools, which I'll need to uncompress the LAZ files that are on the Virginia Lidar website.

Installation of the lidar2dems software is made easy by an installation script,  After downloading the script, run:
$ chmod +x
$ sudo ./
The script took almost two hours to complete on my VirtualBox VM, but it completed without error.

It did not, however, install many of the LAStools utilities, especially laszip, as I had hoped.


To get laszip, I went to and downloaded  Then:
$ unzip
 which created a LAStools directory with several subdirectories, including a bin subdirectory that contained windows .exe binaries and _README.txt text files for many LASzip utilities, including laszip.exe.  Next I:
$ cd LAStools
$ make
This created the following Linux binaries in the bin directory:
las2las  lasdiff   lasinfo   lasprecision  txt2las
las2txt  lasindex  lasmerge  laszip
To test if laszip works, I grabbed a LAZ file for downtown Leesburg:
$ wget
This got me the file 18STJ7733.laz. Then I ran:
$ ls -l
total 16624
-rw-rw-r-- 1 user user 17019618 Aug 17  2012 18STJ7733.laz
$ laszip 18STJ7733.laz
$ ls -l
total 149256
-rw-rw-r-- 1 jelkner jelkner 135810762 Aug  3 12:05 18STJ7733.las
-rw-rw-r-- 1 jelkner jelkner  17019618 Aug 17  2012 18STJ7733.laz
So it appears to have uncompressed the LAZ file into a much larger (almost 8x larger) LAS file.

Given that it is free software with an LGPL license, I don't understand why someone in the Open Source GIS community hasn't made this installation much simpler on Ubuntu yet.  For now, I've made a modest contribution toward that goal by creating a page with the Ubuntu binaries on my Open Book Project site:

Friday, July 31, 2015

Manual of Airborne Topographic Lidar - Chapter 2 Acronym Glossary

One of the most challenging things about reading this text is the incredible number of acronyms used.  I created the following glossary of acronyms to enable me to make sense of what I read in chapter two.

Glossary of Acronyms

  • ALS - Airborne Laser Scanning
  • AM - Amplitude modulation 
  • AOBD - Acousto-optic beam deflector
  • AOL - Airborne Oceanic Lidar
  • AOI - Angle of incidence
  • APD - Avalanche photodiode
  • ARD - Automatic rendezvous and docking 
  • ARV - Autonomous robotic vision
  • ASC - Advanced Scientific Concepts (company) 
  • ASPRS -  American Society for Photogrammetry and Remote Sensing
  • ATM - Airborne Topographic Mapper 
  • ATR - Automatic target recognition 
  • BATC - Ball Aerospace and Technologies Corporation
  • CCD - Charge-coupled device 
  • CIR - Color-infrared
  • CMOS - Complementary metal–oxide–semiconductor
  • COTS - Commercial off-the-shelf
  • CW - Continuous wave
  • DEM - Digital elevation model 
  • DESDynI - Deformation, Ecosystem Structure, and Dynamics of Ice
  • DHS SBIR - Department of Homeland Security Small Business Innovation Research program
  • DSM - Digital surface model
  • DSP - Digital signal processing
  • DTM - Digital terrain model
  • EDME - Electronic Distance Measuring Equipment
  • ESFL - Electronically Steerable Flash Lidar
  • FDC - Frequency-to-distance conversion
  • FFPA - Flash-focal-plane-array
  • FM - Frequency modulation
  • FOR - Field of regard
  • FOV - Field of view
  • FPGA - Field programmable gate array
  • FW - Full waveform
  • FWHM - Full width half maximum
  • GmAPD - Geiger Mode Avalanche Photodiode
  • GPS - Global positing system
  • HLDA - Autonomous hazard detection and avoidance 
  • ICP - Iterative closest point
  • IFOV - Instantaneous field of view 
  • IIP - Instrument Incubator Program
  • IMU - Inertial measurement unit
  • InGaAs - Indium-gallium-arsenide
  • INS - Inertial navigation system 
  • ISD - Integrated scanner, detector
  • ISDT - Integrated scanner, detector, telescope
  • ISR - Intelligence, Surveillance, and Reconnaissance
  • ISS - International space station
  • LAPF - Laser Altimeter Processing Facility 
  • LAS -  Laser (file format) developed by ASPRS
  • LA-WASM -  large-aperture wide-angle scanning mirror
  • LIF - Laser induced florescence 
  • LINS - Laser Inertial Measurement system
  • LRF - Laser range finder
  • LVIS - Laser Vegetation Imaging Sensor
  • MASER - Microwave amplification by stimulated emission of radiation
  • MIR - Mid range infrared
  • NASA - National Aeronautics and Space Administration
  • NAVOCEANO - Naval Oceanographic Office
  • NOAA - National Oceanic and Atmospheric Administration
  • OMCMFNC - Organic Mine Countermeasures Future Naval Capabilities
  • PRF - Pulse repetition frequency
  • PROXOPS - Proximity operations
  • PRR - Pulse repetition rate 
  • RAMS - Real-time Aerial Mapping System
  • RASCAL - RAster SCanning Airborne Laser Altimeter 
  • RF - Radio frequency
  • RGB - Red-green-blue
  • ROAR - Rapid Overt Airborne Reconnaissance
  • ROIC - Readout integrated circuit
  • RVS - Raytheon Vision Systems (company)
  • SLICER - Scanning Lidar Imager of Canopies by Echo Recovery
  • SNR - Signal-to-noise ratio 
  • SWaP - Size, weight, and power
  • SWIR - Short wave infrared
  • TMFL - Topographic mapping flash lidar
  • TOF - Time-of-flight
  • TRN - Terrain-relative navigation 
  • UAV - Unmanned aerial vehicles 
  • UPS - Uninterruptabke Power Supply
  • UV - Ultraviolet
  • VIS/IR - Visible / infrared
  • VTUAV - Vertical takeoff unmanned aerial vehicle

Additional Resources

Tuesday, July 21, 2015


In a previous post I described how I setup a virtual machine and installed Node.js on it.  The goal was to help a group of Summer student interns get started with both Node.js and Firefox OS app development by working on a web application created by a previous student intern.

That plan didn't work out well, as the application in question did not provide the kind of on ramp to the technologies we wanted to learn that I had hoped it would. Since Node.js is new to both the student interns and me, we all decided to look at the tutorials on NodeSchool, particularly learnyounode.

With node, and thus npm already installed, I added the learnyounode materials to my virtual machine with:
$ sudo npm install -g learnyounode
$ learnyounode
brings up a menu:
with links to instructions for each of the 13 exercises in the tutorial.

I spent the morning working through the first 6 exercises, and created a github repository with my solutions.  I'm setting myself the additional goal of learning to write good clean JavaScript while I'm at it, so I installed JSLint:
$ sudo npm install -g jslint
And made every effort to minimize the number of warnings jslint reports on each of my solutions.  Using JSLint is almost like having an automated Douglas Crockford to look over your coding style, so I feel I'm in pretty good hands as I learn.

Learnyounode is a fabulous tutorial.  It moves you quickly and efficiently toward learning key ideas of important Node.js programming concepts by having you complete a series of well designed exercises.

I'm out of time for today, but the main concept I think I have finally begun to wrap my head around from the first 6 exercises is the concept of a callback, which I had to use in my solution to exercise 6.

Thursday, July 16, 2015

Notes from Section 2.1 of Manual of Airborne Topographic Lidar

Chapter 2 of the Manual of Airborne Topographic Lidar, titled "An Overview of ALS Technology", provides a broad introduction to the operating principles and key elements of Lidar, and then discusses the specifics of several existing Lidar systems. Since it contains 90 dense pages of information, I'll break my notes on it into several posts.

Operating Principles

  • LASER - Light Amplified by Stimulated Emission of Radiation.
  • Lidar works with light in electromagnetic spectrum from ultraviolet through infrared.
  • Wavelengths are chosen with regard to range-performance, atmospheric or water absorption, eye-safety, and reflectivity of lased materials.
  • Lasers are nearly monochromatic (they emit light within an extremely narrow band of wavelengths).
  • The emitted wavelength is a function of the material in which the light is stimulated, so wavelength vs. material tables can be produced, like this one:
  • When laser light strikes a surface, parts of the light are:
    1. transmitted
    2. absorbed
    3. reflected
  • The distribution of these 3 parts is a function of the lased surface and the wavelength of the beam (see graph below for example distributions of common earth surfaces).
  •  All lidar systems are comprised of 2 parts:
    1. transmitter emitting laser pulses
    2. optical receiver detecting backscattered pulse
    1. The illustration below shows a simplified model of the laser ranging principle:

      which is represented by the equation:
      Δt = 2R/c
      Δt is the time between when the laser pulse is emitted and when its backscattered pulse reaches the receiver. c is the speed of light, which is approximately 299700 km/s in the earth's atmosphere. Solving for the range, R, we get:
      R = Δt·c/2
      so the range is half the product of the elapsed time and the speed of light.
    2. The square wave idealization is not a very realistic representation of an actual laser pulse, which is much more closely modeled by a Gaussian wave like this:
      The text discusses in some detail how this waveform is analyzed to determine the pulse. One common way this is done is referred to as full width half maximum (FWHM), which is the pulse duration at 50% or more of its peak intensity. Two other common measurements are 1/e and 1/e2 (shown in the above illustration) of peak intensity.

    Additional Resources

      Wednesday, July 15, 2015

      Learning Leaflet 1

      Leaflet is a free software JavaScript library for creating web mapping applications.  Since the Photovoltaic Viability Map project uses it, it's time I learned how to use it as well.

      Using the following three resources:
      I created a web page with a map displaying the area around where I live:

      Which rendered in a browser produces this:

      Lessons Learned

      From The Basics page on the SWITCH2OSM site, I learned that there are two components to a web based map application:
      1. A database of 'tiles' which are rendered (usually in 256x256 pixel sections) side by side to make the map.
      2. A JavaScript API for viewing the database tiles.
      While the OpenStreetMap (OSM) tiles are open data, and are thus free to use, the serving of them is not. Given the cost of maintaining and serving 46 gigabytes of data to the whole world, this is understandable, but I never thought about it before. OSM has a usage policy for their tile servers, which also suggests alternative tile servers for their tiles.

      So to setup a web map application that uses OSM tiles, one has two choices:
      1. Download the OSM database and generate the tiles yourself.
      2. Use a third-party supplier (some charge, some don't).
      The German company Geofabrik hosts the OSM database separated by regions, countries and states (within the US) in North America here. So instead of the 46 gigabytes of the entire OSM database, I could get just the 211 megabytes for the state of Virginia here.

      A content delivery network (CDN) is a distributed system of servers across the Internet for the purpose of rapid delivery of content.  The quick start guide on the leaflet site uses the Cloudfare CDN to serve leaflet. The pv-viability-map uses MapQuest's CDN as a tile server.  For my first leaflet map, I used both of these together.

      This is JavaScript that (using leaflet) does the work of rendering the map:
      var map ='map').setView([38.8605579, -77.1166921], 15);

      L.tileLayer('http://otile{s}{z}/{x}/{y}.jpg', {
          attribution: 'Portions Courtesy NASA/JPL-Caltech and U.S. Depart. of Agriculture, Farm Service Agency. Tiles Courtesy of <a href="" target="_blank">MapQuest</a> <img src="">',
          maxZoom: 18,
          subdomains: '1234'
      The tileLayer method does the job retrieving and rendering the tiles selected by the setView method on the Map object.  The string substitution variables s, x, y, and z contain:
      Sequence of available subdomains for the CDN. Set to '1234' in this case.
      The longitude.
      The latitude.
      The zoom level.
      In this source listing I changed what had been 'sat' in the url to 'map' (colored green in the source listing for identification), which changed the rendered map to:
       I experimented with the zoom level and maxZoom settings.  Changing the zoom level to 18 and setting maxZoom to 19, the map looked like this:
      and allowed the + button to be clicked once, after which it became grayed out.  The highest usable zoom level seems to be 19 anyway, since when I tried 20 the map became solid gray. A zoom level of 1 shows the map of the entire planet.