; docformat = 'rst' ; ; NAME: ; cgChangeMapProjection ; ; PURPOSE: ; This function warps a map projected image from one map projection to another, using ; Map_Proj_Image to do the warping. ; ;******************************************************************************************; ; ; ; Copyright (c) 2012, by Fanning Software Consulting, Inc. All rights reserved. ; ; ; ; Redistribution and use in source and binary forms, with or without ; ; modification, are permitted provided that the following conditions are met: ; ; ; ; * Redistributions of source code must retain the above copyright ; ; notice, this list of conditions and the following disclaimer. ; ; * Redistributions in binary form must reproduce the above copyright ; ; notice, this list of conditions and the following disclaimer in the ; ; documentation and/or other materials provided with the distribution. ; ; * Neither the name of Fanning Software Consulting, Inc. nor the names of its ; ; contributors may be used to endorse or promote products derived from this ; ; software without specific prior written permission. ; ; ; ; THIS SOFTWARE IS PROVIDED BY FANNING SOFTWARE CONSULTING, INC. ''AS IS'' AND ANY ; ; EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED WARRANTIES ; ; OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT ; ; SHALL FANNING SOFTWARE CONSULTING, INC. BE LIABLE FOR ANY DIRECT, INDIRECT, ; ; INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED ; ; TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR PROFITS; ; ; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ; ; ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT ; ; (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS ; ; SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. ; ;******************************************************************************************; ; ;+ ; This function warps a map projected image from one map projection to another, using ; Map_Proj_Image to do the warping. Useful in general, it is used specifically to ; warp map projected images into the Equirectangular map projection with a WGS 84 ellipsoid ; that is used by Google Earth. See the program `cgImage2KML` to see how it is used. ; ; :Categories: ; Map Projections, Utilities ; ; :Returns: ; An image warped into the output map projection is returned. ; ; :Params: ; image: in, required ; Any 2D or true-color image with or without an alpha channel. This image will ; be warped to the output map projection (`mapOut`). ; mapin: in, required, type=object ; A map coordinate object (cgMap) describing the map projection and coordinates of the ; input image that is to be warped. ; ; :Keywords: ; boundary: in, optional, type=dblarr ; A four-element array specifying the Cartesian (XY) coordinates of the input ; image range, in the form [xmin, ymin, xmax, ymax]. If this parameter is not ; present, the boundary will be obtained from the `MapIn` map coordinate object. ; On output, this variable will contain the Cartesian boundary of the warped image. ; bilinear: in, optional, type=boolean, default=0 ; Set this keyword to warp the image with bilinear interpolation. The default is ; to do nearest neighbor interpolation. ; latlonbox: out, optional, type=fltarr ; A four-element float array containing the boundaries of the warped image in ; the [north, south, east, west] form preferred by Google Earth. Listed as degrees. ; mapout: in, optional, type=object ; A map coordinate object (cgMap) describing the map projection and coordinates of the ; warped image. This image will be warped into this map projection. The [XY]Range of ; this object will be set to the `XYRange` of the output image. ; mask: out, optional ; Set this keyword equal to a named variable that will contain a byte array of the same ; dimensions as the output image, containing a mask of the “good” values. Values in the output ; image that were set to `Missing` (that is, the values were off the map) will have a mask value ; of 0, while all other mask values will be 1. ; missing: in, optional, type=integer, default=0 ; Set this keyword equal to an integer value to be used for pixels that fall outside of ; the valid map coordinates. The default value is 0. ; xyrange: out, optional, type=dblarr ; The Cartesian (XY) coordinates associated with the output image. These are the map image ; boundaries of the output image. ; ; :Examples: ; To prepare a GeoTiff file for creating a KML overlay on Google Earth:: ; netObject = Obj_New('IDLnetURL') ; url = 'http://www.idlcoyote.com/data/AF03sep15b.n16-VIg.tif' ; returnName = netObject -> Get(URL=url, FILENAME=AF03sep15b.n16-VIg.tif') ; Obj_Destroy, netObject ; map = cgGeoMap('AF03sep15b.n16-VIg.tif') ; googleMap = Obj_New('cgMap', 'Equirectangular', Ellipsoid='WGS 84') ; warpedImage = cgChangeMapProjection(image, map, MAPOUT=googleMap) ; ; :Author: ; FANNING SOFTWARE CONSULTING:: ; David W. Fanning ; 1645 Sheely Drive ; Fort Collins, CO 80526 USA ; Phone: 970-221-0438 ; E-mail: david@idlcoyote.com ; Coyote's Guide to IDL Programming: http://www.idlcoyote.com ; ; :History: ; Change History:: ; Written, 30 October 2012 by David W. Fanning. ; Fixed a problem with a TRANSPOSE command for true-color images. Bad logic. 4 Jan 2012. DWF. ; Fixed a problem in which Map_Proj_Image was not handling true-color image warping correctly ; and also a problem with handling missing values passed from cgImage2KML correctly. 20 Feb 2013. DWF. ; Forgot to specify the default value for MISSING. 25 Feb 2014. DWF. ; ; :Copyright: ; Copyright (c) 2012, Fanning Software Consulting, Inc. ;- FUNCTION cgChangeMapProjection, image, mapIn, $ BOUNDARY=boundary, $ BILINEAR=bilinear, $ LATLONBOX=latlonbox, $ MAPOUT=mapOut, $ MASK=mask, $ MISSING=missing, $ XYRANGE=xyrange Compile_Opt idl2 ; Error handling. ON_Error, 2 ; An image and an input map coordinate object are required. IF N_Elements(image) EQ 0 THEN Message, 'An input image is required.' IF N_Elements(mapIn) EQ 0 THEN Message, 'An input map coordinate object is required.' ; Assign default values. IF N_Elements(boundary) EQ 0 THEN mapIn -> GetProperty, BOUNDARY=boundary bilinear = Keyword_Set(bilinear) IF N_Elements(mapOut) EQ 0 THEN BEGIN mapOut = Obj_New('cgMap', 'Equirectangular', Ellipsoid='WGS 84') ENDIF ; If the missing value is a color, change it to a color triple for true-color images. IF N_Elements(missing) NE 0 THEN BEGIN IF (Size(image, /N_Dimensions) GT 2) THEN BEGIN IF Size(missing, /TNAME) EQ 'STRING' THEN missing = Reform(cgColor(missing, /Triple)) ENDIF ENDIF ELSE BEGIN IF (Size(image, /N_Dimensions) GT 2) THEN missing = BytArr(3) ELSE missing = 0B ENDELSE ; Gather information about the input image. Do this differently if this is a 2D ; image as opposed to a true-color image. dims = Image_Dimensions(image, XIndex=xindex, YIndex=yindex, TRUEINDEX=trueindex, $ XSize=xsize, YSize=ysize, ALPHACHANNEL=alphachannel) IF trueindex EQ -1 THEN BEGIN warped = Map_Proj_Image(image, boundary, Image_Structure=mapIn->GetMapStruct(), $ Map_Structure=mapOut->GetMapStruct(), UVRANGE=xyrange, Missing=missing[0], MASK=mask) ; Update the range values in the output map coordinate object. mapOut -> SetProperty, XRANGE=xyrange[[0,2]], YRANGE=xyrange[[1,3]] ; Return information of interest to user. mapOut -> GetProperty, LATLONBOX=latlonbox, BOUNDARY=boundary ENDIF ELSE BEGIN IF trueIndex LT 2 THEN img = Transpose(image, [xindex, yindex, trueindex]) frames = alphachannel ? 4 : 3 warped = Make_Array(xsize, ysize, frames, TYPE=Size(image, /TYPE)) mapInStruct = mapIn->GetMapStruct() mapOutStruct = mapOut->GetMapStruct() ; Warp the three image planes separately. FOR j=0,frames-1 DO BEGIN warped[*,*,j] = Map_Proj_Image(img[*,*,j], boundary, Image_Structure=mapInStruct, $ Map_Structure=mapOutStruct, UVRANGE=xyrange, Missing=missing[j], Mask=mask) ENDFOR ; Handle alpha channel separately. IF frames EQ 4 THEN BEGIN warped[*,*,3] = Map_Proj_Image(img[*,*,3], boundary, Image_Structure=mapInStruct, $ Map_Structure=mapOutStruct, UVRANGE=xyrange) ENDIF ; Update the range values in the output map coordinate object. mapOut -> SetProperty, XRANGE=xyrange[[0,2]], YRANGE=xyrange[[1,3]] ; Return information of interest to user. mapOut -> GetProperty, LATLONBOX=latlonbox, BOUNDARY=boundary ENDELSE RETURN, warped END