; $Id: //depot/idl/IDL_71/idldir/lib/map_grid.pro#1 $ ; ; Copyright (c) 1996-2009, ITT Visual Information Solutions. All ; rights reserved. Unauthorized reproduction is prohibited. ;+ ; NAME: ; FIXED_MAP_GRID ; ; PURPOSE: ; The MAP_GRID procedure draws the graticule of parallels and meridians, ; according to the specifications established by MAP_SET. MAP_SET must be called ; before MAP_GRID to establish the projection type, the center of the ; projection, polar rotation and geographical limits. ; ; CATEGORY: ; Mapping. ; ; CALLING SEQUENCE: ; MAP_GRID ; ; INPUTS: ; NONE ; ; OPTIONAL INPUTS: ; NONE ; ; KEYWORD PARAMETERS: ; ; ; BOX_AXES: Surround the map window with a "box" style axes with ; annotations, outside the box, where the parallels intersect the ; sides, and the meridians intersect the bottom and top edges of the ; box. The border of the box is drawn in alternating foreground and ; background colors, with color changes at each intersection with ; a parallel or meridian. This keyword determines the thickness of ; the box's border, in millimeters. If LABEL is not explicitly ; specified, it defaults to 1 when this keyword is present. If this ; feature is selected, be sure to leave enough room around the map ; window for the annotation, usually by specifying the XMARGIN and ; YMARGIN keywords to MAP_SET. See the example below. ; CHARSIZE: The size of the characters used for the labels. The default is 1. ; COLOR: The color index for the grid lines. ;FILL_HORIZON: Fills the current map_horizon. ; HORIZON: Draws the current map horizon. ; INCREMENT: Determines the spacing between graticle points. ; GLINESTYLE: If set, the line style used to draw the grid of parallels and ; meridians. See the IDL User's Guide for options. The default is ; dotted. ; GLINETHICK: The thickness of the grid lines. Default is 1. ; LABEL: Set this keyword to label the parallels and meridians with their ; corresponding latitudes and longitudes. Setting this keyword to ; an integer will cause every LABEL gridline to be labeled (i.e, ; if LABEL=3 then every third gridline will be labeled). The ; starting point for determining which gridlines are labeled is the ; minimum latitude or longitude (-180 to 180), unless the LATS or ; LONS keyword is set to a single value. In this case, the starting ; point is the value of LATS or LONS. ; LATALIGN: This keyword controls the alignment of the text baseline for ; latitude labels. A value of 0.0 left justifies the label, ; 1.0 right justifies it, and 0.5 centers it. ; LATDEL: The spacing in degrees between parallels of latitude in the grid. ; If not set, a suitable value is determined from the current map ; projection. ; LATS: The desired latitudes to be drawn (and optionally labeled). If ; this keyword is omitted, appropriate latitudes will be generated ; with consideration of the optional LATDEL keyword. If this ; keyword is set to a single value, drawn (and optionally labeled) ; latitudes will be automatically generated as ; [...,LATS-LATDEL,LATS,LATS+LATDEL,...] over the extent of the ; map. The single valued LATS is taken to be the starting point ; for labelling (See the LABEL keyword). ; LATLAB: The longitude at which to place latitude labels. The default is ; the center longitude on the map. ; LATNAMES: An array specifing the names to be used for the latitude labels. ; By default, this array is automatically generated in units of ; degrees. This array can be a string array or numeric, but should ; not be of mixed type. When LATNAMES is specified, LATS must also ; be specified, but the number of elments in the two arrays need not ; be equal. Extra LATNAMES are ignored, while LATNAMES not supplied ; are automatically reported in degrees. LATNAMES can be used when ; LATS is set to a single value. It this case, the LATS used will ; start at the specified latitude and progress northward, wrapping ; around the globe if appropriate. Caution should be used when ; using LATNAMES in conjunction with a single LATS, since the ; number of visible latitude gridlines is dependent on many factors. ; LONALIGN: This keyword controls the alignment of the text baseline for ; longitude labels. A value of 0.0 left justifies the label, ; 1.0 right justifies it, and 0.5 centers it. ; LONDEL: The spacing in degrees between meridians of longitude in the grid. ; If not set, a suitable value is determined from the current map ; projection. ; LONLAB: The latitude at which to place longitude labels. The default is ; the center latitude on the map. ; LONS: The desired longitudes to be drawn (and optionally labeled). If ; this keyword is omitted, appropriate longitudes will be generated ; with consideration of the optional LONDEL keyword. If this ; keyword is set to a single value, drawn (and optionally labeled) ; longitudes will be automatically generated as ; [...,LONS-LONDEL,LONS,LONS+LONDEL,...] over the extent of the map. ; The single valued LONS is taken to be the starting point for ; labeling (See the LABEL keyword). ; LONNAMES: An array specifing the names to be used for the longitude labels. ; By default, this array is automatically generated in units of ; degrees. This array can be a string array or numeric, but should ; not be of mixed type. When LONNAMES is specified, LATS must also ; be specified, but the number of elments in the two arrays need not ; be equal. Extra LONNAMES are ignored, while LATNAMES not supplied ; are automatically reported in degrees. LONNAMES can be used when ; LONS is set to a single value. It this case, the LONS used will ; start at the specified longitude and progress eastward, wrapping ; around the globe if appropriate. Caution should be used when ; using LONNAMES in conjunction with a single LONS, since the number ; of visible longitude gridlines is dependent on many factors. ; ; MAP_STRUCTURE: Set this keyword to a !MAP structure, as returned from ; MAP_PROJ_INIT. If this keyword is set then the gridlines are passed ; through MAP_PROJ_FORWARD and the resulting lines are drawn using ; Cartesian coordinates. If this keyword is not set then it is assumed ; that MAP_SET has been called to set up a map projection. ; ; NO_GRID: Set this keyword if you only want labels but not gridlines. ; ;ORIENTATION: Specifies the clockwise angle in degrees from horizontal ; of the text baseline that the labels should be rotated. Note, ; that the rotation used in MAP_SET is also clockwise, but XYOUTS ; is normally counterclockwise. ; ; T3D: Set this keyword to indicate that the generalized transformation ; matrix in !P.T is to be used. If not present, the user-supplied ; coordinates are simply scaled to screen coordinates. ; ; ZVALUE: Sets the Z coordinate, in normalized coordinates in the ; range of 0 to 1, at which to output the continents. ; ; Note - This keyword has effect only if keyword T3D is set and the ; transformation is stored in !P.T ; ;OUTPUTS: ; Draws gridlines and labels on the current map. ; ; EXAMPLE: ; map_set,10,20,23.5,/cont,/ortho,/horizon ; establish a projection ; lats=[-65,-52,-35,-20,-2.59,15,27.6,35,45,55,75,89] ; ; choose the parallels to draw ; map_grid,lats=lats,londel=20,lons=-13,label=2,lonlab=-2.5,latlab=7 ; ;draw the grid ; Make a map with a grid surrounded by a box style axis: ; map_set, /STEREO, 40, -90,scale=50e6,/CONTINENTS, XMARGIN=3, YMARGIN=3 ; map_grid, /BOX_AXES, COLOR=3, CHARSIZE=1.5 ; ; ; MODIFICATION HISTORY: ; Written by: SVP, May, 23 1996. ; DMS, Feb, 1996 Note that this version plots all gridlines ; unless a 4 element LIMIT was specified to MAP_SET. ; SVP, Nov 1996 Changed over !map1 (new IDL5 maps) ; SVP, May 1996 Map_Grid used to live inside of map_set.pro, now ; it lives here. ; SVP, May 1996 Changed LABEL to determine the skip between labelled ; gridlines. ; Also, added the LATS and LONS keywords to give complete ; control over where the labels go. ; SVP, Sept 1996 Added LATNAMS,LONAMES, and ORIENTATION keywords ; DMS, Nov, 1996 Rev 2 of maps. ; DMS, Jul, 1997 Added Box Axes ; CT, Jan 2004: Added MAP_STRUCTURE keyword. ; Renamed the function FIXED_MAP_GRID with fixed at lines 555-556 and ; 664-665 as described here: http://www.dfanning.com/map_tips/missinggrid.html. ; This is required to get proper grid output in some applications, but ; it breaks other code. ;- function map_grid_check_range, xy ; ; Determine if the xy label point is inside the plot coordinate axes. ; If not, return a single value, which is NOT drawn on the plot. This ; solves the problem of non-Hershey fonts not being clipped to the plot ; area. x = xy[0] y = xy[1] if x lt Min(!X.CRange) then return, xy[0] if x gt Max(!X.CRange) then return, xy[0] if y lt Min(!Y.CRange) then return, xy[0] if y gt Max(!Y.CRange) then return, xy[0] return, xy end function fixed_map_grid_incr, span ; ; Determine LONDEL or LATDEL if not specified ; COMPILE_OPT hidden, IDL2 IF span eq 0 THEN return, 45. ipow = 0 t = abs(span) < 450. WHILE t lt 5 DO BEGIN t = t * 10 & ipow = ipow +1 & ENDWHILE increments = [ 1., 2., 4., 5., 10., 15., 30., 45.] i = 0 WHILE t gt (increments[i] * 10) DO i = i + 1 t = increments[i] / 10^ipow return, span ge 0 ? t : -t end ;------------------------------------------------------------------------- ; Find the point on the line between points c0 and c1, expressed in ; DEVICE coordinates, where the longitude (Icoord = 0) or latitude ; (Icoord = 1) is equal to Gwant. If the segment between c0 and c1 ; doesn't intersect the given meridan or parallel, or either endpoint ; is not mappable, return NaN. ; Otherwise, return the device coordinate, X if Icoord = 0, or Y if ; Icoord = 1, of the intersection. ; Function fixed_map_grid_SOLVE, c0, c1, Icoord, Gwant, $ MAP_STRUCTURE=mapStruct, RECURSIVE=recursive COMPILE_OPT hidden, IDL2 hasMap = N_TAGS(mapStruct) gt 0 p0 = CONVERT_COORD(c0, /DEVICE, /TO_DATA) p1 = CONVERT_COORD(c1, /DEVICE, /TO_DATA) if (hasMap) then begin ; Convert from UV to latlon p0 = MAP_PROJ_INVERSE(p0[0], p0[1], MAP_STRUCTURE=mapStruct) p1 = MAP_PROJ_INVERSE(p1[0], p1[1], MAP_STRUCTURE=mapStruct) endif p0 = p0[Icoord] p1 = p1[Icoord] ; Not mappable or zero interval. if ~finite(p0) || ~finite(p1) || (p1 eq p0) then $ return, !values.f_nan if (Icoord eq 0) && (p0 gt p1) then begin ;Wrap if we cross dateline if gwant ge 0 then p1 += 360. $ else p0 -= 360. endif t = (Gwant - p0) / (p1-p0) if (t lt 0.0) || (t gt 1.0) then $ return, !values.f_nan low = 0.0 high = 1.0 tol = 1.0e-5 del = c1 - c0 while abs(high-low) gt tol do begin ;Binary chopping method t = (low + high) / 2. c = c0 + t * del p = CONVERT_COORD(c, /DEVICE, /TO_DATA) if (hasMap) then $ ; Convert from UV to latlon p = MAP_PROJ_INVERSE(p[0], p[1], MAP_STRUCTURE=mapStruct) p = p[Icoord] if (~FINITE(p)) then $ return, !values.f_nan if (Icoord eq 0) then begin ;Check for dateline? if p lt p0 then p += 360. $ ;Wrap? P should be in interval p0-p1. else if p gt p1 then p -= 360. endif if (Gwant-p0) * (Gwant - p) gt 0.0 then begin ;In same interval as p0 : low low = t p0 = p endif else high = t ;Keep low, and fcn at low = p0. endwhile return, c[Icoord] end ;------------------------------------------------------------------------- PRO fixed_map_grid, LABEL=label, LATDEL = latdel, NO_GRID=no_grid, $ LONDEL=londel, GLINESTYLE=glinestyle, GLINETHICK=glinethick, $ LONLAB=lonlab, LATLAB=latlab, LONALIGN=lonalign, $ LATALIGN=latalign, LONNAMES=lonnames, LATNAMES=latnames, $ LATS=lats, LONS=lons, ZVALUE=zvalue, $ COLOR=color, T3D=t3d, CHARSIZE=charsize, ORIENTATION=orientation, $ HORIZON=horizon, FILL_HORIZON=fill_horizon, _EXTRA=extra, $ INCREMENT=increment, CLIP_TEXT=clip_text, BOX_AXES=box_axes, $ MAP_STRUCTURE=mapStruct, FORMAT=format, $ WHOLE_MAP=obsolete_keyword compile_opt idl2 ON_ERROR, 2 hasMap = N_TAGS(mapStruct) gt 0 if ((!x.type NE 3) && ~hasMap) THEN $ message, 'fixed_map_grid---Current ploting device must have mapping coordinates' ; Put a grid on a previously established map projection. ; ; no grid? - in case someone wants just to put labels no_grid = keyword_set(no_grid) ; if Label = n, then Labels are added every n gridlines ; If box_axes is set, and LABEL isn't explicitly specified, set label. ; nlabel = (N_ELEMENTS(label) gt 0) ? $ FIX(ABS(label[0])) : KEYWORD_SET(box_axes) have_lons = n_elements(lons) gt 0 have_lats = n_elements(lats) gt 0 if n_elements(zvalue) eq 0 THEN zvalue = 0 ; CLIP_TEXT (default value = 1) = 1 to clip text within the map area, ; 0 to not clip text. noclip = (N_ELEMENTS(clip_text) gt 0) ? ~KEYWORD_SET(clip_text) : 0 ; Append the graphics keywords: if n_elements(t3d) then map_struct_append, extra,'T3D',t3d if n_elements(color) then map_struct_append, extra, 'COLOR',color if n_elements(glinestyle) eq 0 then glinestyle = 1 ;Default = dotted map_struct_append, extra, 'LINESTYLE', glinestyle ;Append it if n_elements(glinethick) then map_struct_append, extra,'THICK',glinethick if n_elements(charsize) then map_struct_append, extra,'CHARSIZE', charsize ;Orientation is reversed & conflicts w/box_axes if n_elements(orientation) and (keyword_set(box_axes) eq 0) then $ map_struct_append, extra,'ORIENTATION', -1 * orientation ; ; The gridlines can be specified by ; ; 1) an array of lats and/or lons ; 2) a single lats or lons which is taken to be the center ; or 'for sure' lat or lon with gridlines every latdel or londel from it ; 3) automatically calculated if lats or lons are not specified. ; ; ; Require that LATS be specified when LATNAMES is ALSO SPECIFIED ; if (n_elements(latnames) gt 0) and n_elements(lats) le 1 then $ message,'fixed_map_grid---The LATNAMES keyword MUST be used in conjuction '+$ 'with the LATS keyword.' if n_elements(lonnames) gt 0 and have_lons eq 0 then $ message,'fixed_map_grid---The LONNAMES keyword MUST be used in conjuction '+$ 'with the LONS keyword.' ; Get lat/lon ranges from !MAP. Did MAP_SET specify 4 element limit? map = hasMap ? mapStruct : !MAP if n_elements(lats) gt 1 then latmin = min(lats, max=latmax) $ else if map.ll_box[0] ne map.ll_box[2] then begin latmin = map.ll_box[0] latmax = map.ll_box[2] endif else begin latmin = -90 latmax = 90 endelse if have_lons then begin ;Lons directly specified? lonmin = lons[0] lonmax = lons[n_elements(lons)-1] endif else if (map.ll_box[1] ne map.ll_box[3]) and $ ;Lon limit specified? (latmax lt 90.) and (latmin gt -90.) then begin ; and poles not included lonmin = map.ll_box[1] lonmax = map.ll_box[3] ;Copy limits endif else begin ;If not, use entire globe lonmin = -180 lonmax = 180 endelse IF lonmax le lonmin THEN lonmax = lonmax + 360. ;Default grid spacings... IF n_elements(latdel) eq 0 THEN begin latdel = fixed_map_grid_incr(latmax - latmin) latd = 1 endif else latd = latdel IF n_elements(londel) eq 0 THEN begin londel = fixed_map_grid_incr(lonmax - lonmin) lond = 1 endif else lond = londel ; IF the deltas are < 1, ; do not convert the limits into integers IF abs(latmax - latmin) gt 5. and latd ge 1 THEN BEGIN ;Make range integers latmin = float(floor(latmin)) latmax = ceil(latmax) ENDIF IF abs(lonmax - lonmin) gt 5 and lond ge 1 THEN BEGIN ;Integerize long spans lonmin = float(floor(lonmin)) lonmax = ceil(lonmax) ENDIF ; Where we label things... IF N_Elements(Latlab) eq 0 THEN Latlab = (lonmin + lonmax)/2 IF N_ELements(LonLab) eq 0 THEN LonLab = (latmin +latmax)/2 IF n_elements(latalign) eq 0 THEN latalign = .5 ;Text alignment of lat labels IF n_elements(lonalign) eq 0 THEN lonalign = .5 ;Text alignment of lon labels ; Is this a cylindrical proj? if (hasMap) then $ is_cyl = 0 $ else $ map_proj_info, iproj, CYLINDRICAL=is_cyl, /CURRENT if keyword_set(increment) then step = increment $ else step = 4 < (latmax - latmin)/10. len = long(float((latmax-latmin)) / float(step) + 1.0) ; Clip to avoid roundoff errors which can cause the latitude to exceed ; 90 degs by a very small amount. lati = (float(latmax-latmin) / (len-1.)) * findgen(len) + latmin > (-90) < 90 ; This fudge avoids curved meridians at the poles because of the split planes if is_cyl and map.p0lat eq 0 then begin del = 2.0e-2 if lati[0] eq -90 then lati[0] = del - 90. if lati[len-1] eq 90 then lati[len-1] = 90. - del endif ;Compute longit distance between points for latitude lines. step = 4 < (lonmax - lonmin)/10. ;At most 4 degrees len = (lonmax-lonmin)/step + 1 loni = findgen(len) * step + lonmin IF (loni[len-1] NE lonmax) THEN loni = [loni, lonmax] ; ; Determine the number of lons and the lon array ; if n_elements(lons) eq 0 then begin n_lons = 1+fix((lonmax-lonmin) / londel) longitudes = lonmin - (lonmin mod londel) + findgen(n_lons) * londel endif else if n_elements(lons) eq 1 then begin i0 = ceil((lonmin - lons[0]) / float(londel)) ;First tick i1 = floor((lonmax - lons[0]) / float(londel)) ;Last tick n_lons = i1 - i0 + 1 > 1 longitudes = (findgen(n_lons) + i0) * londel + lons[0] endif else begin n_lons=n_elements(lons) longitudes=lons endelse ; ; Determine the number of lats and the lat array ; if n_elements(lats) eq 0 then begin lat0 = latmin - (latmin mod float(latdel)) ;1st lat for grid n_lats = 1 + fix((latmax-lat0)/float(latdel)) latitudes = lat0 + findgen(n_lats)*latdel endif else if n_elements(lats) eq 1 then begin i0 = ceil((latmin - lats[0]) / float(latdel)) ;First tick i1 = floor((latmax - lats[0]) / float(latdel)) ;Last tick n_lats = i1 - i0 + 1 > 1 latitudes = (findgen(n_lats) + i0) * latdel + lats[0] endif else begin n_lats=n_elements(lats) latitudes=lats endelse ; ; Build the Latitude/Longitude Label Flags ; lon_label = bytarr(n_lons) lat_label = bytarr(n_lats) if nlabel ne 0 then begin if n_elements(lons) eq 1 then begin ; Ensure center is set and then go out index=where(longitudes eq lons[0]) for i=(index[0] > 0) mod nlabel, n_lons-1, nlabel do lon_label[i] = 1 endif else begin for i=0, n_lons-1, nlabel do lon_label[i] = 1 endelse if n_elements(lats) eq 1 then begin ; Make sure the center one is set ; and go out from there index=where(latitudes eq lats[0], count) for i=(index[0] > 0) mod nlabel, n_lats-1, nlabel do lat_label[i] = 1 endif else begin ; Start with latmin and label each nlabel point for i=0, n_lats-1, nlabel do lat_label[i] = 1 endelse endif ; ; Dont repeat 180 labelling if the projection is cylindrical or ; polar and both 180 and -180 are present. This can be defeated by using ; LONS=-180 ; if is_cyl or (abs(map.p0lat) eq 90) then begin id_180 = where(longitudes eq 180,count) id_m180 = where(longitudes eq -180,mcount) if count gt 0 and mcount gt 0 then begin if n_elements(lons) eq 1 then begin if lons[0] eq -180 then lon_label[id_180]=0 endif else lon_label[id_m180]=0 endif endif n = n_lons > n_lats ; latlontxt = strarr(n, 2) if keyword_set(box_axes) then begin ;Draw a Box legend? box_thick = box_axes * 0.1 ;From mm to Thickness in cm dc = !d.y_ch_size ;Char height to draw if n_elements(charsize) then dc = dc * charsize xw = !x.window * !d.x_size ;Window coords in x & y yw = !y.window * !d.y_size ; xww and yww = corners of the uv_range that is mappable. If NOBORDER ; was set for MAP_SET, this is the same as the window coords (xw,yw), ; otherwise, this rectangle is smaller than the window rectangle. ; Fudge factor for window to ensure that the edges are mappable. del = [1,-1]* 0.01 xww = (map.uv_box[[0,2]] * !x.s[1] + !x.s[0]) * !d.x_size + del yww = (map.uv_box[[1,3]] * !y.s[1] + !y.s[0]) * !d.y_size + del bdel = box_thick * !d.y_px_cm ;Thickness of box in device units if n_elements(color) eq 0 then bcolor = !p.color $ ;Box color else bcolor = color xp = xw[0] - [0,bdel, bdel,0] ;X & Y polygon coords for outer box yp = yw[0] - [0,0,bdel,bdel] ;Draw the outline of the box plots, xw[[0,1,1,0,0]], yw[[0,0,1,1,0]], /DEVICE, COLOR=bcolor plots, xw[[0,1,1,0,0]]+[-bdel, bdel, bdel, -bdel, -bdel], $ yw[[0,0,1,1,0]]+[-bdel, -bdel, bdel, bdel, -bdel], /DEVICE, COLOR=bcolor ychar = [yw[0]-bdel-dc, yw[1]+bdel+dc/4.] xchar = [xw[0] - bdel - dc/4., xw[1]+bdel+dc/4.] boxpos = replicate(!values.f_nan, n, 2,2) ;Device coordinates for box annotations. Go in to avoid edges of map & ;border which are frought with singularities. Also to avoid effects ;of MAP_SET,/NOBORDER. For box axes to be annotated, all the edges of the ;map rectangle must be mappable. ; endif else box_thick = 0 ; Do the horizon if specified. if keyword_set(horizon) then map_horizon, _EXTRA=e if keyword_set(fill_horizon) then map_horizon, _EXTRA=e, /FILL ; ; ****************** Draw/Label the meridians ****************** ; FOR i=0,n_lons-1 DO BEGIN lon=longitudes[i] lon2 = (lon lt -180) ? (lon + 360) : $ ((lon gt 180) ? (lon - 360) : lon) ; This block of code draws longitude lines that are at the center + or ; - 180 degrees, at center + or - (180-eps) to ensure that the grid ; appears on the correct side. Its really not necessary if people ; would use the /HORIZON keyword, but they don't. if is_cyl then begin del = lon - map.p0lon while del gt 180 do del -= 360. while del lt -180 do del += 360. if abs(del) eq 180 then $ lon -= 1.0e-5 * ((del ge 0) ? 1 : -1) ;fudge it (sign(1.0e-5, del)) endif IF (~no_grid) THEN begin if (hasMap) then begin ; Make sure to clear out variable in case all points are clipped. polylines = -1 ; uv = MAP_PROJ_FORWARD(REPLICATE(lon, N_ELEMENTS(lati)), lati, $ ; MAP_STRUCTURE=mapStruct, POLYLINES=polylines) uv = MAP_PROJ_FORWARD(REPLICATE(lon, N_ELEMENTS(lati)), lati, $ MAP_STRUCTURE=mapStruct) polylines = [N_Elements(lati), Indgen(N_Elements(lati))] index = 0L npoly = N_ELEMENTS(polylines) ; Loop thru our polylines connectivity array. while (index lt npoly) do begin nline = polylines[index] if (nline eq -1) then $ break if (nline gt 0) then begin indices = polylines[index + 1 : index + nline] PLOTS, REFORM(uv[0,indices]), REFORM(uv[1,indices]), $ zvalue, $ NOCLIP=1, $ _EXTRA=extra endif index += nline + 1 endwhile endif else begin PLOTS, lon, lati, zvalue, NOCLIP=0, _EXTRA=extra endelse endif IF N_Elements(format) EQ 0 THEN BEGIN fmt = (lon2 ne long(lon2)) ? '(f7.2)' : '(i4)' ENDIF ELSE BEGIN IF format EQ "" THEN fmt = (lon2 ne long(lon2)) ? '(f7.2)' : '(i4)' ELSE fmt = format ENDELSE IF lon_label[i] THEN BEGIN IF i lt n_elements(lonnames) then begin ;User specified label? IF (reverse(size(lonnames[i])))[1] eq 7 then $ ;String? lonname=lonnames[i] else $ lonname=strtrim(string(lonnames[i], FORMAT=fmt),2) endif else $ lonname=strtrim(string(lon2, format=fmt),2) latlontxt[i,0] = lonname if (box_thick eq 0) then begin xy = 0 if (hasMap) then begin ; We need to convert from latlon to UV ourself. uv = MAP_PROJ_FORWARD(lon, LonLab, MAP_STRUCTURE=mapStruct) if (FINITE(uv[0]) && FINITE(uv[1])) then begin xy = uv[0:1] xy = map_grid_check_range(xy) endif endif else begin if (noclip eq 1) || map_point_valid(lon, LonLab) then $ xy = [lon, LonLab] endelse if (N_ELEMENTS(xy) eq 2) then $ XYOUTS, xy[0], xy[1], lonname, ALIGNMENT=lonalign, $ NOCLIP=noclip, Z=zvalue, _EXTRA=extra endif ENDIF if box_thick ne 0 then begin dy = (yw[1] - yw[0]) * 0.01 ;1% of the height for j=0,1 do begin ;Save longitude crossings, try for edge... k = 0 ; Try to find longitude crossings. If it doesn't cross exactly ; at the edge, try going in until it crosses and is valid. while ~finite(boxpos[i,j,0]) && abs(k) lt 3 do begin boxpos[i, j, 0] = fixed_map_grid_solve( $ [xww[0], yw[j]+k*dy], $ [xww[1], yw[j]+k*dy], 0, lon, $ MAP_STRUCTURE=mapStruct) k = k + (j ? -1 : 1) endwhile endfor endif ENDFOR ; ; Draw/Label the parallels of latitude ****************** ; FOR i=0,n_lats-1 DO BEGIN lat=latitudes[i] IF N_Elements(format) EQ 0 THEN BEGIN fmt = (lat ne long(lat)) ? '(f7.2)' : '(i4)' ENDIF ELSE BEGIN IF format EQ "" THEN fmt = (lat ne long(lat)) ? '(f7.2)' : '(i4)' ELSE fmt = format ENDELSE IF lat_label[i] THEN BEGIN IF i lt n_elements(latnames) then begin ;User specified latname? IF (reverse(size(latnames[i])))[1] eq 7 then $ latname=latnames[i] else $ latname=strtrim(string(latnames[i],format=fmt),2) endif else latname=strtrim(string(lat, format=fmt),2) latlontxt[i, 1] = latname if (box_thick eq 0) then begin xy = 0 if (hasMap) then begin ; We need to convert from latlon to UV ourself. uv = MAP_PROJ_FORWARD(latlab, lat, MAP_STRUCTURE=mapStruct) if (FINITE(uv[0]) && FINITE(uv[1])) then begin xy = uv[0:1] xy = map_grid_check_range(xy) endif endif else begin if (noclip eq 1) || map_point_valid(latlab, lat) then $ xy = [latlab, lat] endelse if (N_ELEMENTS(xy) eq 2) then $ XYOUTS, xy[0], xy[1], latname, $ alignment=latalign, NOCLIP=noclip, Z=zvalue, _EXTRA=extra endif ENDIF if (~no_grid && (ABS(lat) ne 90)) then begin if (hasMap) then begin ; Make sure to clear out variable in case all points are clipped. polylines = -1 ; uv = MAP_PROJ_FORWARD(loni, REPLICATE(lat, N_ELEMENTS(loni)), $ ; MAP_STRUCTURE=mapStruct, POLYLINES=polylines) uv = MAP_PROJ_FORWARD(loni, REPLICATE(lat, N_ELEMENTS(loni)), $ MAP_STRUCTURE=mapStruct) polylines = [N_Elements(loni), Indgen(N_Elements(loni))] index = 0L npoly = N_ELEMENTS(polylines) ; Loop thru our polylines connectivity array. while (index lt npoly) do begin nline = polylines[index] if (nline eq -1) then $ break if (nline gt 0) then begin indices = polylines[index + 1 : index + nline] PLOTS, REFORM(uv[0,indices]), REFORM(uv[1,indices]), $ zvalue, $ NOCLIP=0, $ _EXTRA=extra endif index += nline + 1 endwhile endif else $ PLOTS, loni, lat, zvalue, NOCLIP=0, _EXTRA=extra endif if box_thick ne 0 then begin for j=0,1 do begin ;Save latitude crossings ; Start at edge and try for an intersection. ; If that doesn't work, go in some. k = 0 dx = (xw[1] - xw[0]) * 0.01 while ~finite(boxpos[i,j,1]) && abs(k) lt 3 do begin boxpos[i, j, 1] = fixed_map_grid_solve( $ [xw[j]+dx*k, yww[0]], $ [xw[j]+dx*k, yww[1]], 1, lat, $ MAP_STRUCTURE=mapStruct) k = k + (j ? -1 : 1) endwhile endfor endif endfor ; ******************************** Do the box axes ************************** if box_thick ne 0 then begin for iaxis=0,1 do for j=0,1 do begin v = boxpos[*,j,iaxis] ;Values along axes good = where(finite(v), count) ;Ignore bad values if (count eq 0) then $ ;Anything there? continue dy = iaxis eq 1 v = v[good] ;Remove unmappable elements subs = sort(v) ;Sort the axis crossings v = v[subs] ;Sorted Y values vtext = (latlontxt[good,iaxis])[subs] v0 = ([xw[0], yw[0]])[iaxis] ;Starting value on axis xp0 = xp + j * (xw[1]-xw[0] + bdel) ;Polygon X coords yp0 = yp + j * (yw[1]-yw[0] + bdel) ;Y coords xychar = [xchar[j], ychar[j]] ;Char position for i=0, count-1 do begin ;Draw each item z = v[i] ;Axis crossing value if iaxis eq 0 then begin xp0 = (i eq (count-1) && (i and 1) && ~vtext[i]) ? $ [v0, xw[1], xw[1], v0] : [v0, z, z, v0] endif else yp0 = [v0, v0, z, z] if (i and 1) then $ polyfill, xp0, yp0, /DEVICE, COLOR=bcolor xychar[iaxis] = z if strlen(vtext[i]) gt 0 then begin xyouts, xychar[0], xychar[1], vtext[i], $ ORIENTATION=dy * (90-180*j), $ ALIGN=0.5, CLIP=0, Z=zvalue, /DEVICE, _EXTRA=extra endif v0 = z endfor ;Fill to the end of the axis if i and 1 then begin if iaxis eq 0 then xp0 = [v0, xw[1], xw[1], v0] $ else yp0 = [v0, v0, yw[1], yw[1]] polyfill, xp0, yp0, /DEVICE, COLOR=bcolor endif endfor endif ; box_thick end