For a good portion of my life I have been a space nerd and I have always enjoyed looking at images of other planets in our solar system, nebulae, and other galaxies. Recently there was a new addition to the group of photographed celestial objects: Pluto. On July 14th, 2015 the New Horizons probe flew within 12,600 km of the surface of Pluto and took some awesome images. With these new images, I thought it would be a neat idea to try and mosaic a few of the images together using ENVI to see if it could be done.
In order to get some images that would be suitable for mosaicking, I needed to make sure that they were taken in quick succession so that they would have a good amount of overlap with each other. This allows ENVI to be able to register the images to one another without too much trouble. To find these images I got data from New Horizon's instrument LORRI at the following link:
It was really important when selecting images to mosaic that the images needed to have similar acquisition times. After looking through the available images from LORRI, I found 13 images from the acquisition time from 11:36:00 UTC to 11:36:36 UTC. It is great that the images of Pluto are available, but unfortunately there was no metadata associated with the data to accurately geolocate the images. To work around this missing metadata, I decided that I would just create a dummy spatial reference centered around a reference image. For my base image, I chose to use the picture taken at 11:36:00 UTC and gave this image a basic WGS-84 coordinate system with fake pixel sizes of 1 meter by 1 meter. At this point it is important to note that my fake coordinate system and pixel sizes likely do not correspond at all to the actual coordinate system or pixel sizes, I just needed to create a spatial reference that I could use to register and mosaic the images using ENVI.
With my reference image and fake spatial reference, the next step was to figure out how to get the images registered properly to one another. After doing some investigation with the Image Registration Workflow, I realized that in order to get the best results I needed to have the images placed as close to their actual relative positions as possible. To do this, I looked at my reference image and first image to be registered in the series (reference image had the acquisition time of 2015-07-14 11:36:00 UTC and the first image to be registered was acquired at 2015-07-14 11:36:03 UTC).
After looking closely at the two images and tweaking things by hand a bit, I found that the approximate pixel offset of the second image with respect to the first was about [-580, -18] pixels in the x and y directions. This mean that the first image in series after the base image had an offset of [-580, -18] pixels and the second image had an offset of [-1160, -36] pixels. At this point I decided to use the ENVI API to programmatically apply this offset to each image in the series. A big reason for this choice to use the ENVI API is that I wouldn't have to step through the Image Registration Workflow for every image in my series and, if I changed the offset at all, it would be quick to apply to a different image.
After I re-projected each image based on my pixel offset the images were relatively close at this point. All that was left was to apply color balancing (through a custom histogram matching function), register the images to one another, and mosaic them together. Here is a quick summary of the steps that I performed using the ENVI API:
Apply color matching to reference image -> Use task GenerateTiePointsByCrossCorrelation -> Use task FilterTiePointsByGlobalTransform -> Use task ImageToImageRegistration -> Use task BuildMosaicRaster.
Below is the output mosaic and the IDL code used to create the image.
;+
; :Author: Zach Norman
;
; :Description:
; Small example program to register and mosaic some Pluto images from the New Horizon's
; LORRI instrument. Save this code as pluto_mosaicking.pro and to run you can press 'Run'
; in the IDL Workbench. See the output section below for where to find the results of the
; program and details on the results.
;
; The code is written to automatically download the images of Pluto to a directory called images
; that will be created in the same location that this .pro file will reside. If there are any
; issues downloading the files then you may have to download them from the web by hand. This
; could be necessary if the URLs change for the images from NASA. You may need to
; make some additional changes as well if you download the images by hand.
;
; Specifically, you may have to change some of the parameters in the procedure pluto_mosaicking
; so that they are current for your machine. If you download the Pluto images from
; http://pluto.jhuapl.edu/soc/Pluto-Encounter/index.php then you will need to specify
; the value of "imagedir" in the pluto_mosaicking procedure.
;
; If you want to or have to download the images by hand, there are 13 images which have
; acquisition times ranging from 2015-07-14 11:36:00 UTC to 2015-07-14 11:36:36 UTC. As
; is, the code is written so that you place this .pro file next to the directory which
; contains the images.
;
;
; :Output:
; Three files will be generated in the same directory as this .pro code and they will be
; named: pluto_mosaic, pluto_mosaic.hdr, pluto_mosaic.tif.
;
; The first file is an ENVI format image that can be read into ENVI with the header file (.hdr)
; and the third file is a TIFF images so that you can view the mosaic with a general photo viewer.
; Note that the ENVI header file contains fake projection information. I created a dummy spatial
; reference for each image since there was no metadata provided with the LORRI images.
;
;
;-
;small function for the histogram matching between images
;this is useful because some of the images have bad automatic scaling from NASA and
;they appear very dark
function histogram_matching, match, base
compile_opt idl2
maxval = (max(match))>(max(base))
;calculate the histograms and cumulatide densities
match_histogram = Histogram(match, Min=0, Max=maxval, Binsize=1)
match_cdf = TOTAL(match_histogram, /CUMULATIVE) / N_ELEMENTS(match_histogram)
base_histogram = Histogram(base, Min=0, Max=maxval, Binsize=1)
base_cdf = TOTAL(base_histogram, /CUMULATIVE) / N_ELEMENTS(base_histogram)
;create the transform function z to switch our data space in the match image
z = bytscl(base_histogram)
FOR j=0,n_elements(z)-1 DO BEGIN
i = Where(match_cdf LT base_cdf[j], count)
IF count GT 0 THEN z[j] = i[-1] ELSE z[j]=0
ENDFOR
;do the histogram matching
matched = z[base]
return, matched
end
;method for ENVI to close all data sources and
pro ENVI::reset
compile_opt idl2
e = envi(/current)
if (e eq !NULL) then return
if (e.WIDGET_ID ne 0) then begin
;reset our views
views = e.getview(/all)
foreach view,views do view.close
endif
;close all open rasters
opendata = (e.data).Get(/raster,/vector, /series)
foreach raster, opendata do raster.close
end
pro pluto_mosaicking
compile_opt idl2
;start ENVI and make it a new state!
e = envi(/current)
if (e eq !NULL) then e = envi(/headless) else e.reset
;approximate pixel offset for each image
px_off = [-580, -18]
;directories for running the code
thisdir = FILE_DIRNAME( (((scope_traceback(/STRUCTURE))[-1]).FILENAME))
;imagedir is the location of the original images
imagedir = thisdir + path_sep() + 'images'
;copydir is the location of copies of the images from imagedir that have a dummy spatial reference
copydir = thisdir + path_sep() + 'image_copies'
;warpeddir is where the registered images go
warpeddir = thisdir + path_sep() + 'images_warped'
;if imagedir does not exist, then make it and download the files to that directory
;if the URLs change then these will need to be modified by hand or downloaded by hand from the following URL:
; http://pluto.jhuapl.edu/soc/Pluto-Encounter/index.php
if ~file_test(imagedir) then begin
print, 'Did not find the directory ' + imagedir
print, string(9b) + 'Creating directory and downloading images...'
FILE_MKDIR, imagedir
urls = ['http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179679_0x636_sci_3.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179682_0x636_sci_3.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179685_0x636_sci_2.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179688_0x636_sci_4.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179691_0x636_sci_4.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179694_0x636_sci_4.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179697_0x636_sci_4.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179700_0x636_sci_4.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179703_0x636_sci_4.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179706_0x636_sci_4.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179709_0x636_sci_3.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179712_0x636_sci_3.jpg',$
'http://pluto.jhuapl.edu/soc/Pluto-Encounter/data/pluto/level2/lor/jpeg/029917/lor_0299179715_0x636_sci_3.jpg']
;use the IDLnetURL object to download each image
urlobj = IDLnetURL()
for i=0, n_elements(urls)-1 do begin
outfile = imagedir + path_sep() + strmid(urls[i],strpos(urls[i],'/',/REVERSE_SEARCH)+1)
void = urlobj.Get(FILENAME = outfile,$
URL = urls[i])
print, string(9b) + string(9b) + 'Downloaded image ' + strtrim(i+1,2) + ' of ' + strtrim(n_elements(urls),2)
endfor
urlobj.cleanup
print, string (9b) + 'Downloaded all images!'
endif
;make our output directories if they don't exist
if ~file_test(copydir) then FILE_MKDIR, copydir
if ~file_test(warpeddir) then FILE_MKDIR, warpeddir
;go to imagedir and look for the Pluto images
cd, imagedir, CURRENT = first_dir
files = file_search('*.jpg')
;check to see if we already made copies fo the images
;create a fake spatial reference
psize = 1d*(180/!CONST.R_EARTH);pixel size of 1 meter
spatialRef = ENVIStandardRasterSpatialRef( $
COORD_SYS_CODE=4326, /GEOGCS, $
PIXEL_SIZE=[psize,psize], TIE_POINT_PIXEL=[0.0D,0.0D], $
TIE_POINT_MAP=[0D,0D])
;make sure to sort the files so that we are arranged by time
files = files[sort(files)]
;copy all of the files and create ENVI rasters witha base spatial reference
reference = e.openraster(files[0], SPATIALREF_OVERRIDE=spatialref)
reference_data = reference.GetData()
;export original raster to copy raster
outfile = copydir + path_sep() + file_basename(files[0], '.jpg')
if file_test(outfile) then FILE_DELETE, outfile
reference.export, outfile, 'ENVI'
;perform histogram matching on all images and export them to the copydir
for i=1, n_elements(files)-1 do begin
spatialRef = ENVIStandardRasterSpatialRef( $
COORD_SYS_CODE=4326, /GEOGCS, $
PIXEL_SIZE=[psize,psize], TIE_POINT_PIXEL=[0.0D,0.0D], $
TIE_POINT_MAP=i*psize*px_off)
warp = e.openraster(files[i], SPATIALREF_OVERRIDE=spatialRef)
warp_data = warp.GetData()
warp_match = histogram_matching(reference_data, warp_data)
outfile = copydir + path_sep() + file_basename(files[i], '.jpg')
if file_test(outfile) then FILE_DELETE, outfile
warp_copy = ENVIRaster(warp_data, INHERITS_FROM=warp, URI = outfile)
warp.close
warp_copy.save
warp_copy.close
endfor
reference.close
;go to the copy directory and register our images to one another
cd, copydir
files = file_basename(files, '.jpg')
reference = e.openraster(files[0])
;copy the reference image to the copy directory and open it in ENVI
outfile = warpeddir + path_sep() + files[0]
if file_test(outfile) then FILE_DELETE, outfile
reference.export,outfile, 'ENVI'
reference.close
reference = e.openraster(outfile)
;perform image registration on all images
print, 'Found ' + strtrim(files.length,2) + ' images to register!'
for i=1, n_elements(files)-1 do begin
warp_this = e.openraster(files[i])
; Get the auto tie point generation task from the catalog of ENVITasks
TieTask = ENVITask('GenerateTiePointsByCrossCorrelation')
TieTask.INPUT_RASTER1 = reference
TieTask.INPUT_RASTER2 = warp_this
TieTask.Execute
TieTask_TiePoints = TieTask.OUTPUT_TIEPOINTS
; Get the tie point filter task from the catalog of ENVITasks
FilterTask = ENVITask('FilterTiePointsByGlobalTransform')
FilterTask.INPUT_TIEPOINTS = TieTask_TiePoints
FilterTask.Execute
FilterTask_TiePoints = FilterTask.OUTPUT_TIEPOINTS
; Get the image-to-image registration task from the catalog of ENVITasks
RegistrationTask = ENVITask('ImageToImageRegistration')
RegistrationTask.INPUT_TIEPOINTS = FilterTask_TiePoints
RegistrationTask.WARPING = 'RST'
outfile = warpeddir + path_sep() + files[i]
if file_test(outfile) then FILE_DELETE, outfile
RegistrationTask.OUTPUT_RASTER_URI = warpeddir + path_sep() + files[i]
RegistrationTask.Execute
Reference.close
Reference = RegistrationTask.OUTPUT_RASTER
;close the raster that was warped
warp_this.close
;delete tie points files
FILE_DELETE, TieTask.OUTPUT_TIEPOINTS_URI
FILE_DELETE, FilterTask.OUTPUT_TIEPOINTS_URI
print, string(9b) + 'Completed ' + strtrim(i,2) + ' of ' + strtrim(files.length-1,2) + ' image registrations'
endfor
print, 'Registered all images!'
print, ''
;close the output raster
reference.close
;got to the warped directory and open each warped raster for the mosaic
cd, warpeddir
scenes = objarr(n_elements(files))
for i=0, n_elements(files)-1 do scenes[i] = e.openraster(files[i])
print, 'Generating mosaic...'
; Get the task from the catalog of ENVITasks
Task = ENVITask('BuildMosaicRaster')
; Define inputs
Task.INPUT_RASTERS = scenes
Task.COLOR_MATCHING_METHOD = 'Histogram Matching'
Task.COLOR_MATCHING_STATISTICS = 'Entire Scene'
Task.FEATHERING_METHOD = 'Edge'
Task.FEATHERING_DISTANCE = make_array(n_elements(files), /INTEGER, VALUE = 40)
; Define outputs
outfile = thisdir + path_sep() + 'pluto_mosaic'
if file_test(outfile) then file_delete, outfile
Task.OUTPUT_RASTER_URI = outfile
; Run the task
Task.Execute
;export the raster as a TIF image
outfile = thisdir + path_sep() + 'pluto_mosaic.tif'
if file_test(outfile) then file_delete, outfile
(Task.OUTPUT_RASTER).export, outfile, 'TIFF'
(Task.OUTPUT_RASTER).close
;return to first directory
cd, first_dir
print, 'Mosaic complete!'
;close ENVI
e.close
end
To run the code, simply create a text file called "pluto_mosaicking.pro" and copy/paste the code into the file. When you run the program, IDL will automatically download the Pluto images, create directories, and produce the mosaicked output in the same directory as the pluto_mosaicking.pro file. If there is an issue while downloading the files, you may have to manually download them from NASA at https://pluto.jhuapl.edu/soc/Pluto-Encounter/. Additional instructions are included at the top of the code for how to run the program along with comments throughout the example to explain the steps taken.