qgis-icon_smallThanks to the work of Volker Fröhlich and other Fedora/EPEL packagers I was able to create RPM packages of QGIS 2.10 Pisa for Fedora 21, Centos 7, and Scientific Linux 7 using the great COPR platform.

Repo: https://copr.fedoraproject.org/coprs/neteler/QGIS-2.10-Pisa/

The following packages can now be installed and tested on epel-7-x86_64 (Centos 7, Scientific Linux 7, etc.), and Fedora-21-x86_64:

  • qgis 2.10.1
  • qgis-debuginfo 2.10.1
  • qgis-devel 2.10.1
  • qgis-grass 2.10.1
  • qgis-python 2.10.1
  • qgis-server 2.10.1

Installation instructions (run as “root” user or use “sudo”):

su

# EPEL7:
yum install epel-release
yum update
wget -O /etc/yum.repos.d/qgis-2-10-epel-7.repo https://copr.fedorainfracloud.org/coprs/neteler/QGIS-2.10-Pisa/repo/epel-7/neteler-QGIS-2.10-Pisa-epel-7.repo
yum update
yum install qgis qgis-grass qgis-python

# Fedora 21:
dnf copr enable neteler/QGIS-2.10-Pisa
dnf update
dnf install qgis qgis-grass qgis-python

Enjoy!

GRASS GIS and dockerSometimes, we developers get reports via mailing list that this & that would not work on whatever operating system. Now what? Should we be so kind and install the operating system in question in order to reproduce the problem? Too much work… but nowadays it has become much easier to perform such tests without having the need to install a full virtual machine – thanks to docker.

Disclaimer: I don’t know much about docker yet, so take the code below with a grain of salt!

In my case I usually work on Fedora or Scientific Linux based systems. In order to quickly (i.e. roughly 10 min of automated installation on my slow laptop) try out issues of GRASS GIS 7 on e.g., Ubuntu, I can run all my tests in docker installed on my Fedora box:

# we need to run stuff as root user
su
# Fedora 21: install docker 
yum -y docker-io

# Fedora 22: install docker
dnf -y install docker

# enable service
systemctl start docker
systemctl enable docker

Now we have a running docker environment. Since we want to exchange data (e.g. GIS data) with the docker container later, we prepare a shared directory beforehand:

# we'll later map /home/neteler/data/docker_tmp to /tmp within the docker container
mkdir /home/neteler/data/docker_tmp

Now we can start to install a Ubuntu docker image (may be “any” image, here we use “Ubuntu trusty” in our example). We will share the X11 display in order to be able to use the GUI as well:

# enable X11 forwarding
xhost +local:docker

# search for available docker images
docker search trusty

# fetch docker image from internet, establish shared directory and display redirect
# and launch the container along with a shell
docker run -v /data/docker_tmp:/tmp:rw -v /tmp/.X11-unix:/tmp/.X11-unix \
       -e uid=$(id -u) -e gid=$(id -g) -e DISPLAY=unix$DISPLAY \
       --name grass70trusty -i -t corbinu/docker-trusty /bin/bash

In almost no time we reach the command line of this minimalistic Ubuntu container which will carry the name “grass70trusty” in our case (btw: read more about Working with Docker Images):

root@8e0f233c3d68:/# 
# now we register the Ubuntu-GIS repos and get GRASS GIS 7.0
add-apt-repository ppa:ubuntugis/ubuntugis-unstable
add-apt-repository ppa:grass/grass-stable
apt-get update
apt-get install grass7

This will take a while (the remaining 9 minutes or so of the overall 10 minutes).

Since I like cursor support on the command line, I launch (again?) the bash in the container session:

root@8e0f233c3d68:/# bash
# yes, we are in Ubuntu here
root@8e0f233c3d68:/# cat /etc/issue

Now we can start to use GRASS GIS 7, even with its graphical user interface from inside the docker container:

# create a directory for our data, it is mapped to /home/neteler/data/docker_tmp/
# on the host machine 
root@8e0f233c3d68:/# mkdir /tmp/grassdata
# create a new LatLong location from EPSG code
# (or copy a location into /home/neteler/data/docker_tmp/)
root@8e0f233c3d68:/# grass70 -c epsg:4326 ~/grassdata/latlong_wgs84
# generate some data to play with
root@8e0f233c3d68:/# v.random n=30 output=random30
# start the GUI manually (since we didn't start GRASS GIS right away with it before)
root@8e0f233c3d68:/# g.gui

Indeed, the GUI comes up as expected!

GRASS GIS 7 GUI in docker container

GRASS GIS 7 GUI in docker container

You may now perform all tests, bugfixes, whatever you like and leave the GRASS GIS session as usual.
To get out of the docker session:

root@8e0f233c3d68:/# exit    # leave the extra bash shell
root@8e0f233c3d68:/# exit    # leave docker session

# disable docker connections to the X server
[root@oboe neteler]# xhost -local:docker

To restart this session later again, you will call it with the name which we have earlier assigned:

[root@oboe neteler]# docker ps -a
# ... you should see "grass70trusty" in the output in the right column

# we are lazy and automate the start a bit
[root@oboe neteler]# GRASSDOCKER_ID=`docker ps -a | grep grass70trusty | cut -d' ' -f1`
[root@oboe neteler]# echo $GRASSDOCKER_ID 
[root@oboe neteler]# xhost +local:docker
[root@oboe neteler]# docker start -a -i $GRASSDOCKER_ID

### ... and so on as described above.

Enjoy.

GRASS GIS logoAfter months of development a first release candidate of GRASS GIS 6.4.5 is now available. This is a stability release of the GRASS GIS 6 line.

Source code download:
https://grass.osgeo.org/grass64/source/
https://grass.osgeo.org/grass64/source/grass-6.4.5RC1.tar.gz

Binaries download:
https://grass.osgeo.org/download/software/#g64x

To get the GRASS GIS 6.4.5RC1 source code directly from SVN:
svn checkout https://svn.osgeo.org/grass/grass/tags/release_20150406_grass_6_4_5RC1/

Key improvements:
Key improvements of the GRASS GIS 6.4.5RC1 release include stability fixes (esp. vector library), some fixes for wxPython3 support, some module fixes, and more message translations.

See also our detailed announcement:
https://trac.osgeo.org/grass/wiki/Release/6.4.5RC1-News

First time users should explore the first steps tutorial after installation:
https://grasswiki.osgeo.org/wiki/Quick_wxGUI_tutorial

Release candidate management at
https://trac.osgeo.org/grass/wiki/Grass6Planning

Please join us in testing this release candidate for the final release.

Consider to donate pizza or beer for the next GRASS GIS Community Sprint (following the FOSS4G Europe 2015 in Como):
https://grass.osgeo.org/donations/

Thanks to all contributors!

About GRASS GIS

The Geographic Resources Analysis Support System (https://grass.osgeo.org), commonly referred to as GRASS GIS, is an Open Source Geographic Information System providing powerful raster, vector and geospatial processing capabilities in a single integrated software suite. GRASS GIS includes tools for spatial modeling, visualization of raster and vector data, management and analysis of geospatial data, and the processing of satellite and aerial imagery. It also provides the capability to produce sophisticated presentation graphics and hardcopy maps. GRASS GIS has been translated into about twenty languages and supports a huge array of data formats. It can be used either as a stand-alone application or as backend for other software packages such as QGIS and R geostatistics. It is distributed freely under the terms of the GNU General Public License (GPL). GRASS GIS is a founding member of the Open Source Geospatial Foundation (OSGeo).

The GRASS Development Team, April 2015

The GRASS GIS Development team has announced the release of the new major version GRASS GIS 7.0.0. This version provides many new functionalities including spatio-temporal database support, image segmentation, estimation of evapotranspiration and emissivity from satellite imagery, automatic line vertex densification during reprojection, more LIDAR support and a strongly improved graphical user interface experience. GRASS GIS 7.0.0 also offers significantly improved performance for many raster and vector modules: “Many processes that would take hours now take less than a minute, even on my small laptop!” explains Markus Neteler, the coordinator of the development team composed of academics and GIS professionals from around the world. The software is available for Linux, MS-Windows, Mac OSX and other operating systems.

Detailed announcement and software download:
https://grass.osgeo.org/news/42/15/GRASS-GIS-7-0-0/

About GRASS GIS
The Geographic Resources Analysis Support System https://grass.osgeo.org/, commonly referred to as GRASS GIS, is an open source Geographic Information System providing powerful raster, vector and geospatial processing capabilities in a single integrated software suite. GRASS GIS includes tools for spatial modeling, visualization of raster and vector data, management and analysis of geospatial data, and the processing of satellite and aerial imagery. It also provides the capability to produce sophisticated presentation graphics and hardcopy maps. GRASS GIS has been translated into about twenty languages and supports a huge array of data formats. It can be used either as a stand-alone application or as backend for other software packages such as QGIS and R geostatistics. It is distributed freely under the terms of the GNU General Public License (GPL). GRASS GIS is a founding member of the Open Source Geospatial Foundation (OSGeo).

GRASS GIS 7 just got better: When reprojecting vector data, now automated vertex densification is applied. This reduces the reprojection error for long lines (or polygon boundaries). The needed improvement has been kindly added in v.proj by Markus Metz.

Example

As an (extreme?) example, we generate a box in LatLong/WGS84 (EPSG: 4326) which is of 10 degree side length (see below for screenshot and at bottom for SHAPE file download of this “box” map):

[neteler@oboe ~]$ grass70 ~/grassdata/latlong/grass7/
# for the ease of generating the box, set computational region:
g.region n=60 s=40 w=0 e=30 res=10 -p
projection: 3 (Latitude-Longitude)
zone:       0
datum:      wgs84
ellipsoid:  wgs84
north:      60N
south:      40N
west:       0
east:       30E
nsres:      10
ewres:      10
rows:       2
cols:       3
cells:      6
# generate the box according to current computational region:
v.in.region box_latlong_10deg
exit

Next we start GRASS GIS in a metric projection, here the EU LAEA:

# EPSG 3035, metric EU LAEA:
grass70 ~/grassdata/europe_laea/user1/
GRASS 7.0.0svn (europe_laea): >

Now we first reproject the map the “traditional way” (no vertex densification as in most GIS, here enforced by smax=0):

v.proj box_latlong_10deg out=box_latlong_10deg_no_densification
location=latlong mapset=grass7 smax=0

Then we do a second reprojection with new automated vertex densification (here we use the default values for smax which is a 10km vertex distance in the reprojected map by default):

v.proj box_latlong_10deg out=box_latlong_10deg_yes_densification
location=latlong mapset=grass7

Eventually we can compare both reprojected maps:

g.region vect=box_latlong_10deg_no_densification

# compare:
d.mon wx0
d.vect box_latlong_10deg_no_densification color=red
d.vect box_latlong_10deg_yes_densification color=green fill_color=none
Comparison of the reprojection of a 10 degree large LatLong box to the metric EU LAEA (EPSG 3035): before in red and new in green. The grid is based on WGS84 at 5 degree spacing.

Comparison of the reprojection of a 10 degree large LatLong box to the metric EU LAEA (EPSG 3035): before in red and new in green. The grid is based on WGS84 at 5 degree spacing.

The result shows how nicely the projection is now performed in GRASS GIS 7: the red line output is old, the green color line is the new output (its area filling is not shown).

Consider to benchmark this with other GIS… the reprojected map should not become a simple trapezoid.

Sample dataset download

Download of box_latlong_10deg.shp for own tests (1kB).

The beautiful days in early November 2014 allowed to get some nice views of the Trentino (Northern Italy) – thanks to Landsat 8 and NASA’s open data policy:

Landsat 8: Northern Italy 1 Nov 2014
Landsat 8: Northern Italy 1 Nov 2014

Trento captured by Landsat8
Trento captured by Landsat8

Landsat 8: San Michele - 1 Nov 2014
Landsat 8: San Michele – 1 Nov 2014

The beauty of the landscape but also the human impact (landscape and condensation trails of airplanes) are clearly visible.

All data were processed in GRASS GIS 7 and pansharpened with i.fusion.hpf written by Nikos Alexandris.

Do you also sometimes get maps which contain zero (0) rather than NULL (no data) in some parts of the map? This can be easily solved with “floodfilling”, even in a GIS.

My original map looks like this (here, Trentino elevation model):

The light blue parts should be no data (NULL) rather than zero (0)...

Now what? In a paint software we would simply use bucket fill but what about GIS data? Well, we can do something similar using “clumping”. It requires a bit of computational time but works perfectly, even for large DEMs, e.g., all Italy at 20m resolution. Using the open source software GRASS GIS 7, we can compute all “clumps” (that are many for a floating point DEM!):

# first we set the computational region to the raster map:
g.region rast=pat_DTM_2008_derived_2m -p
r.clump pat_DTM_2008_derived_2m out=pat_DTM_2008_derived_2m_clump

The resulting clump map produced by r.clump is nicely colorized:

Clumped map derived from DEM (generated with r.clump)

As we can see, the area of interest (province) is now surrounded by three clumps. With a simple map algebra statement (r.mapcalc or GUI calculator) we can create a MASK by assigning these outer boundary clumps to NULL and the other “good” clumps to 1:

r.mapcalc "no_data_mask = if(pat_DTM_2008_derived_2m_clump == 264485050 || \
  pat_DTM_2008_derived_2m_clump == 197926480 || \
  pat_DTM_2008_derived_2m_clump == 3, null(), 1)"

This mask map looks like this:

Mask map from all clumps except for the large outer clumps

We now activate this MASK and generate a copy of the original map into a new map name by using map algebra again (this just keeps the data matched by the MASK). Eventually we remove the MASK and verify the result:

# apply the mask
r.mask no_data_mask
# generate a copy of the DEM, filter on the fly
r.mapcalc "pat_DTM_2008_derived_2m_fixed = pat_DTM_2008_derived_2m"
# assign a nice color table
r.colors pat_DTM_2008_derived_2m_fixed color=srtmplus
# remove the MASK
r.mask -r

And the final DEM is now properly cleaned up in terms of NULL values (no data):

DEM cleaned up for no data

Enjoy.

brainscan1Last year (2013) I “enjoyed” a brain CT scan in order to identify a post-surgery issue – luckily nothing found. Being in Italy, like all patients I received a CD-ROM with the scan data on it: so, something to play with! In this article I’ll show how to easily turn 2D scan data into a volumetric (voxel) visualization.

The CT scan data come in a DICOM format which ImageMagick is able to read and convert. Knowing that, we furthermore need the open source software packages GRASS GIS 7 and Paraview to get the job done.

First of all, we create a new XY (unprojected) GRASS location to import the data into:

# create a new, empty location (or use the Location wizard):
grass70 -c ~/grassdata/brain_ct

We now start GRASS GIS 7 with that location. After mounting the CD-ROM we navigate into the image directory therein. The directory name depends on the type of CT scanner which was used in the hospital. The file name suffix may be .IMA.

Now we count the number of images, convert and import them into GRASS GIS:

# list and count
LIST=`ls -1 *.IMA`
MAX=`echo $LIST | wc -w`

# import into XY location:
curr=1
for i in $LIST ; do

# pretty print the numbers to 000X for easier looping:
curr=`echo $curr | awk ‘{printf “%04d\n”, $1}’`
convert “$i” brain.$curr.png
r.in.gdal in=brain.$curr.png out=brain.$curr
r.null brain.$curr setnull=0
rm -f brain.$curr.png
curr=`expr $curr + 1`

done

At this point all CT slices are imported in an ordered way. For extra fun, we can animate the 2D slices in g.gui.animation:

Animation of brain scan slices
(click to enlarge)

# enter in one line:
g.gui.animation rast=`g.mlist -e rast separator=comma pattern=”brain*”`

The tool allows to export as animated GIF or AVI:

Animation of brain scan slices (click to enlarge)

Now it is time to generate a volume:

# first count number of available layers
g.mlist rast pat=”brain*” | wc -l

# now set 3D region to number of available layers (as number of depths)
g.region rast=brain.0003 b=1 t=$MAX -p3

At this point the computational region is properly defined to our 3D raster space. Time to convert the 2D slices into voxels by stacking them on top of each other:

# convert 2D slices to 3D slices:
r.to.rast3 `g.mlist rast pat=”brain*” sep=,` out=brain_vol

We can now look at the volume with GRASS GIS’ wxNVIZ or preferably the extremely powerful Paraview. The latter requires an export of the volume to VTK format:

# fetch some environment variables
eval `g.gisenv -s`
# export GRASS voxels to VTK 3D as 3D points, with scaled z values:
SCALE=2
g.message “Exporting to VTK format, scale factor: $SCALE”
r3.out.vtk brain_vol dp=2 elevscale=$SCALE \
output=${PREFIX}_${MAPSET}_brain_vol_scaled${SCALE}.vtk -p

Eventually we can open this new VTK file in Paraview for visual exploration:

# show as volume
# In Paraview: Properties: Apply; Display Repres: volume; etc.
paraview –data=brain_s1_vol_scaled2.vtk

markus_brain_ct_scan3 markus_brain_ct_scan4 markus_brain_ct_scan2

 

 

 

 

 

 

 

 

 

 

 

 

Fairly easy!

BTW: I have a scan of my non-smoker lungs as well :-)

Drowning in too many maps? Have some fun exploring fascinating geometries of changing landscapes in Space Time Cube and creating 2D and 3D animations from time series of geospatial data. Learn about the new capabilities for spatio-temporal data handling in GRASS GIS 7 (https://grass.osgeo.org/grass7/) and explore various techniques for dynamic visualizations.

First, we will introduce you to GRASS GIS 7, including its spatio-temporal capabilities and you will learn how to manage and analyze geospatial data time series. Then, we will explore new tools for visualization of spatio-temporal data. You will create both 2D and 3D dynamic visualizations directly in GRASS GIS 7. Additionally, we will explain the Space Time Cube concept using various applications based on raster and vector data time series. You will learn to manage and visualize data in space time cubes (voxel models). No prior knowledge of GRASS GIS is necessary, we will cover the basics needed for the workshop. All relevant material including an overview of the tools and hands-on practical instructions along with the sample data sets will be available on-line. And, by the way, GRASS GIS is a free and open source geographic information system (GIS) used for geospatial data management, analysis, modeling, image processing, and visualization which runs on Linux, MS Windows, Mac OS X and other systems.

Presenters: Vaclav Petras, Anna Petrasova, Helena Mitasova, Markus Neteler

When: FOSS4G 2014, Sept 8th-13th 2014, Portland, OR, USA

Register at: https://2014.foss4g.org/schedule/workshops/#wshop-526

banner_lidar

[toc]

LiDAR point cloud data are commonly delivered in the ASPRS LAS format. The format is supported by libLAS, a BSD-licensed C++ library for reading/writing these data. GRASS GIS 7 supports the LAS format directly when built against libLAS (as the case for most binary packages being available for download).

In this exercise we will import a sample LAS data set covering a tiny area close to Raleigh, NC (USA), belonging to the North Carolina sample data set. Sample LAS data download: https://grass.osgeo.org/sampledata/north_carolina/  (25MB).

For a full exercise, we will, however, assume that no GRASS GIS location is ready so far (so: newbies are welcome!) and create a new one initially.

1. Having the LAS file: now what?

In the first place, check the metadata of the LAS file using the lasinfo command (comes with libLAS; here only parts of the output shown):

lasinfo points.las
---------------------------------------------------------
  Header Summary
---------------------------------------------------------
  Version:  1.2
  Source ID:  0
  Reserved:  0
  Project ID/GUID:  '00000000-0000-0000-0000-000000000000'
  System ID:  'libLAS'
  Generating Software:  'libLAS 1.2'
[...]
  Spatial Reference:
None
[...]
---------------------------------------------------------
  Point Inspection Summary
---------------------------------------------------------
  Header Point Count: 1287775
  Actual Point Count: 1287775

  Minimum and Maximum Attributes (min,max)
---------------------------------------------------------
  Min X, Y, Z:   6066629.86, 2190053.45, -3.60
  Max X, Y, Z:   6070237.92, 2193507.74, 906.00
  Bounding Box:  6066629.86, 2190053.45, 6070237.92, 2193507.74
  Time:  0.000000, 0.000000
  Return Number:  1, 3
  Return Count:  1, 7
  Flightline Edge:  0, 0
  Intensity:  0, 256
  Scan Direction Flag:  0, 0
  Scan Angle Rank:  0, 0
  Classification:  2, 7
  Point Source Id:  0, 0
  User Data:  0, 0
  Minimum Color (RGB):  0 0 0 
  Maximum Color (RGB):  0 0 0 

  Number of Points by Return
---------------------------------------------------------
  (1) 1225886  (2) 61430  (3) 459

  Number of Returns by Pulse
---------------------------------------------------------
  (1) 30877  (2) 153  (5) 1225886  (6) 30706  (7) 153

  Point Classifications
---------------------------------------------------------
  647337 Ground (2) 
  639673 Low Vegetation (3) 
  740 Building (6) 
  25 Low Point (noise) (7) 
  -------------------------------------------------------
  0 withheld
  0 keypoint
  0 synthetic
  -------------------------------------------------------

We see: no spatial reference system indicated!
Luckily we know from here that the projection is  NAD83(HARN) / North Carolina, LCC 2SP metric, EPSG code 3358). Furthermore we see:

  • Number of Points by Return: 3 (i.e., first, mid, last)
  • Point Classifications: the points are already classified as “Ground” (class 2),  “Low Vegetation” (3),  “Building” (6), and Low Point (noise) (class 7). Something to play with later.

Time to create a GRASS GIS location and import the LAS file.

2. Creating a GRASS GIS location for the LAS file

Since we know the EPSG code of the projection, that’s an easy task. Please note that GRASS GIS can generate locations directly from SHAPE files (with .prj file), GeoTIFF and more.

We fire up GRASS GIS 7 and open the Location Wizard:

grass7_startup

In the Location Wizard, we first define a name for the new location:

grass7_loc_wizard1

We select the “EPSG code” method for creating a new location:

grass7_loc_wizard2

You can search for “North Carolina” and select the EPSG code 3358 from the list:

grass7_loc_wizard3

Next summary should show up as follows (be sure to have the metric projection shown!):

grass7_loc_wizard4

With “Finish” you reach this notification (indeed, nothing to change! It is already fine):

grass7_loc_wizard5

Since we want to import the LAS file, no need to manually define any region extent here – just say “No”:

grass7_loc_wizard6

While we could import the data also into the PERMANENT mapset, we prefer to create an own mapset “lasdata” for our LAS data (once you reach hundreds of maps to manage, you will be happy about the concept of mapsets):

grass7_loc_wizard7

Voilà, we get back to the initial startup screen and can now start our GRASS GIS session with our “nc_nad83_lcc” location and “lasdata” mapset within the location: “Start GRASS”!

grass7_loc_wizard8

3. Import of the LAS file

When creating a new location from a GeoTIFF or SHAPE file (or other GDAL supported format), then the data set is imported right away. This is not the case for LAS files, also due to the fact that we can directly apply binning statistics during import of the LAS file (e.g. percentiles, min or max) and create a raster surface from the points right away rather than importing them as vector points.

3. a) Creating a raster surface from LAS during import

The LAS import into a raster surface is available through r.in.lidar:

grass7_las_import1

First the LAS file needs to be selected and an output file name specified (in this example, we want to extract the 95th percentile as binning method, hence a reasonable map name):

grass7_las_import2

In the “Statistic” tab, we select the “percentile” method along with 95 as value:

grass7_las_import3

In the “Optional” tab we activate to extend the computational region from the LAS file and, since the spatial reference system metadata are lacking from the LAS file, also “override dataset projection” to use that predefined in the location. Finally, we define 5m as desired raster resolution for the resulting raster map:

grass7_las_import4

Upon conpletion of the import/binning, the new raster elevation map is shown after zooming to the map (r.in.lidar -e … restores upon completion the previous region settings, hence we may have to zoom):

grass7_las_import5

Now we can start to analyze or visually explore the imported LAS file.

4. Visual LiDAR data exploration

Using the wxNVIZ 3D viewer, we can easily fly through the new DEM. Switch in the Map Display to “3D view” (1). Note that the default rendering is initially done at low visual resolution for speed reasons. You can switch to “Rotate mode” as well to easily navigate with the mouse. In the “Data” tab (2) you can increase the visual resolution (3) to obtain a crisp view:

grass7_las_viz1

 

Now all kinds of analysis steps may follow.

Outlook

For true LiDAR processing as points, see the following GRASS GIS 7 modules: v.in.lidar (for point import), v.lidar.correction, v.lidar.edgedetection, v.lidar.growing, v.outlier, and v.surf.bspline.