28

Sep

2016

The Amazing Race!

Author: Zachary Norman

It wasn't so long ago that the IDL-Python bridge was introduced to IDL 8.5. It was with this new version, that I got my first experience programming with Python and testing the IDL-Python bridge. Through the past year it has been exciting to see the new changes and improvements that have become a part of the bridge.



Some of these new features include:

-Better error catching with the IDL-Python bridge

-Enhanced Jupyter notebook that allows for the development of full IDL programs and Python code in the same environment

-Improved support for variables passing back and forth


With all the time I have spent with Python, I have always wondered what some of the advantages are between Python and IDL. One thing that I have commonly heard several engineers say was that IDL was much faster than Python. For this blog, I decided to put that to the test and see how Python and IDL really compared to one another.



Before talking about the test, I do just want to explain things a bit about how it was set up and some potential caveats about the processing times that will be shown. With the tests I created, I did my best to choose tests that were comparable between IDL and Python. Since I'm no expert at Python, there very well may have been other methods that could be faster than what I will show. Most of the pieces I included in the test are things I found easily by doing a web search - meaning that most of the approaches I used were the most common programming methods that people are likely using. This shows how much faster IDL might be than a general program than something that someone might write in Python.





The test:

Here is what was actually tested between IDL and Python with an array of [10000,10000] or 10000*10000 elements

-Array creation time
-Type conversion times
-Index array creation times (i.e. [0,1,2,3,4...,n-1])

-Incrementing array values of all elements by 1

-Complex math expression with array (exact equation: sqrt(sin(arr*arr)))

-Single threaded for IDL and multithreaded

-Array element access times (i.e. setting y = arr[i])

-Simple image processing filter times (filters: sobel, roberts, prewitt)




The results:

Average array creation time (seconds):

    Python : 0.213000 +/- 0.00953933

    IDL    : 0.0936666 +/- 0.0155028

Total time (seconds):

    Python : 0.639000

    IDL    : 0.281000

Python/IDL time ratio: 2.27402

Average array data type conversion time (seconds):

     Python : 0.171333 +/- 0.0155028

    IDL    : 0.0730000 +/- 0.00866031

Total time (seconds):

     Python : 0.514000

     IDL    : 0.219000

Python/IDL time ratio: 2.34703





Average index array creation time (seconds):

     Python : 0.229000 +/- 0.00866031

     IDL    : 0.124667 +/- 0.0160104

Total time (seconds):

    Python : 0.687000

    IDL    : 0.374000

Python/IDL time ratio: 1.83690





Average increasing array value time (seconds):

     Python : 0.0933333 +/- 0.000577446

     IDL    : 0.0313334 +/- 0.000577377

Total time (seconds):

    Python : 0.280000

    IDL    : 0.0940001

Python/IDL time ratio: 2.97872





Average complex math statements (1 thread) time (seconds):

    Python : 6.36967 +/- 0.0645319

    IDL    : 8.34667 +/- 0.0155028

Total time (seconds):

    Python : 19.1090

    IDL    : 25.0400

Python/IDL time ratio: 0.763139





Average complex math statements (8 thread) time (seconds):

    Python : 6.34400 +/- 0.0321871

    IDL    : 1.93933 +/- 0.00923762

Total time (seconds):

    Python : 19.0320

    IDL    : 5.81800

Python/IDL time ratio: 3.27123





Average loop through array element time (seconds):

    Python : 11.5290 +/- NaN

    IDL    : 3.29100 +/- NaN

Total time (seconds):

    Python : 11.5290

    IDL    : 3.29100

Python/IDL time ratio: 3.50319





Average image processing routines time (seconds):

    Python : 15.3660 +/- 0.0829635

    IDL    : 1.39900 +/- 0.0238955

Total time (seconds):

    Python : 46.0980

    IDL    : 4.19700

Python/IDL time ratio: 10.9836







Conclusion:



In short, IDL significantly outperformed Python will all the speed tests apart from the complex math statement. However, IDL has access to built in multithreading for large arrays and, with multithreading enabled, IDL outperforms Python significantly when using all available cores.



Below is the IDL code used to compare the processing speed of IDL and Python. To use it you will need a few Python modules which can be found at the beginning of the procedure "python_test". 

IDL Code:
pro printTimes, pythontimes, idltimes, TIMED = timed
  compile_opt idl2
  
  ;check if we have a name for the process to print with the mean
  if ~keyword_set(timed) then begin
    add = ' time (seconds):'
  endif else begin
    add = timed + ' time (seconds):'
  endelse
  
  print, 'Average ' + add 
  print, '  Python : ' + strtrim(pythontimes.mean(),2) + ' +/- ' + strtrim(stddev(pythontimes),2)
  print, '  IDL    : ' + strtrim(idltimes.mean(),2) + ' +/- ' + strtrim(stddev(idltimes),2)
  print, 'Total time (seconds):'
  print, '  Python : ' + strtrim(total(pythontimes),2)
  print, '  IDL    : ' + strtrim(total(idltimes),2)
  print, 'Python/IDL time ratio: ' + strtrim(total(pythontimes)/total(idltimes),2)
  
end

pro python_test
  compile_opt idl2
  
  
  ;initialize Python
  >>> 'import numpy as np'
  >>> 'import time'
  >>> 'import math'
  >>> 'from idlpy import *'
  >>> 'from skimage.filters import roberts, sobel, prewitt'
  
  ;number of times we want to run each test
  nloops = 3
  
  ;array dimensions for iterator tests
  dims = 10000
  ;>>>'dims = ' + strtrim(dims,2)
  python.dims = dims
  
  ;array dimensions for filter tests
  dims_filter = 10000
  python.dims_filter = dims_filter
  
  ;initialize arrays to hol dpython times
  pythontimes = fltarr(nloops)
  idltimes = fltarr(nloops)
  
  ;test array creation in Python and IDL
  for i=0,nloops-1 do begin
    tStart = systime(/seconds)
    >>> 'arr = np.ones((dims,dims))'
    tEndPython = systime(/seconds)
    arr = fltarr(dims,dims)
    tEndIDL = systime(/seconds)
    
    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes, idltimes, TIMED = 'array creation'

  ;check type conversion times
  for i=0,nloops-1 do begin
    tStart = systime(/seconds)
    >>> 'arr2 = arr.astype(int)'
    tEndPython = systime(/seconds)
    arr2 = fix(arr)
    tEndIDL = systime(/seconds)
    
    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes, idltimes, TIMED = 'array data type conversion'

  ;check index array creation times
  for i=0,nloops-1 do begin
    tStart = systime(/seconds)
    >>> 'arr = np.arange(0, dims*dims, dtype=long)'
    tEndPython = systime(/seconds)
    arr = lindgen(dims*dims)
    tEndIDL = systime(/seconds)

    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes, idltimes, TIMED = 'index array creation'

  ;check adding values to array
  for i=0,nloops-1 do begin
    tStart = systime(/seconds)
    >>> 'arr += 1'
    tEndPython = systime(/seconds)
    arr += 1
    tEndIDL = systime(/seconds)

    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes, idltimes, TIMED = 'increasing array value'

  ;check complex math expressions with a single thread
  pref_set, 'IDL_CPU_TPOOL_NTHREADS', 1, /commit
  for i=0,nloops-1 do begin
    tStart = systime(/seconds)
    >>> 'y = np.sin(arr*arr)**.5'
    tEndPython = systime(/seconds)
    y = sqrt(sin(arr*arr))
    tEndIDL = systime(/seconds)

    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes, idltimes, TIMED = 'complex math statements (' + strtrim(!CPU.TPOOL_NTHREADS,2) + ' thread)'
  pref_set, 'IDL_CPU_TPOOL_NTHREADS', /default, /commit
  
  ;check complex math expressions with all threads
  for i=0,nloops-1 do begin
    tStart = systime(/seconds)
    >>> 'y = np.sin(arr*arr)**.5'
    tEndPython = systime(/seconds)
    y = sqrt(sin(arr*arr))
    tEndIDL = systime(/seconds)

    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes, idltimes, TIMED = 'complex math statements (' + strtrim(!CPU.TPOOL_NTHREADS,2) + ' thread)'

  ;check array element access times
  nhere = 1
  for i=0,nhere-1 do begin
    tStart = systime(/seconds)
    >>>'for x in np.nditer(arr):\n    y=x'
    tEndPython = systime(/seconds)
    foreach x, arr do y = x
    tEndIDL = systime(/seconds)

    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes[0:nhere-1], idltimes[0:nhere-1], TIMED = 'loop through array element'

  ;check times for image processing
  im_dat = lonarr(dims_filter, dims_filter)
  ;set a square in the middle of the image to 1
  im_dat[.4*dims_filter:.6*dims_filter,.4*dims_filter:.6*dims_filter] = 1
  ;send the array to Python as well
  python.im_dat = im_dat
  
  for i=0,nloops-1 do begin
    tStart = systime(/seconds)
    >>> 'edge_sobel = sobel(im_dat)'
    >>> 'edge_roberts = roberts(im_dat)'
    >>> 'edge_prewitt = prewitt(im_dat)'
    tEndPython = systime(/seconds)
    edge_sobel = sobel(im_dat)
    edge_roberts = roberts(im_dat)
    edge_prewitt = prewitt(im_dat)
    tEndIDL = systime(/seconds)

    ;save the times
    pythontimes[i] = tEndPython - tStart
    idltimes[i] = tEndIDL - tEndPython
  endfor
  print
  printTimes, pythontimes, idltimes, TIMED = 'image processing routines'

  stop
end

Comments (0) Number of views (1856) Article rating: 5.0

Categories: IDL Blog | IDL Data Point

Tags:

23

Sep

2016

Coming Soon in ENVI 5.4: ENVITask Style Sheets

Author: Brian Griglak

Since the early days of ENVITasks in ENVI 5.2, there has been what we call the Dynamic UI.  Some of the toolbox items that wrap ENVITasks use it, but you can use it at any time to set some or all of the parameters of your task using the ENVIUI::SelectTaskParameters() method.  When you call this method, it will inspect the task to get the NAME, TYPE, and VALUE properties of each parameter, and construct a dialog with widgets for each.  If the VALUE is !NULL, then the widget will be empty, but if it has a valid value the widget will display that.  We have a set of heuristics that try to select the best widget to use for each TYPE, but it may not choose the widget you desire for a particular parameter.

To address this, we are introducing task style sheets in the upcoming release of ENVI 5.4.  The idea is borrowed from HTML and CSS, which separate content from presentation.  In the ENVITask context, the ENVITask template defines the content, and the style sheet defines the presentation logic.  These style sheets are nowhere as sophisticated as CSS can be, starting with the lack of cascading, but they don’t need to be for our purposes.  The style sheet is a JSON file or in-memory Hash that allows you to hide a couple of the widgets at the bottom of the dialog, and specify which widget to use for each parameter.  You don’t need to specify every parameter, those that are not in the style sheet will fall back to the heuristics.

Let’s make this concrete using the example from the help page:

  e = ENVI()
  Task = ENVITask('NNDiffusePanSharpening')
  ret = e.UI.SelectTaskParameters(Task)

The resulting pop-up dialog looks like this:

The “Pixel Size Ratio” parameter has to be an integer, so we used the edit widget with the increment/decrement buttons next to it.  But maybe you want to use a wheel instead; this is where a style sheet comes into play.

I can change the “Pixel Size Ratio” parameter to use a wheel widget, like  this:

  style = Hash("schema", "envitask_3.0", $
               "parameters", List(Hash("name", "pixel_size_ratio", $
                                       "type", "IDLWheel_UI")))
  ret = e.UI.SelectTaskParameters(Task, STYLE_SHEET=style)

One key point is that the “parameters” value must be a List of Hashes, even for a single parameter.  The equivalent JSON file is this:

{
    "schema": "envitask_3.0",
    "parameters": [
        {
            "name": "pixel_size_ratio",
            "type": "IDLWheel_UI"
        }
    ]
}

The dialog now looks like this:

All of the possible values for “type” are well documented, under the heading “User Interface Elements”.  Some of the widgets allow you to specify extra properties to customize the behavior.  The IDLWheel_UI class, for example, supports specification of MIN and MAX properties.  If the parameter has these properties set, they will be used, but as you can see in the case of NNDifusePanSharpening it doesn’t.  But we can specify them in the style sheet if we want to add them.  We do this by putting all these properties into a Hash under the key “keywords”:

  style = Hash("schema", "envitask_3.0", $
               "parameters", List(Hash("name", "pixel_size_ratio", $
                                       "type", "IDLWheel_UI", $
                                       "keywords", Hash("min", 1, $
                                                        "max", 9))))
  ret = e.UI.SelectTaskParameters(Task, STYLE_SHEET=style)

The equivalent JSON file is this:

{
    "schema": "envitask_3.0",
    "parameters": [
        {
            "name": "pixel_size_ratio",
            "type": "IDLWheel_UI",
            "keywords": {
                "max": 9,
                "min": 1
            }
        }
    ]
}

The dialog now looks like this:

In all the examples, I have manually specified the STYLE_SHEET keyword in the call to SelectTaskParameters().  If you want the style sheet to be automatically used any time you call this method with your task, then create the JSON file named <task name>.style and put it in the same folder as the .task file.

Comments (0) Number of views (1285) Article rating: 5.0

Categories: IDL Blog | IDL Data Point

Tags:

20

Sep

2016

Using Attributes to Improve Image Classification Accuracy

Author: Jason Wolfe

In this article I will show an example of how you can improve the accuracy of a supervised classification by considering different attributes.

Supervised classification typically involves a multi-band image and ground-truth data for training. The classification is thus restricted to spectral information only; for example, radiance or reflectance values in various wavelengths. However, adding different types of data to the source image gives a classifier more information to consider when assigning pixels to classes.

Attributes are unique characteristics that can help distinguish between different objects in an image. Examples include elevation, texture, and saturation. In ENVI you can create a georeferenced layer stack  image where each band is a separate attribute. Then you can classify the attribute image using training data.

Creating an attribute image is an example of data fusion in remote sensing. This is the process of combining data from multiple sources to produce a dataset that contains more detailed information than each of the individual sources.

Think about what attributes will best distinguish between the different classes. Spectral indices are easy to create and can help identify different features. Vegetation indices are a good addition to an attribute image if different vegetation types will be classified.

Here are some other attributes to experiment with:

  • Texture measures (occurrence and co-occurrence)
  • Synthetic aperture rader (SAR) data
  • Color transforms (hue, saturation, intensity, lightness)

Including elevation data with your source imagery helps to identify features with noticeable height such as buildings and trees. If point clouds are available for the study area, you can use ENVI LiDAR to create digital surface model (DSM) and digital elevation model (DEM) images. Use Band Math to subtract the DEM from the DSM to create a height image, where pixels represent absolute heights in meters. In the following example, the brightest pixels represent trees. You can also see building outlines.

Example height image at 1-meter spatial resolution

Locating coincident point cloud data for your study area (preferably around the same date) can be a challenge. Luckily I found the perfect set of coincident data for a classification experiment.

Example

The National Ecological Observatory Network (NEON) provides field measurements and airborne remote sensing datasets for terrestrial and aquatic study sites throughout the world. Their Airborne Observation Platform (AOP) carries an imaging spectrometer with 428 narrow spectral bands extending from 380 to 2510 nanometers with a spectral sampling of 5 nanometers. Also onboard is a full-waveform LiDAR sensor and a high-resolution red/green/blue (RGB) camera. For a given observation site, the hyperspectral data, LiDAR data, and high-resolution orthophotography are available for the same sampling period.

I acquired a sample of NEON data near Grand Junction, Colorado from July 2013. My goal was to create an urban land-use classification map using an attribute image and training data. For experimentation and to reduce processing time, I only extracted the RGB and near-infrared bands from the hyperspectral dataset and created a multispectral image with these bands. I used ENVI LiDAR to extract a DEM and DSM from the point clouds. Then I created a height image whose pixels lined up exactly with the multispectral image.

I created an Enhanced Vegetation Index (EVI) image using the Spectral Indices tool. Finally, I combined the RGB/NIR images, the relative height image, and the EVI image into a single attribute image using the Layer Stacking tool.

Next, I used the Region of Interest (ROI) tool to collect samples of pixels from the multispectral image that I knew represented five dominant classes in the scene: Water, Asphalt, Concrete, Grass, and Trees. I used a NEON orthophoto to help verify the different land-cover types.

I ran seven different supervised classifiers with the multispectral image, then again with the attribute image. Here are some examples of Maximum Likelihood classifier results:

Maximum Likelihood classification result from the NEON multispectral image

Notice how the classifier assigned some Building pixels to Asphalt throughout the image. 

Here is an improved result using the attribute image. Those pixels are correctly classified as Building now.

Maximum Likelihood classification result using the attribute image

A confusion matrix and accuracy metrics can help verify the accuracy of the classification.

Confusion matrix calculated from the attribute image classification

The following table shows how the Overall Accuracy value is higher with the attribute image when using different supervised classifiers:

  Multispectral image  Attribute Image 
 Mahalanobis Distance  72.5  83.8
 Minimum Distance  57.69  95.22
 Maximum Likelihood  91.3  98.91
 Parallelepiped  57.18  95.8
 Spectral Angle Mapper  54.97  61.79
 Spectral Information Divergence  62.1  66.93
 Support Vector Machine  85.03  99.16

The accuracy of a supervised classification depends on the quality of your training data as well as a good selection of attributes. In some cases, too many attributes added to a multispectral image can make the classification result worse, so you should experiment with what works best for your study area. Also, some classifiers work better than others when considering different spatial and spectral attributes. Finally, you may need to normalize the different data layers in the attribute image if their pixel values have a wide variation.

New tutorials will be available in ENVI 5.4 for creating attribute images and for defining training data for supervised classification.

Resources

National Ecological Observatory Network. 2016. Data accessed on 28 July 2016. Available on-line at http://www.neonscience.org from National Ecological Observatory Network, Boulder, CO, USA.

Comments (0) Number of views (2288) Article rating: 5.0

Categories: ENVI Blog | Imagery Speaks

Tags:

13

Sep

2016

Managing FMV with Jagwire and Passive Analytics

Author: Albert Spratley

The rapid growth of unmanned aerial vehicles (UAVs) and payloads has resulted in an ever growing deluge of data that has to be archived, sorted, analyzed, and distributed to consumers across the defense, agriculture, and utility markets. In many cases, especially in the case of full motion video (FMV), a single flight can result in several hours of data that has to be viewed and analyzed. Often only a small fraction of that data is useful for analysis purposes. For larger UAV fleets, with multiple, simultaneous missions, substantial resources are required to perform the analysis. The resources required to analyze these data products increases cost proportionally to the amount of data collected.

For systems that have adopted the use of properly formatted metadata, we can attempt to filter this glut of data by analyzing patterns and attempting to infer some operator intent based on domain knowledge. For example, identifying temporal “pauses” for the sensor center field of view may indicate an area or point of interest for further analysis. Circular patterns in the sensor center field of view could indicate the inspection of a building, object, or structure of significance. Smooth “pans” during the video or “sweeping” motions across the ground can infer a collection aimed at covering an area on the ground.

Harris Geospatial Solutions has designed and prototyped algorithms capable of identifying these useful segments of video by analyzing the metadata embedded within the video stream. These “passive analytics” run in real time, during the UAV flight, and identify sub-sections of video that are far more likely to be useful in a more detailed analysis. By dynamically detecting, and setting aside these sub-clips of video, the burden of first-phase analysis can be greatly reduced, allowing the user to focus their analytical and dissemination resources on meeting the challenges of their market space rather than wading through a sea of irrelevant data.

Comments (0) Number of views (1857) Article rating: No rating

Categories: ENVI Blog | Imagery Speaks

Tags:

1

Sep

2016

App Update Check - Coming soon in IDL 8.6 and ENVI 5.4

Author: Jim Pendleton

First, a dose of nostalgia.

In the early days of commercial IDL, new drops of the application occurred when new features were added and as bugs were fixed for important customers. It was not unusual for there to be a couple different releases per month, a release cycle concept not dissimilar to that of most of today's smart phone apps.

The new releases were signaled through an update to a README-type of text file that was available on a node on a computer network called DECNET, shared by many universities and national labs in the days before the internet as we know it today.

At the organization where I worked, I had a nightly batch job that would alert me to the availability of a new version of IDL. It was a highlight of my week when new features appeared in IDL. Yes, I had that sort of a life back then. Frequent updates, from a customer point of view, kept the buzz alive with the products.  

Here we are about 30 years later and this functionality to check for updates is now built into the upcoming releases of both IDL 8.6 and ENVI 5.4. By default, IDL will check once a week automatically for updates when a new session is started.

Additionally, with the click of a single "Help/Check for Updates..." menu button in the IDL Workbench or from the ENVI UI, you can ensure that you, too, are enjoying the freshest and best we have to offer.

If your version is current, you will receive a dialog assuring you all is well.

To check that your ENVI version is current, find the menu item in the main application user interface.

You will need an internet connection to the commercial world for the status check to complete successfully.

Extra "80s-throwback" credit if you use a packet sniffer to find the syntax of the request and write a script that will automatically launch a WGET query each time you start a new IDL session.

Does this new feature signal a return to more frequent application updates? Have we gone back to the future? Only time will tell.

Comments (0) Number of views (2460) Article rating: No rating

Categories: IDL Blog | IDL Data Point

Tags:

12345678910Last

MONTHLY ARCHIVE

«March 2017»
SunMonTueWedThuFriSat
2627281234
567891011
12131415161718
19202122232425
2627282930311
2345678

MOST POPULAR POSTS

AUTHORS

Authors

Authors

Authors

Authors

Authors

Authors

Authors

Authors

Authors

Authors

GUEST AUTHORS

Authors

Authors

Authors

Authors

Authors

Authors

Authors



© 2017 Exelis Visual Information Solutions, Inc., a subsidiary of Harris Corporation