;------------------------------------------------------------- ;+ ; NAME: ; GAMAP ; ; PURPOSE: ; Menu-driven user interface for creating plots with ; the GAMAP package subroutines. The actual data retrieval ; and plotting is donre with ctm_plot.pro. This routine ; mainly collects all user requests and passes them on to ; ctm_plot. ; ; CATEGORY: ; Plotting / CTM Tools ; ; CALLING SEQUENCE: ; GAMAP, [ DiagN [, Keywords ] ] ; ; INPUTS: ; DIAGN -> A diagnostic number or category to restrict ; the record selection (default is: use all). ; ; KEYWORD PARAMETERS: ; General keywords: ; ----------------------------------------------------------------- ; FILENAME -> CTM output file name. Default is to display a ; pickfile dialog and let the user select. You can have ; wildcards ('*', '?') in your filename which restricts ; the file selection. ; ; /NOFILE -> Don't query for filename but display all records that ; have already been loaded. This can save you a couple of ; mouse clicks when you want to create several plots with ; data from one file, and it also useful when you want ; to plot data from 'external' files that were converted ; with ctm_make_datainfo. If a filename is given or no ; data were loaded, the file selection dialog will appear ; anyhow. ; ; /RESETDEFAULTS -> If set, will reset all GAMAP values to their ; defaults. ; ; /HELP -> Displays a help page. ; ; Keywords to restrict the number of records displayed for selection: ; ----------------------------------------------------------------- ; TRACER -> A tracer number to restrict record selection ; ; TAU0 -> Time value (at beginning of record) ; ; DATE -> 6 digit date (e.g. 940101) at the beginning of the ; output record (this gets translated into a TAU0 ; value via the function nymd2tau). You can specify more ; than one date at a time using a vector (e.g. [940101, 940301]). ; For the GISS model(s), you also have to specify /GISS in ; order to get correct tau values. The time is assumed to ; be 00 GMT. For other times use the TAU0 keyword as ; TAU0=nymd2tau(dates,times). ; ; Keywords defining output options (these override defaults in ; gamap.defaults) ; ----------------------------------------------------------------- ; /PS -> If set, will directly send output to the 'idl.ps' file. ; If not set, GAMAP will prompt the user whether to create ; the 'idl.ps' file. ; ; OUTFILENAME -> Name of file to send PostScript output to. ; ; /NOTIMESTAMP -> Do not include a user ID and time stamp ; on the postscript plot. Unnecessary if the TIMESTAMP value ; in gamap.defaults is set to 0. ; ; /DO_GIF -> If set, GAMAP will save animation frames as GIF ; files. If not set, GAMAP will prompt the user whether ; to save animation frames to GIF files. DO_GIF overrides ; the default setting in "gamap.defaults". ; ; GIFFILENAME -> Name of GIF file to save animation frames to. ; If the token %N% is used in GIFFILENAME, then GAMAP ; will replace %N% with the actual frame number. If ; GIFFILENAME is not set, or if DO_GIF is set to *QUERY in ; "gamap.defaults", GAMAP will prompt user to supply ; GIFFILENAME. ; ; /DO_MPEG -> If set, GAMAP will save animation frames as GIF ; files. If not set, GAMAP will prompt the user whether ; to save animation frames to GIF files. DO_MPEG overrides ; the default setting in "gamap.defaults". ; ; MPEGFILENAME -> Name of MPEG file to save animation frames to. ; If MPEGFILENAME is not set, or if DO_MPEG is set to *QUERY in ; "gamap.defaults", GAMAP will prompt user to supply ; MPEGFILENAME. ; ; _EXTRA=e -> Picks up extra keywords for CTM_PLOT, etc... ; ; OUTPUTS: ; ; SUBROUTINES: ; Internal subroutines: ; ------------------------------------------------------ ; GAMAP_CheckDataBlockConsistency (function) ; GAMAP_FindNearestEdges (function) ; GAMAP_GetDataBlockRanges (function) ; GAMAP_GetDefaultRanges ; GAMAP_AutoYRange (function) ; GAMAP_PrintDimInfo ; GAMAP_QueryAnimationOptions ; GAMAP_QueryAverageOrTotal ; GAMAP_QueryPostScriptOptions ; GAMAP_SelectDataBlocks (function) ; GAMAP_SelectPlotType ; GAMAP_StoreGridInfo ; GAMAP_UserRangeEntry (function) ; ; ; Also uses external subroutines: ; ------------------------------------------------------ ; CHOICE (function) CLOSE_DEVICE ; CTM_PLOT OPEN_DEVICE ; STRREPL (function) STRWHERE (function) ; DEFAULT_RANGE (function) CHKSTRU (function) ; REPLACE_TOKEN (function) CTM_GRID (function) ; MAKE_SELECTION(function) ; ; REQUIREMENTS: ; Uses GAMAP package subroutines. ; ; NOTES: ; (1) GAMAP can read ASCII punch files with GEOS or GISS, model II ; diagnostics, binary punch files (as defined in Jan 1999), ; and GEOS-CTM binary restart files. Binary punch files ; are processed much faster and allow "windowing" of output. ; ; (2) For pixel plots, GAMAP can only plot cylindrical maps, ; with rectangular projections. Arbitrary map projections ; are possible with any type of contour plot ( *** hasn't been ; tested ***). ; ; (3) GAMAP forces map ranges to coincide with the grid box ; edges, so that the map and pixel plot will be aligned. ; Each "pixel" size corresponds to one full grid box. ; For grids with half-polar boxes, it is therefore recommended ; not to plot the polar latitudes, since those boxes will ; show up as full-size boxes and shift the rest of the plot ; accordingly. ; ; (4) When the user selects multiple data blocks, GAMAP will produce ; a multi-panel plot if !p.multi indicates more than one panel ; on the screen (use the MULTIPANEL procedure to turn it on). ; If you plot only one panel per screen, GAMAP will automatically ; start XInterAnimate to present your own movie to you. Be ; aware that XInterAnimate is limited by your system resources. ; With default window sizes, we can usually display no more than ; 30 frames. ; ; (5) Animation frames can also be written to GIF or MPEG ; files. Defaults can be set in "gamap.defaults", or ; specified via the command line. You can also save individual ; GAMAP plots as GIF files. If you want to animate them later ; (e.g. using ULead's GIF-Animator), make sure to specify the ; RANGE keyword to get identical color schemes (or use contours). ; ; (6) The GAMAP authors wish to point out that it is still relatively ; expensive to produce color plots on the printer. We encourage ; you to try out contour plots and make a test print on a black ; and white printer before you make a color print. ; ; (7) See the documentation in "gamap.pdf" for more information. ; ; EXAMPLES: ; (1) ; GAMAP ; ; operates fully interactively ; ; MULTIPANEL,nplots=6 ; turn on multi-panel environment ; GAMAP ; ; same as above, but produce multi-panel plots with ; ; 6 plots per page ; ; (2) ; GAMAP, 'IJ-AVG-$', tra=4 ; ; Will create a CO (tracer=4) plot for the ND45 diagnostic. ; ; GAMAP will display dialog pickfile box and will scan the ; ; file for all records with ND45 and tracer 4. Those will be ; ; displayed and the user can then select a record to be plotted. ; ; (3) ; GAMAP, [ 45, 28 ], tra=[2,4], date=[940601, 940801], $ ; FileName='~bmy/terra/CTM4/ctm.pch',/ps ; ; In this example the file is fully specified, hence no file ; ; selection dialog will be displayed. GAMAP scans the file ; ; for all records of 'IJ-AVG-$' and 'BIOBSRCE' and tracers ; ; 2 (OX) and 4 (NOX) and it seelcts only those records that ; ; begin on the first of JUNE or AUGUST 1994. Because the ps ; ; flag is set, the output will be directed to the postscript ; ; file 'idl.ps' without first being displayed on the screen. ; ; ; MODIFICATION HISTORY: ; mgs, 12 Nov 1998: VERSION 1.00 ; bmy, 16 Nov 1998: - added defaults for LAT, LEV, LON, PTYPE ; - now prompts for PS ; - now prompts user for /PS output ; bmy, 17 Nov 1998: - now call DEFAULT_RANGE to ensure that ; that LAT, LON, LEVEL have two elements, ; even if there is only one unique value. ; - now uses N_UNIQ.PRO to test for the number ; of unique elements in LON, LAT, and LEVEL. ; mgs, 17 Nov 1998: - finishing touches for first release. ; - added NOFILE keyword ; - added plot type b/w contours ; mgs, 18 Nov 1998: - added timestamp as default when closing ; postscript files ; bmy, 08 Jan 1999: - Will also prompt for totaling (if ; averaging is not selected) ; bmy, 13 Jan 1999: - now prompt user for OUTFILENAME ; bmy, 15 Jan 1999: VERSION 1.02 ; - add support for 3-D data slices ; - clean up user interface so that the user ; menu of plotting options is only invoked ; when plotting a 2-D map. ; bmy, 19 Jan 1999: - added binary flag masking ; - added defaults for averaging and selection ; - improved echoback of information to user ; bmy, 20 Jan 1999: - prompts user again if data block selection ; or averaging selection is out of range ; - fixed bug: now default data block ; selection is saved. ; - Reset PS to 0 and OUTFILENAME to '' if we ; are plotting a 0-D or 3-D data block ; - updated comments ; mgs, 21 Jan 1999: - dimensional information now in subroutine ; - improved binary masking ; - added several Quit options ; - Postscript options now controlled from ; gamap.defaults ; - removed NoTimeStamp keyword (now set in ; gamap.defaults) ; bmy, 12 Feb 1999: VERSION 1.03 ; - now works with data blocks that are ; sub-regions of the globe ; - added functions GAMAP_GetDataBlockRanges ; GAMAP_SelectDataBlock, and ; GAMAP_QueryAverageOrTotal ; - updated comments ; bmy, 17 Feb 1999: VERSION 1.20 ; - Replace DATAINFO.OFFSET by DATAINFO.FIRST, ; which contains the I, J, L indices of ; the first grid box ; - Animation facility added ; - added functions GAMAP_GetModelInfo, ; GAMAP_CheckDataBlockConsistency, ; GAMAP_SelectPlotType, and ; GAMAP_QueryPostScriptOptions. ; - Also renamed function GAMAP_SelectDataBlock to ; GAMAP_SelectDataBlocks, since one can now ; select multiple data blocks ; bmy, 18 Feb 1999: - added /RESETDEFAULTS keyword ; - removed /ANIMATE keyword ; bmy, 19 Feb 1999: - added NOAUTOYRANGE keyword ; - added function GAMAP_AutoYRange ; - added GIFFILENAME keyword ; - added GIF_SAV to common block SAVEVALUES ; - call REPLACE_TOKEN to replace token text ; in DEFAULTGIFFILENAME ; bmy, 22 Feb 1999: - added more animation options ; - added DO_GIF, DO_MPEG, DO_ANIMATE, and ; MPEGFILENAME keywords ; - added GAMAP_QueryAnimationOptions routine ; bmy, 23 Feb 1999: - small bug fixes ; bmy, 04 Mar 1999: - added internal routines GAMAP_FindNearestEdges ; and GAMAP_GetDefaultRanges ; - now force lat/lon ranges to coincide with ; grid box edges ; - warn user if lat range contains half-polar ; boxes, since TVIMAGE will treat them as ; whole boxes and the map overlay will be ; inaccurate! ; bmy, 05 Mar 1999: - Clean up FILEINFO/DATAINFO matching process ; - renamed/reorganized internal functions\ ; bmy, 20 Mar 1999: - bug fix for default ranges (may need more ; fixing later on) ; mgs, 22 Mar 1999: - added ALREADY_PS flag for multi-panel use ; - animation now only if !p.multi does not ; have more than 1 panel to display ; mgs, 23 Mar 1999: - improved comments and examples ; - removed unnecessary function MatchFileInfo... ; (easier with make_selection) ; - changed all "string booleans" to booleans ; - Do_Animation now an entirely local variable ; mgs, 25 Mar 1999: - few minor bug fixes ; - improved handling of default ranges ; - detect out of range in record selection ; ;- ; Copyright (C) 1998, 1999, Martin Schultz and Bob Yantosca, ; Harvard University ; This software is provided as is without any warranty ; whatsoever. It may be freely used, copied or distributed ; for non-commercial purposes. This copyright notice must be ; kept with any copy of this software. If this software shall ; be used commercially or sold as part of a larger package, ; please contact the author to arrange payment. ; Bugs and comments should be directed to mgs@io.harvard.edu ; or bmy@io.harvard.edu with subject "IDL routine gamap" ;------------------------------------------------------------- function GAMAP_CheckDataBlockConsistency, D, FileInfo, FileIndex ;================================================================= ; GAMAP_CheckDataBlockConsistency checks data blocks selected ; for animation to make sure they meet the following criteria: ; ; (1) All data blocks have the same dimensions ; Since old ASCII punch files do not contain dimensional ; information we need to supplement it by using the grid ; information ; (2) All data blocks start with the same grid box ; (3) All data blocks were produced by the ame model with ; identical resolution. ;================================================================= ;----------------------------------------------------------------- ; Index array of data blocks that match criteria ;----------------------------------------------------------------- GoodDataBlocks = IntArr( N_Elements( D ) ) ;----------------------------------------------------------------- ; Get the MODELINFO structure corresponding to D[0] ;----------------------------------------------------------------- RefModelInfo = FileInfo[ FileIndex[0] ].ModelInfo ;----------------------------------------------------------------- ; Get the max # of elements in D[0].DIM, D[0].FIRST, ; REFMODELINFO.NAME, and REFMODELINFO.RESOLUTION. ; The equivalence tests will produce a match for each ; element of these fields. ;----------------------------------------------------------------- DimSize = N_Elements( D[0].Dim ) FirstSize = N_Elements( D[0].First ) NameSize = N_Elements( RefModelInfo.Name ) ResSize = N_Elements( RefModelInfo.Resolution ) ; ### debug ; print,'### DIMSIZE, FIRSTSIZE, NAMESIZE, RESSIZE' ; print,DIMSIZE, FIRSTSIZE, NAMESIZE, RESSIZE ; print,'D[0].dim = ',D[0].Dim ;----------------------------------------------------------------- ; Loop over the records of D, skipping D[0] ;----------------------------------------------------------------- for N = 0, N_Elements( D ) - 1 do begin ; Get the MODELINFO structure for this data block ThisModelInfo = FileInfo[ FileIndex[N] ].ModelInfo ThisGridInfo = * (FileInfo[ FileIndex[N] ].GridInfo ) if ( FileInfo[ FileIndex[N] ].FileType eq 0 ) then begin if ( Total( D[N].Dim ) le 0 ) then begin D[N].Dim[0] = ThisGridInfo.IMX D[N].Dim[1] = ThisGridInfo.JMX D[N].Dim[2] = ThisModelInfo.NTROP print,'### CheckConsistency: dimensions set from grid info !' endif if ( Total( D[N].First ) eq 0 ) then D[N].First = [1, 1, 1] endif ; Check dimensions of this data block if (D[N].Dim[2] lt 0) then D[N].Dim[2] = max(D.Dim[2]) Ind0 = Fix( Total( D[N].Dim eq D[0].Dim ) ) ; Check the indices of the first grid box of this data block Ind1 = Fix( Total( D[N].First eq D[0].First ) ) ; Check model name against that of first data block Ind2 = Fix( ThisModelInfo.Name eq RefModelInfo.Name ) ; Check model resolution against that of first data block Ind3 = Fix( Total( ThisModelInfo.Resolution eq $ RefModelInfo.Resolution ) ) ; Require dimension, offset, model name, model res to match if ( Ind0 eq DimSize AND Ind1 eq FirstSize AND $ Ind2 eq NameSize AND Ind3 eq ResSize ) then GoodDataBlocks[N] = 1 ; ### debug ; print,'### ind0,ind1,ind2,ind3' ; print,ind0,ind1,ind2,ind3 ; print,'D[N].dim = ',D[N].Dim endfor ;-------------------------------------------------------------------- ; Return data blocks that match the criteria ;-------------------------------------------------------------------- Ind0 = Where( GoodDataBlocks eq 1 ) D = D[ Ind0 ] ;-------------------------------------------------------------------- ; Print warning messages for those data blocks that ; did not meet the criteria and were eliminated. ; Will get more sophisticated later... ;-------------------------------------------------------------------- Ind1 = Where( GoodDataBlocks eq 0 ) if ( Ind1[0] ge 0 ) then begin for N = 0, N_Elements( Ind1 ) - 1 do begin S = StrTrim( String( N, Format='(i10)' ), 2 ) S = 'Data Block ' + S + ' has been eliminated from animation...' Message, S, /Info, /NoName endfor endif ; Successful Return! return, 1 end ;---------------------------------------------------------------------------- function GAMAP_FindNearestEdges, $ F, Range, Latitude=Latitude, Longitude=Longitude ;================================================================= ; Function GAMAP_FindNearestEdges (bmy, 3/4/99) converts an ; input range for lat or lon to the nearest grid box edge. ; This function is called for one individual data block. ; LATITUDE and LONGITUDE are boolean keywords that determine ; whether to return latitude or longitude edges. ; ; mgs, 23 Mar 1999: - de-hardwired half polar check and streamlined ; operation. ; bmy, 25 Mar 1999: - bug fix (gt instead of ge) ;================================================================= ; Get structures from FILEINFO ModelInfo = F.ModelInfo GridInfo = *( F.GridInfo ) ; Grid box edge arrays for latitude or longitude Longitude = Keyword_Set( Longitude ) Latitude = Keyword_Set( Latitude ) if ( Longitude ) then Edges = GridInfo.XEDGE if ( Latitude ) then Edges = GridInfo.YEDGE NEdges = n_elements(Edges) ; Half-Polar box warning if ( Latitude AND ModelInfo.HalfPolar ) then begin if ( Range[0] lt Edges[1] OR Range[1] gt Edges[NEdges-2] ) then begin S = StrTrim( Replicate( Byte( '=' ), 70 ), 2 ) print Message, S, /Info, /NoName Message, 'WARNING!!! TVIMAGE will display Half-Polar Boxes' + $ ' as full boxes', /Info, /NoName Message, 'You might want to exclude the poles from your lat range!',$ /Info, /NoName Message, S, /Info, /NoName print endif endif ; return nearest left/bottom edge Ind = where(Edges le Range[0]) if (Ind[0] ge 0) then $ LowerEdge = Edges[Ind[n_elements(Ind)-1]] $ else $ LowerEdge = Edges[0] ; can't go beyond first edge Ind = where(Edges ge Range[1]) if (Ind[0] ge 0) then $ UpperEdge = Edges[Ind[0]] $ else $ UpperEdge = Edges[NEdges-1] ; can't go beyond last edge ; *** ??? do we need to take extra care of longitude ranges ; *** ??? if we span the dateline ??? ; *** commented out old version (mgs, 23 Mar 1999) ;* ; Work array...add Range to the Edges array ;* Arr = [ Edges, Range ] ;* ; Sort array and find unique values ;* Arr = Arr( Uniq( Arr, Sort( Arr ) ) ) ;* ; If RANGE[0] is an lower edge, then return it ;* ; Otherwise, locate the edge closest to RANGE[0] ;* Ind = Where( Edges eq Range[0] ) ;* if ( Ind[0] ge 0 ) then begin ;* LowerEdge = Edges[ Ind ] ;* endif else begin ;* Ind = Where( Arr eq Range[0] ) ;* LowerEdge = Arr[ Ind-1 ] ;* endelse ;* ; If RANGE[1] is an upper edge, then return it ;* ; Otherwise, locate the edge closest to RANGE[1] ;* Ind = Where( Edges eq Range[1] ) ;* if ( Ind[0] ge 0 ) then begin ;* UpperEdge = Edges[ Ind ] ;* endif else begin ;* Ind = Where( Arr eq Range[1] ) ;* UpperEdge = Arr[ Ind+1 ] ;* endelse ; Return the nearest grid box edges return, [ LowerEdge, UpperEdge ] end ;---------------------------------------------------------------------------- function GAMAP_GetDataBlockRanges, D, F, Lon_Range, Lat_Range, Lev_Range ;================================================================= ; GAMAP_GETDATABLOCKRANGES computes the extent of the data block ; in longitude, latitude, and altitude (bmy, 3/5/99 ; ; NOTES: ; (1) D and F are the DATAINFO and FILEINFO structures that ; correspond to theh first data block. Thus we pass ; D[0] and FILEINFO[ FILEINDEX[0] ] from GAMAP.PRO. ; ; (2) **** This function should be moved to ctm_retrieve_data at ; **** some point, storing the coordinate information in the ; **** datainfo structure as soon as we read the file header ; (mgs, 03/23/99) ; ; (3) I don't like clipping to -88 .. 88 here! If the data does ; contain the polar latitudes, the associated range info should ; reflect this!! ;================================================================= ; Common blocks @gamap_cmn ;----------------------------------------------------------------- ; Get the MODELINFO and GRIDINFO fields that correspond to D[0] ;----------------------------------------------------------------- ModelInfo = F.ModelInfo GridInfo = *( F.GridInfo ) ;----------------------------------------------------------------- ; For ASCII punch files (FILEINFO.FILETYPE = 0), the dimensions ; don't become visible until after each data block is read in. ; If D[0].DIM or D[0].FIRST are all zero, then get the dimensions ; from the GRIDINFO structure. ;----------------------------------------------------------------- if ( F.FileType eq 0 ) then begin if ( Total( D.Dim ) le 0 ) then begin D.Dim[0] = GridInfo.IMX D.Dim[1] = GridInfo.JMX D.Dim[2] = 1 ;ModelInfo.NTROP print,'### GetDataBlockRanges: dimensions set from grid info!' endif if ( Total( D.First ) eq 0 ) then D.First = [1, 1, 1] endif ;----------------------------------------------------------------- ; Longitude extent of the current data block ; Take care of longitudes that wrap around the date line ; ; NOTE: Convert from FORTRAN notation to IDL notation by ; subtracting 1 from D[0].FIRST and D[0].DIM ;----------------------------------------------------------------- I0 = D.First[0] - 1 I1 = ( D.First[0] + D.Dim[0] - 2 ) mod GridInfo.IMX Lon_Range = [ GridInfo.XEdge[I0], GridInfo.XEdge[I1+1] ] ;----------------------------------------------------------------- ; Latitude extent of the current data block ; Truncate lat range at -88 and +88, since TVIMAGE can't handle ; half polar boxes -- it can't print half a pixel!! (bmy, 3/19/99) ;----------------------------------------------------------------- I0 = D.First[1] - 1 I1 = D.First[1] + D.Dim[1] - 2 Lat_Range = [ GridInfo.YEdge[I0], GridInfo.YEdge[I1+1] ] If ( Lat_Range[0] lt -89 ) then Lat_Range[0] = GridInfo.YEdge[I0+1] If ( Lat_Range[1] gt 89 ) then Lat_Range[1] = GridInfo.YEdge[I1 ] ;----------------------------------------------------------------- ; Vertical extent of the current data block ; Set the range correctly for data blocks that only have 1 level ;----------------------------------------------------------------- if (ChkStru(GridInfo,'LMX')) then begin ; stop ; ######### LevInd = Indgen( GridInfo.LMX ) + 1 I0 = D.First[2] - 1 I1 = D.First[2] + D.Dim[2] - 2 if ( I1 ge 0 ) then begin Lev_Range = [ LevInd[I0], LevInd[I1] ] endif else begin Lev_Range = [ LevInd[I0], LevInd[I0] ] endelse endif else begin ; 2D gridinfo structure LevInd = 1 Lev_Range = [ 1, 1 ] endelse ;----------------------------------------------------------------- ; Successful return! ;----------------------------------------------------------------- return, 1 end ;---------------------------------------------------------------------------- function GAMAP_AutoYRange, $ D, AvMode, TotalMode, Lon, Lat, Level ;================================================================= ; Function GAMAP_AutoYRange computes the absolute ; minimum and maximum values for all of the animation frames. ; This is needed so that the colorbars of the individual ; animation frames will all have the same range. (bmy, 2/19/99) ; This routine is also called for multi-panel plots if keyword ; /AUTORANGE is set. ;================================================================= ; Common blocks (DEBUG flag) @gamap_cmn ; Local arrays to hold min and max values. MinData = FltArr( N_Elements( D ) ) MaxData = FltArr( N_Elements( D ) ) ; Loop over all selected data blocks... for N = 0, N_Elements( D ) - 1 do begin Success = CTM_Get_DataBlock( Data, D[N].Category, $ ILun=D[N].Ilun, $ Tracer=D[N].Tracer, $ Tau0=D[N].Tau0, $ Average=AvMode, Total=TotalMode, $ Lev=Level, Lat=Lat, $ Lon=Lon, /NoPrint ) ; For each data block found, compute its max and min ; values, and store them in the MAXDATA and MINDATA arrays. if ( Success ) then begin TmpMin = Min( Data, Max=TmpMax ) MinData[ N ] = TmpMin MaxData[ N ] = TmpMax ; Set status of that record to read D[N].status = 1 ;**** Debug output (bmy, 2/19/99) if ( DEBUG ) then begin print, '### GAMAP.PRO: N, MINDATA, MAXDATA: ', N, TmpMin, TmpMax endif endif endfor ; AUTOYRANGE contains the "global" min and max values AutoYRange = [ Min( MinData ), Max( MaxData ) ] ; Undefine the data array so that we don't take up resources Undefine, Data ; Return the "global" Min and Max values return, AutoYRange end ;---------------------------------------------------------------------------- pro GAMAP_PrintDimInfo,Arr,Ind,NDim=NDim,Message=MStr ;================================================================= ; GAMAP_PRINTDIMINFO prints the current number of dimensions and ; their labels. ; ; Arr contains flags as Arr = [ NLon, NLat, NLevel ] ; Nxxx = 0 : dimension was averaged/totaled ; = 1 : dimension contains single value ; = 2 : is a multi element vector ; ; Ind is a selection array which dimensions to print by name ; (default: all dimensions with Arr=2) ; ; NDim returns the number of "active" dimensions ;================================================================= DimStr = [ 'Longitude', 'Latitude', 'Altitude' ] DimInd = Where( Arr gt 1 ) if (n_elements(Ind) eq 0) then Ind = DimInd NDim = n_elements(DimInd) * ( DimInd[0] ge 0 ) if (n_elements(Mstr) eq 0) then MStr = '' $ else MStr = '. ' + MStr MesgStr = 'Selected data is ' + $ StrTrim( String( NDim, format='(i4)'), 2 ) + '-D' if ( NDim gt 0 ) then begin MesgStr = MesgStr + MStr + ' [ ' + $ string(DimStr[Ind],format='(3(A,:,","))' ) + ' ].' endif else begin MesgStr = MesgStr + '.' endelse if (not !QUIET) then print Message, MesgStr, /Info, /NoName return end ;---------------------------------------------------------------------------- pro GAMAP_QueryAnimationOptions, NDim, $ Do_GIF, GIFFileName, Do_MPEG, $ MPEGFileName, Quit=Quit ;================================================================= ; Subroutine GAMAP_QueryAnimationOptions queries the user for ; animation options. (bmy, 2/22/99) ; ; mgs, 23 Mar 1999: - now uses yesno function and is streamlined ; - doesn't taper with animation status any more ; - Do_XX are always returned as boolean values ; although they may have been passed in as strings ;================================================================= FORWARD_FUNCTION yesno ; Common blocks common SaveValues, Sel_Sav, Lev_Sav, Lon_Sav, $ Lat_Sav, PType_Sav, Avg_Sav Quit = 0 ;----------------------------------------------------------------- ; if NDim = 0 then no GIF or MPEG feasible ; If DO_GIF = *QUERY, query user whether or not to save ; animation frames to GIF output ;----------------------------------------------------------------- if (NDim eq 0) then begin Do_GIF = 0 Do_MPEG = 0 return endif if ( string(Do_GIF) eq '*QUERY' ) then begin print Do_GIF = yesno('Save frames to GIF files?',default=0,/QUIT) if (Do_GIF lt 0) then begin Quit = 1 return endif ;* LStr = '' ;* LStrDefault = GIF_Sav ;* Read, LStr, $ ;* Prompt='Save frames to GIF files? (Y/N) (default : ' + $ ;* LStrDefault + ', Q=Quit) >> ' ;* ;* if ( LStr eq '' ) then LStr = LStrDefault ;* LStr = StrMid( StrUpCase( StrTrim( LStr, 2 ) ), 0, 1 ) ;* ;* ; QUIT if LSTR = Q ;* if ( LStr eq 'Q' ) then begin ;* Quit = 1 ;* return ;* endif ;* ;* ; Force LSTR to be Y or N ;* if ( LStr ne 'Y' ) then LStr = 'N' ;* if ( LStr eq 'Y' ) then Do_GIF = '1' ;* endif else Do_GIF = fix(Do_GIF) ;----------------------------------------------------------------- ; Query user for the GIF file name. The GIF file name may ; contain a replaceable token for the frame number. ;----------------------------------------------------------------- if ( Do_GIF eq 1 ) then begin if ( StrUpCase( GIFFileName ) eq '*QUERY' ) then begin LStr = '' LStrDefault = 'frame%NNN%.gif' Read, LStr, Prompt='GIF File Name (default : ' + $ LStrDefault + ', Q=Quit) >> ' LStr = StrTrim( LStr, 2 ) ; Quit if 'Q' if ( StrUpCase( LStr ) eq 'Q' ) then begin Quit = 1 return endif if ( LStr eq '' ) then LStr = LStrDefault GIFFileName = LStr endif endif ;----------------------------------------------------------------- ; If 0-D or 3-D , no MPEG possible ; If DO_MPEG = *QUERY, query user whether or not to save ; animation frames to MPEG output ;----------------------------------------------------------------- if ( NDim eq 3 ) then begin ; 0-D handled above Do_MPEG = 0 return endif if ( string(Do_MPEG) eq '*QUERY' ) then begin Print Do_MPEG = yesno('Save frames to MPEG file?',default=0,/QUIT) if (Do_MPEG lt 0) then begin Quit = 1 return endif endif else Do_MPEG = fix(Do_MPEG) ;----------------------------------------------------------------- ; Query user for the MPEG file name. ;----------------------------------------------------------------- if ( Do_MPEG eq 1 ) then begin if ( StrUpCase( MPEGFileName ) eq '*QUERY') then begin LStr = '' LStrDefault = 'result.mpeg' Read, LStr, Prompt='MPEG File Name (default : ' + $ LStrDefault + ', Q=Quit) >> ' LStr = StrTrim( LStr, 2 ) ; Quit if 'Q' if ( StrUpCase( LStr ) eq 'Q' ) then begin Quit = 1 return endif if ( LStr eq '' ) then LStr = LStrDefault MPEGFileName = LStr endif endif return end ;---------------------------------------------------------------------------- pro GAMAP_QueryAverageOrTotal, Arr, NDim, AvMode, TotalMode, Quit=Quit ;==================================================================== ; Procedure GAMAP_QueryAverageOrTotal queries the user as to ; averaging or totaling options (bmy, 2/12/99) ;==================================================================== ; GAMAP common blocks @gamap_cmn common SaveValues, Sel_Sav, Lev_Sav, Lon_Sav, $ Lat_Sav, PType_Sav, Avg_Sav ; Initialize variables AvMode = 0 TotalMode = 0 ;-------------------------------------------------------------------- ; If at least one dimension contains more than one value, ; query user for averaging (or totaling) ;-------------------------------------------------------------------- if ( NDim ge 1 ) then begin LStr = '' LStrDefault = StrTrim( String( Avg_Sav, Format='(i4)' ), 2 ) if (NDim eq 3) then Str0 = 'Visualize 3D' $ else Str0 = 'No averaging' Print Repeat_AvMode_Selection: print,'Do you want to average or total the data?' Read, LStr, $ prompt='(0='+Str0+', 1=lon, 2=lat, 4=alt, ' + $ '8=Total, Q=Quit, Default=' + LStrDefault + ') >> ' LStr = StrTrim(Lstr,2) if ( StrUpcase(StrMid(LStr,0,1)) eq 'Q' ) then begin Quit = 1 return endif if ( LStr eq '' ) then LStr = LStrDefault ;----------------------------------------------------------------- ; Ask user again if AVMODE is out of range ;----------------------------------------------------------------- AvMode = Fix( LStr ) if ( AvMode lt 0 OR AvMode gt 15 ) then begin Message, 'Selection must be >= 0 and <= 15. Choose again...', $ /Info, /NoName goto, Repeat_AvMode_Selection endif ; Save for next time... Avg_Sav = AvMode ;----------------------------------------------------------------- ; Binary masking. Make sure that the dimension(s) selected ; for averaging or totaling each have more than one point. ; Reduce NDim by one for each dimension that is to be ; averaged or totaled. ; ; If a dimension contains a single point, then we cannot ; average or total along that dimension. Use a Hexadecimal bit ; mask for a bitwise AND operation to clear AVMODE and ; TOTALMODE. This will assure that CTM_PLOT will select the ; proper plot labels. ; ; A refresher course: ; Decimal: 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 ; Hexadecimal: 0 1 2 3 4 5 6 7 8 9 A B C D E F ; ; See the IDL manual for info about the Bitwise AND operation. ;----------------------------------------------------------------- ; Clear averaging flag for all dimensions with only one value mask = 'F8'x + 4*(Arr[2] gt 1) + 2*(Arr[1] gt 1) + (Arr[0] gt 1) ;**** Debug output (bmy, 2/19/99) if ( DEBUG ) then print,'### AVERAGING MASK = ',mask AvMode = (AvMode AND mask) ;----------------------------------------------------------------- ; Set dimensions to be averaged to zero and reduce NDim ; Ind contains 1 for each dimenson *not* to be averaged ;----------------------------------------------------------------- Ind = ( (AvMode AND [1,2,4]) eq 0 ) Arr = Arr * Ind ; sets averaged dimensions to zero ; Settings for total TotalMode = 0 if ( (AvMode AND 8) gt 0 ) then begin TotalMode = ( AvMode AND 'F7'x ) AvMode = 0 endif ;----------------------------------------------------------------- ; Echo back information about which dimensions will be ; averaged or totaled by GAMAP (those for which ARR[i] = 0) ; NDIM is adjusted in GAMAP_PrintDimInfo ;----------------------------------------------------------------- Ind = Where( Arr eq 0, Count ) if (Count gt 0) then begin if (AvMode gt 0) then $ Gamap_PrintDimInfo,Arr,Ind,NDim=NDim, $ Message='Data averaged over' $ else if (TotalMode gt 0 ) then $ Gamap_PrintDimInfo,Arr,Ind,NDim=NDim, $ Message='Data totaled over' endif endif end ;---------------------------------------------------------------------------- pro GAMAP_QueryPostScriptOptions, $ Do_PS, TimeStamp, OutFileName ;================================================================= ; Subroutine GAMAP_QueryPostScriptOptions queries the user for ; PostScript options (bmy, 2/16/99) ; ; mgs, 23 Mar 1999: - removed all animation related parameters. ; This routine is only called if no animation is active. ; - removed "double" keywords (Bob's bad day ;-) ; - removed PScrLabel, since it was equivalent ; to the default in close_device ; - Do_PS will always return numeric value. -1 ; if user shall be queried after plot ;================================================================= FORWARD_FUNCTION yesno ;----------------------------------------------------------------- ; Set or query postscript options if Do_PS ;----------------------------------------------------------------- if ( Do_PS eq string(1) ) then begin ;-------------------------------------------------------------- ; Query user for PostScript file name ;-------------------------------------------------------------- if (strupcase(OutfileName) eq '*QUERY') then begin LStr = '' LStrDefault = 'idl.ps' Read, LStr, Prompt='Postscript File Name (default : ' + $ LStrDefault + ') >> ' LStr = StrTrim(LStr,2) if ( LStr eq '' ) then LStr = LStrDefault OutFileName = LStr endif ;-------------------------------------------------------------- ; Query User to add a time stamp to the PostScript plot ;-------------------------------------------------------------- if (strupcase(TimeStamp) eq '*QUERY') then begin TimeStamp = yesno('Add time stamp to plot?',default=1) endif else $ TimeStamp = fix(TimeStamp) ; either '0' or '1' endif if ( StrUpCase(Do_PS) eq '*QUERY' ) then $ Do_PS = -1 $ else $ Do_PS = fix(Do_PS) return end ;---------------------------------------------------------------------------- function GAMAP_InterpreteSel, S ;================================================================= ; Interprete a selection string of the form N1,N2-N3,... ;================================================================= ;----------------------------------------------------------------- ; First try to seperate Sel by blanks. This form can be used ; to list individual records. Then seperate each term by ',' ; which is an alternative method. Each comma (or blank) seperated ; term can then contain a '-' to indicate a range. ;----------------------------------------------------------------- res = -1L S = StrTrim(S,2) on_ioerror,bad_format csel = str_sep(S,',') for j=0,n_elements(csel)-1 do begin dsel = strtrim(str_sep(csel[j],'-'), 2) for i=0, n_elements(dsel)-1 do $ if (strpos(dsel[i], ' ') gt 0) then $ goto, Bad_Format if (n_elements(dsel) eq 1) then $ dsel = [ dsel, dsel ] rec0=-1 rec1=-1 reads,dsel,rec0,rec1 res = [ res, indgen(rec1-rec0+1)+rec0 ] endfor on_ioerror,NULL ; remove dummy if (n_elements(res) gt 1) then res = res[1:*] ; everything seems OK: return res-1 as sel for proper indexing return,res-1 Bad_Format: message,'Bad numerical format in selection string!',/Continue,/NoName message,'Allowed syntax is N[-M][,N2[-M2]][,...] (example: 1,3,6-10)', $ /NoName,/INFO return,res-1 end ;---------------------------------------------------------------------------- function GAMAP_SelectDataBlocks, Max_Sel, Quit=Quit ;================================================================= ; GAMAP_SelectDataBlocks prompts the user to select a a data ; block (for simple plotting) or range of data blocks (for ; (timeseries animation). (bmy, 2/17/99) ;================================================================= ; GAMAP Common blocks common SaveValues, Sel_Sav, Lev_Sav, Lon_Sav, $ Lat_Sav, PType_Sav, Avg_Sav ; Initialize variables Sel = -1 Quit = 0 Repeat_Selection: print ; set default string: make sure it does not exceed allowed range LStrDefault = StrTrim(Sel_Sav,2) saved = GAMAP_InterpreteSel(LStrDefault) if (max(saved) gt Max_Sel) then LStrDefault = '1' LStr = '' Read, LStr, Prompt='Enter data block or data block range (default : ' + $ LStrDefault+', Q=Quit) >> ' LStr = strtrim(LStr,2) if ( LStr eq '' ) then LStr = LStrDefault if ( StrUpcase( LStr ) eq 'Q' ) then begin Quit = 1 return,-1 endif sel = GAMAP_InterpreteSel(LStr) ; already complained about syntax errors if (min(Sel) lt -1) then goto,Repeat_Selection ; complain about selections out of range test = where(sel lt 0 OR sel ge Max_Sel) if (test[0] ge 0) then begin message,'Record number(s) out of range! Try again!', $ /Continue,/NoName goto,Repeat_Selection endif ; everything seems OK: return res-1 as sel for proper indexing print,'Selected records : ',sel+1,format='(A,512I4)' ; Save default for next time and return Sel_Sav = LStr return, Sel end ;---------------------------------------------------------------------------- pro GAMAP_SelectPlotType, $ NDim, NLon, NLat, NLevel, AvMode, TotalMode, PType, Quit=Quit ;================================================================= ; GAMAP_SelectPlotType (bmy, 2/17/99) does the following: ; ; (1) For a 2-D plot (NDIM=2), query the user for plot options ; (2) For 0-D, 1-D, and 3-D plots, print informational message ; ; NOTES: PTYPE is returned as -1 unless otherwise assigned here. ; ; If QUIT=1 is returned to the calling program, then this ; denotes that the user aborted the selection process, ; and wishes to quit GAMAP. ;================================================================= ; Common blocks common SaveValues, Sel_Sav, Lev_Sav, Lon_Sav, $ Lat_Sav, PType_Sav, Avg_Sav ; Initialize variables Quit = 0 ;-------------------------------------------------------------------- ; Case Statement for selecting plot type based on # of dimensions ;-------------------------------------------------------------------- case ( NDim ) of ;----------------------------------------------------------------- ; 0 Dimensions: Print value of data ;----------------------------------------------------------------- 0: begin if ( NLon eq 1 AND NLat eq 1 AND NLevel eq 1 ) then begin DimStr = '' endif else begin if ( AvMode gt 0 ) $ then DimStr = 'average ' $ else DimStr = 'total ' endelse DimStr = 'GAMAP will print the ' + DimStr + 'value of the data.' Message, DimStr, /Info, /NoName PType = -1 end ;----------------------------------------------------------------- ; 1 Dimension: Create a line plot ;----------------------------------------------------------------- 1: begin Message, 'GAMAP will create a line plot.', /Info, /NoName Ptype = -1 end ;----------------------------------------------------------------- ; 2 Dimensions: query user for contour plot or pixel plot options ;----------------------------------------------------------------- 2: begin print PType = Choice( ['B/W Contour lines', $ 'Colored contour lines', $ 'Filled contours', $ 'Smooth Pixel Plot', $ 'Coarse Pixel Plot'], $ Title='Select 2-D plot type:', $ Default=Ptype_Sav, $ /BY_INDEX) ;if ( PType lt 0 ) then return ; aborted if ( PType lt 0 ) then begin Quit = 1 return endif ; Store for next iteration Ptype_Sav = Ptype ;**** now done in main program (bmy, 2/16/99) ; set contour colors to all black if ptype eq 0: ;if (PType eq 0) then c_colors = 1 end ;----------------------------------------------------------------- ; 3 Dimensions: Visualize data cube ;----------------------------------------------------------------- 3: begin Message, 'GAMAP will visualize the data cube.', /Info, /NoName PType = -1 end endcase return end ;---------------------------------------------------------------------------- pro GAMAP_StoreGridInfo, FileInfo, FileIndex, Quit=Quit ;==================================================================== ; Procedure GAMAP_StoreGridInfo (bmy, 3/5/99) computes the ; GRIDINFO structure those elements of FILEINFO that correspond to ; elements of DATAINFO (contained in index array FILEIND). ;==================================================================== FORWARD_FUNCTION ChkStru, CTM_Grid ; First find the unique elements of FILEINDEX. ; This prevents re-storing the GRIDINFO structure ; for the same FILEINFO multiple times. F = FileIndex( Uniq( FileIndex, Sort( FileIndex ) ) ) ; Loop over the unique elements of FILEINDEX for N = 0, N_Elements( F ) - 1 do begin ; Get the MODELINFO structure from FILEINFO ModelInfo = FileInfo[ F[ N ] ].ModelInfo ; Check to make sure MODELINFO is a valid structure if ( not ChkStru( ModelInfo, [ 'NAME', 'RESOLUTION' ] ) ) then begin Message, 'Invalid MODELINFO structure!!', /Continue Quit = 1 return endif ; Use MODELINFO to construct the GRIDINFO structure ; if none has been previously defined (**** mgs, 03/30/99) if (not Ptr_Valid( FileInfo[ F[ N ] ].GridInfo ) ) then begin GridInfo = CTM_Grid( ModelInfo ) ; Check to make sure GRIDINFO is a valid structure if ( not ChkStru( GridInfo, [ 'IMX', 'JMX' ] ) ) then begin Message, 'Invalid GRIDINFO structure!!', /Continue Quit = 1 return endif ; Make a new pointer to the GRIDINFO structure in FILEINFO FileInfo[ F[ N ] ].GridInfo = Ptr_New( GridInfo ) endif endfor return end ;---------------------------------------------------------------------------- function GAMAP_UserRangeEntry, $ Name, Default, Quit=Quit, Format=FFormat, _EXTRA=e ;================================================================= ; GAMAP_USERRANGEENTRY prompts the user to enter a value or ; range for a dimension of CTM data (e.g. longitude, latitude, ; levels). If the user hits return, then default values are used. ; ; Now add FORMAT keyword (bmy, 2/19/99) ;================================================================= ; Pass External Functions FORWARD_FUNCTION Default_Range ; Error checking if ( N_Elements( FFormat ) eq 0 ) then FFormat = '(f10.2)' ; Initialize variables Quit = 0 ; Construct default string ; **** Change default string to F10.2 (bmy, 2/19/99) LStrDefault = StrTrim( String( default[0], Format=FFormat ), 2 ) + $ '..' + StrTrim( String( default[1], Format=FFormat ), 2 ) ; Query user LStr = '' print Read, LStr, Prompt='Enter '+name+' or '+name+' range (default : ' + $ LStrDefault+', Q=Quit) >> ' if ( LStr eq '' ) then LStr = LStrDefault if ( StrUpcase(LStr) eq 'Q' ) then begin Quit = 1 return,-1 endif ; Convert input to sorted numeric 2 element vector result = Default_Range( LStr, default, /Limit2, _EXTRA=e ) ; Store for next iteration default = [ result[0], result[1] ] return,result end ;---------------------------------------------------------------------------- pro GAMAP, DiagN, $ FileName=FileName, NoFile=NoFile, $ Tracer=Tracer, Tau0=Tau0, Date=Date, $ PS=Do_PS, TimeStamp=TimeStamp, $ OutFileName=OutFileName, $ Do_GIF=Do_GIF, GIFFileName=GIFFileName, $ Frame0=Frame0, $ Do_MPEG=Do_MPEG, MPEGFileName=MPEGFileName, $ ResetDefaults=ResetDefaults, AutoRange=AutoRange, $ YRange=YRange, $ Help=Help, _EXTRA=e ;================================================================= ; GAMAP: Main program ; Primitive user-interface for ctm_plot ; Allows menu driven interactive plotting ;================================================================= ; Pass internal and external functions FORWARD_FUNCTION ChkStru, $ Choice, $ Make_Selection, $ NYMD2Tau, $ Replace_Token ; Common blocks @gamap_cmn common FirstTimeFlag, FirstTime common SaveValues, Sel_Sav, Lev_Sav, Lon_Sav, $ Lat_Sav, PType_Sav, Avg_Sav ; First-time initializations, etc if ( N_Elements( FirstTime ) eq 0 ) then FirstTime = 1 ResetDefaults = Keyword_Set( ResetDefaults ) if (ResetDefaults) then FirstTime = 1 if (FirstTime) then begin print,'=============================================================' print print,' WELCOME TO G A M A P --- Version 1.30, Mar.1999' print print,' Global Atmospheric Modeling (output) Analysis Package' print,' Martin Schultz and Bob Yantosca, Harvard University' print,'=============================================================' print ; Set defaults for user selections ;**** LEV_SAV, LON_SAV, and LAT_SAV are now reset ;**** using GAMAP_GetDefaultRanges (bmy, 3/4/99) Lev_Sav = [ -999 , 999 ] Lon_Sav = [ -999.0, 999.0 ] Lat_Sav = [ -999.0, 999.0 ] PType_Sav = 3 Avg_Sav = 0 Sel_Sav = '1' FirstTime = 0 ResetDefaults = 1 endif Do_Animation = 0 Do_MultiPanel = 0 if ( Keyword_Set( Help ) ) then begin print,'Usage: gamap [,OPTIONS]' print print,' OPTIONS:' print print,' DIAGN : Diagnostic number or category name' print,' TRACER : Tracer number' print,' TAU0, DATE : Start time of the data record' print,' FILENAME : CTM output file name' print,' /NOFILE : Use previously read records only' print,' /PS : Create postscript file without prompt.' print,' OUTFILENAME : Name of postscript output file' print,' /DO_GIF : Save animation frames as GIF files' print,' GIFFILENAME : Name of GIF output file ' print,' /DO_MPEG : Save animation frames into an MPEG movie' print,' MPEGFILENAME : Name of MPEG output file ' print,' /RESETDEFAULTS: Start GAMAP all over but keep data' print,' YRANGE=[Y0,Y1]: Specify scale' print,' /AUTORANGE : Scale all records equally' print,' /HELP : Display this help page.' print print,'For more information consult the GAMAP reference guide' print,'This is available as a pdf file.' return endif ;================================================================= ; set defaults ; Use defaults from gamap.defaults ; Command line options override defaults ;================================================================= if ( n_elements( diagn ) eq 0 ) then diagn = 0 ; all diagnostics ;* if ( n_elements( filename ) eq 0 AND not keyword_set(NoFile) ) then $ ;* filename = '' ; user selects ;*** still not sure whether this is a correct change: ;*** it forces the default filemask if nothing was read before which we don't want if (keyword_set(NoFile) AND ptr_valid(pGlobalFileInfo) ) then $ undefine,FileName $ else $ if (n_elements(FileName) eq 0) then FileName = '' ; postscript options if (n_elements(Do_PS) gt 0) then $ Do_PS = Keyword_Set( Do_PS ) $ else $ Do_PS = CreatePostscript if (n_elements(TimeStamp) gt 0) then $ TimeStamp = Keyword_Set( TimeStamp ) $ else $ TimeStamp = AddTimeStamp if ( N_Elements( OutFileName ) eq 0 ) $ then OutFileName = DefaultPSFilename ; GIF defaults if (n_elements(Do_GIF) gt 0) then $ Do_GIF = Keyword_Set( Do_GIF ) $ else $ Do_GIF = CreateGIF if ( N_Elements( GIFFileName ) eq 0 ) $ then GIFFileName = DefaultGIFFileName ; MPEG defaults if (n_elements(Do_MPEG) gt 0) then $ Do_MPEG = Keyword_Set( Do_MPEG ) $ else $ Do_MPEG = CreateMPEG if ( N_Elements( MPEGFileName ) eq 0 ) $ then MPEGFileName = DefaultMPEGFileName ; first frame number for GIF files if (n_elements(Frame0) eq 0) then Frame0 = 1 Frame0 = Frame0[0] > 1 ; auto-range: set to -1 if GAMAP shall decide itself ; force to zero if YRange is given if (n_elements(AutoRange) gt 0) then $ AutoRange = Keyword_Set( AutoRange ) $ else $ AutoRange = -1 if (n_elements(YRange) eq 2) then $ AutoRange = 0 $ else $ undefine,YRange ; make sure it doesn't interfere with anything ; see if tau value was passed or date if (n_elements(tau0) eq 0 AND n_elements(date) gt 0) then $ tau0 = nymd2tau(date,GISS=GISS) ;================================================================= ; get all data records (but don't read data yet) ;================================================================= ctm_get_data,datainfo,diagn,file=filename, $ tracer=tracer,tau0=tau0,status=2 if (n_elements(datainfo) eq 0) then begin message,'No matching records found !',/Continue return endif ;================================================================= ; Call CTM_GET_DATA to get the global DATAINFO structure ; List available diagnostics with item number on the screen ;================================================================= ctm_print_datainfo,datainfo,output=r,/NOPRINT ns = string(sindgen(n_elements(r)),format='(i4," : ")') ns[0] = ' ' for i=0,n_elements(r)-1 do print,ns[i],r[i] ;================================================================= ; Prompt user to select a data block or data block range ; If more than one data block is selected GAMAP switches to ; either multipanel or animation mode depending on the value ; of !p.multi. ;================================================================= Sel = GAMAP_SelectDataBlocks( N_Elements( R ), Quit=Quit ) ; ### debug ; print,'### sel=',sel if ( Quit ) then return if (n_elements(Sel) gt 1) then begin if (!p.multi[1]*!p.multi[2] gt 1) then Do_MultiPanel = 1 $ else Do_Animation = 1 endif ;================================================================= ; For the diagnostic selected, access the proper records from ; the DATAINFO structure, and store them in D. Also access ; the global FILEINFO structure from the common block. ;================================================================= D = DataInfo[ Sel ] FileInfo = *( pGlobalFileInfo ) ;================================================================= ; For each element of D, find the corresponding element of ; FILEINFO. Return an index array of matches, FILEINDEX. ; Also store a pointer to the GRIDINFO structure for each of ; the matching elements of FILEINFO. ;================================================================= FileIndex = Make_Selection(FileInfo.Ilun, D.Ilun, /REQUIRED) if ( FileIndex[0] lt 0 ) then return GAMAP_StoreGridInfo, FileInfo, FileIndex, Quit=Quit if ( Quit ) then return ;================================================================= ; For animation, check to make sure that all of the data ; blocks have the same dimensions and model grids ;================================================================= if ( Do_Animation ) then begin Success = GAMAP_CheckDataBlockConsistency( D, FileInfo, FileIndex ) if ( not Success ) then return endif ; get number of final records to plot NRec = N_Elements(D) ;================================================================= ; Compute the "natural" extent of the data block in longitude, ; latitude and vertical levels. Since all elements of D now ; must have the same range, pass the DATAINFO and FILEINFO ; structures corresponding to the first data block. ;================================================================= Success = GAMAP_GetDataBlockRanges( D[0], FileInfo[ FileIndex[0] ], $ TmpLon_Range, TmpLat_Range, TmpLev_Range ) if ( not Success ) then return ;**** Debug output (bmy, 2/19/99) ;* DEBUG=1 if ( Debug ) then begin print, '### GAMAP.PRO: internal Lon_Range: ', TmpLon_Range print, '### GAMAP.PRO: internal Lat_Range: ', TmpLat_Range print, '### GAMAP.PRO: internal Lev_Range: ', TmpLev_Range endif ;* DEBUG=0 ;================================================================= ; Check consistency of XXX_Sav fields for geographical domain ; If there is any change in either LONs or LATs reset both ;================================================================= ; level: just need to crop Lev_Sav if data records contain less if (Lev_Sav[0] lt TmpLev_Range[0]) then Lev_Sav[0] = TmpLev_Range[0] if (Lev_Sav[1] gt TmpLev_Range[1]) then Lev_Sav[1] = TmpLev_Range[1] RangeNeedsChange = 0 ; latitude: likewise if (Lat_Sav[0] lt TmpLat_Range[0]) then RangeNeedsChange = 1 if (Lat_Sav[1] gt TmpLat_Range[1]) then RangeNeedsChange = 1 ; longitude: need to reset if data is WE and ranges are EW or vice versa isEW = Lon_Sav[0] gt Lon_Sav[1] isTmpEW = TmpLon_Range[0] gt TmpLon_Range[1] if (isEW + isTmpEW eq 1) then $ RangeNeedsChange = 1 $ else begin ; check if data block fits in saved ranges if (not isEW AND Lon_Sav[0] lt TmpLon_Range[0]) then RangeNeedsChange = 1 if (not isEW AND Lon_Sav[1] gt TmpLon_Range[1]) then RangeNeedsChange = 1 endelse if (RangeNeedsChange) then begin Lon_Sav[0] = TmpLon_Range[0] Lon_Sav[1] = TmpLon_Range[1] Lat_Sav[0] = TmpLat_Range[0] Lat_Sav[1] = TmpLat_Range[1] endif ;**** Debug output (bmy, 2/19/99) ;* DEBUG=1 if ( Debug ) then begin print, '### GAMAP.PRO: Lon_Sav : ', Lon_Sav print, '### GAMAP.PRO: Lat_Sav : ', Lat_Sav print, '### GAMAP.PRO: Lev_Sav : ', Lev_Sav endif ;* DEBUG=0 ;================================================================= ; select geographical region ; ; Enter level, longitude, and latitude range ; If user presses return, use default range ; ; NOTE: Use floating point format for LON and LAT (bmy, 2/19/99) ;================================================================= Try_Again: ; Get vertical levels Level = fix( GAMAP_UserRangeEntry( 'level', Lev_Sav, Format='(i10)', $ Range=[1,100], Quit=Quit ) ) if ( Quit ) then return Print, Level, Format='(" Vertical Levels: ",2i10 )' ; Get longitudes Lon = GAMAP_UserRangeEntry( 'longitude', Lon_Sav, $ Range=TmpLon_Range[sort(TmpLon_Range)], $ /NOSORT, Quit=Quit ) if ( Quit ) then return ; Reset LON to the nearest grid box longitude edges Lon = GAMAP_FindNearestEdges( FileInfo[ FileIndex[0] ], Lon, /Longitude ) Print, Lon, Format='(" Nearest Grid Box Longitude Edges: ",2f10.2 )' ; Get latitudes...user will be warned about half-polar boxes! Lat = GAMAP_UserRangeEntry( 'latitude', Lat_Sav,$ Range=TmpLat_Range, Quit=Quit) if ( Quit ) then return ; print, '### Lat : ', Lat ; Reset LAT to the nearest grid box latitude edges Lat = GAMAP_FindNearestEdges( FileInfo[ FileIndex[0] ], Lat, /Latitude ) Print, Lat, Format='(" Nearest Grid Box Latitude Edges: ",2f10.2 )' ;================================================================= ; Check the number of unique elements in LON, LAT, and LEVEL. ; Here, NLON, NLAT, NLEVEL will either be 1 or 2 ; ; NOTE the following weakness: If a user enters two different ; values that belong to one grid box, GAMAP will think there ; will be more than one value for for this dimension and ask for ; averaging or totaling. The same may happen if index offsets ; are used (which will be supported in the next version). ;================================================================= NLon = N_Uniq( Lon ) NLat = N_Uniq( Lat ) NLevel = N_Uniq( Level ) ;================================================================= ; DIMSTR contains descriptor strings for each dimension. ; ; ARR is an array of dimension information flags: ; ARR[i] = 0 -> Dimension [i] has been averaged/totaled ; ARR[i] = 1 -> There is a single point along Dimension [i] ; ARR[i] > 1 -> There are multiple points along Dimension [i] ; ; NDIM is the number of dimensions of the plot. ;================================================================= Arr = [ NLon, NLat, NLevel ] GAMAP_PrintDimInfo, Arr, NDim=NDim ;==================================================================== ; If at least one dimension contains more than one value, ; query user for averaging (or totaling). ; ; For safety...if multiple 3-D ata blocks are selected, make sure ; only to vizualize the first data block. Print a message. ;==================================================================== GAMAP_QueryAverageOrTotal, Arr, NDim, AvMode, TotalMode, Quit=Quit if ( Quit ) then return if ( N_Elements( D ) gt 1 AND NDim eq 3 ) then begin Message, 'Multiple 3-D data blocks selected...' + $ 'GAMAP will vizualize only the first data block', /Continue D = D[0] endif ;==================================================================== ; For 0-D, 1-D and 3-D plots, print an informational message. ; For 2-D plots, ask user for type of plot (contour or pixel) ; ; NOTE: PTYPE is now returned as -1 unless redefined ; QUIT = 1 is now the criteria for exiting GAMAP (bmy, 2/16/99) ;==================================================================== GAMAP_SelectPlotType, $ NDim, NLon, NLat, NLevel, AvMode, TotalMode, PType, Quit=Quit ; Selection was aborted if ( Quit ) then return ; Contour Plot Options if ( PType ge 0 ) then begin MContour = ( PType le 1 ) FContour = ( PType eq 2 ) Sample = ( PType eq 4 ) if ( PType eq 0 ) then C_Colors = 1 endif ;================================================================= ; print units in D and query whether unit conversion shall be ; made. A unit conversion will attempt to convert the units of ; *ALL* seelcted records to the new unit. ;================================================================= print if (NRec gt 9) then sx = '...' else sx = '' tmps = string(strtrim(d[0:(nrec-1) < 9].unit,2), sx, format='(9(A,2X),A)') print, 'Units : ', tmps ustr = '' read,ustr,prompt='Enter new unit for all (Q=Quit, default : don''t touch) >> ' ustr = strtrim(ustr, 2) if (strupcase(ustr) eq 'Q') then return DoNot_ConvertUnit = (ustr eq '') print,'$$$ DoNotConvertUnit = ',DoNot_ConvertUnit ;================================================================= ; Get the minimum and maximum values over the entire set of ; animation frames or panels. This is necessary if we are to ; ensure that the color indices in each frame cover the same ; range of data values. This option is forced by setting the ; AutoRange flag to 0 or 1. ;================================================================= ; print,'### Autorange = ',autorange if ( AutoRange lt 0 ) then AutoRange = Do_Animation ; print,'### Autorange = ',autorange if ( AutoRange ) then begin AutoYRange = GAMAP_AutoYRange( D, AvMode, TotalMode, $ Lon, Lat, Level ) ; print,'### after AutoYRange : D.Status = ',D.Status ; print,'### AutoYrange = ',autoYrange ; update datainfo[sel] because data may have been read ; doesn't work but we should do something like this to ; prevent multiple reading of the same records !! ***** #### ; ** if (n_elements(d) eq n_elements(sel)) then DataInfo[Sel] = D endif ;================================================================= ; Prompt user for animation options. ;================================================================= if (not Do_Animation) then Do_MPEG = 0 ; doesn't make sense GAMAP_QueryAnimationOptions, NDim, $ Do_GIF, GIFFileName, Do_MPEG, $ MPEGFileName, Quit=Quit if ( Quit ) then return ;================================================================= ; Get the color table vectors for WRITE_GIF ;================================================================= ; myct,23 ; #### ; if (!d.name eq 'WIN') then $ ; myct,27,ncol=32,bot=20,range=[0.2,0.8] TvLCT, R_Cur, G_Cur, B_Cur, /Get R_Cur[0] = 255 & G_Cur[0] = 255 & B_Cur[0] = 255 R_Cur[1] = 0 & G_Cur[1] = 0 & B_Cur[1] = 0 TVLct, R_Cur, G_Cur, B_Cur ;================================================================= ; Create the plot! DO_PS controls whether or not to produce a ; postscript file: '0'=never, '1'=always, '*QUERY'=ask user ; The /PS keyword overrides the default setting in gamap.defaults. ; No postscript output can be created from 0-D data, 3-D data, ; or timeseries animation. ;================================================================= Do_The_Plot: if (Do_Animation OR NDim eq 0 OR NDim eq 3) then $ Do_PS = 0 $ else $ GAMAP_QueryPostScriptOptions, $ Do_PS, TimeStamp, OutFileName ;============================================================== ; Open output device (postscript file or window) ; _EXTRA options for example: /PORTRAIT, YSIZE, etc. ;============================================================== already_ps = (!D.Name eq 'PS') if (already_ps AND Do_Animation) then begin message,'Cannot do animation on postscript device! Reset to screen.',$ /Continue Close_Device already_ps = 0 endif if (not already_ps) then $ Open_Device, OLD_DEVICE, PS=(Do_PS eq 1), $ FileName=OutFileName, /Color, _EXTRA=e ;================================================================= ; Loop over number of data blocks selected ;================================================================= for N = 0, N_Elements( D ) - 1 do begin ;============================================================== ; Call CTM_PLOT to make the plot for each data block ;============================================================== if ( AutoRange ) then YRange = AutoYRange if (DoNot_ConvertUnit) then ustr = D[N].unit ; if ( AutoRange ) then begin ; If we are not animating, or if we are animating but ; NOAUTOYRANGE=0, then call CTM_PLOT without a YRANGE keyword CTM_Plot, D[N].Category, $ tracer=D[N].tracer, use_datainfo=datainfo, $ Tau0=D[N].Tau0, Ilun=D[N].Ilun, $ lev=level, lon=lon, $ lat=lat, total=totalmode, $ average=avmode, $ contour=mcontour, fcontour=fcontour, $ c_colors=c_colors, $ sample=sample, unit=ustr, $ YRange=YRange, $ _EXTRA=e ; endif else begin ; ; Otherwise call CTM_PLOT with YRANGE=AUTOYRANGE ; ; which was computed above... ; CTM_Plot, D[N].Category, $ ; tracer=D[N].tracer, use_datainfo=datainfo, $ ; Tau0=D[N].Tau0, $ ; lev=level, lon=lon, $ ; lat=lat, total=totalmode, $ ; average=avmode, $ ; contour=mcontour, fcontour=fcontour, $ ; c_colors=c_colors, $ ; sample=sample, unit=ustr, $ ; _EXTRA=e ; endelse ;============================================================== ; Define a byte array to hold animation frames ;============================================================== if ( !D.Name ne 'PS' ) then begin if (!d.name eq 'WIN') then begin ; *** USE CLEANER test !!! *** ; readtrue color image and quantize to 256 bits for gif and XInternimate Image = TVRD(0,0,!D.X_SIZE, !D.Y_SIZE, true=3) print,'###>>:',Image[0,0,*] if (N eq 0) then begin ThisFrame = Color_Quan(Image,3,R_Cur,G_Cur,B_Cur,colors=216, $ CUBE=6) ; /DITHER,get_translation=cttrans) endif else begin ThisFrame = Color_Quan(Image,3,R_Cur,G_Cur,B_Cur,colors=216, $ CUBE=6) ; /DITHER,translation=cttrans) endelse ind = where(R_Cur eq Image[0,0,0] AND G_Cur eq Image[0,0,1] AND B_Cur eq Image[0,0,2]) print,'###>>>>>:',ind[0] R_Cur[0]=0 ; 255 G_Cur[0]=0 ; 255 B_Cur[0]=0 ; 255 R_Cur[1]=255 G_Cur[1]=255 B_Cur[1]=255 endif else $ ThisFrame = TVRD(0,0,!D.X_SIZE, !D.Y_SIZE) help,thisframe ; stop endif ;============================================================== ; If we are producing GIF output, write the current ; frame to a GIF file. Produce 1 file per frame. ;============================================================== if ( Do_GIF ) then begin FrameN = N + Frame0 NumText = StrTrim( String( FrameN, Format='(i5.5)' ), 2 ) TmpFileName = Replace_Token( GIFFileName, '%NNNNN%', NumText ) if (FrameN lt 10000) then begin NumText = StrTrim( String( FrameN, Format='(i4.4)' ), 2 ) TmpFileName = Replace_Token( TmpFileName, '%NNNN%', NumText ) endif if (FrameN lt 1000) then begin NumText = StrTrim( String( FrameN, Format='(i3.3)' ), 2 ) TmpFileName = Replace_Token( TmpFileName, '%NNN%', NumText ) endif if (FrameN lt 100) then begin NumText = StrTrim( String( FrameN, Format='(i2.2)' ), 2 ) TmpFileName = Replace_Token( TmpFileName, '%NN%', NumText ) endif if (FrameN lt 10) then begin NumText = StrTrim( String( FrameN, Format='(i1.1)' ), 2 ) TmpFileName = Replace_Token( TmpFileName, '%N%', NumText ) endif Write_GIF, TmpFileName, ThisFrame, R_Cur, G_Cur, B_Cur, _EXTRA=e endif ;============================================================== ; If we are producing MPEG output, write the current frame ; to a MPEG file. If this is the first frame, then call ; MPEG_OPEN to open the MPEG file and assign it an ID. ;============================================================== if ( Do_MPEG ) then begin if ( N eq 0 ) then begin Dimensions = Size( ThisFrame, /Dimensions ) MpegID = MPEG_Open( Dimensions, FileName=MPEGFileName ) endif MPEG_Put, MpegID, Image=ThisFrame, Frame=N+1 endif ;============================================================== ; If we are doing on-screen animation, save all N animation ; frame into a single byte array for XINTERANIMATE ;============================================================== if ( Do_Animation ) then begin if ( N eq 0 ) then begin Frame = BytArr( !D.X_SIZE, !D.Y_SIZE, N_Elements( D ) ) endif Frame[ *, *, N ] = ThisFrame endif endfor ;============================================================== ; Close the device, with timestamp, if needed. ;============================================================== if (not already_ps) then $ Close_Device, OLD_DEVICE, TimeStamp=TimeStamp, $ FileName=OutFileName, _EXTRA=e ;================================================================= ; To save on resources, undefine THISFRAME and DATA ;================================================================= Undefine, ThisFrame Undefine, Data ;================================================================= ; If we have just plotted a 0-D or 3-D data block, quit GAMAP. ;================================================================= if ( NDim eq 0 OR NDim eq 3 ) then return ;================================================================= ; If we are saving MPEG output, close the MPEG file ;================================================================= if ( Do_MPEG ) then begin MPEG_Save, MpegID, FileName=MPEGFileName MPEG_Close, MpegID endif ;================================================================= ; If we are doing on-screen animation, call XINTERANIMATE ; to display the sequence of frames stored in the FRAME array ;================================================================= if ( Do_Animation ) then begin XInterAnimate, Set=[ !D.X_SIZE, !D.Y_SIZE, N_Elements( D ) ], /Track for N = 0, N_Elements( D ) - 1 do begin XInterAnimate, Frame=N, Image=Frame[ *, *, N ] endfor XInterAnimate, 10, _EXTRA=e endif ;================================================================= ; If DO_PS has the *QUERY value, ask user whether to produce the ; same plot as postscript file. Otherwise quit ;================================================================= if (Do_PS lt 0) then begin Do_PS = yesno('Create postscript file?',default=0) if (Do_PS) then begin Do_GIF = 0 Do_MPEG = 0 goto,Do_The_Plot endif endif return end