;+
;NAME:
; PLOT3D_NEW
;PROCEDURE: plot3d_new, dat , latitude, longitude
;PURPOSE:
; Draws a series of 3d color plots, one plot per energy step.
;INPUT:
; dat: 3d data structure.
; latitude: latitude of center of plot
; longitude: longitude of center of plot
;KEYWORDS:
; EBINS: Specifies energy bins to plot.
; SUM_EBINS: Specifies how many bins to sum, starting with EBINS. If
; SUM_EBINS is a scalar, that number of bins is summed for
; each bin in EBINS. If SUM_EBINS is an array, then each
; array element will hold the number of bins to sum starting
; at each bin in EBINS. In the array case, EBINS and SUM_EBINS
; must have the same number of elements.
; BNCENTER: if > 0 then lat,lon set so that B direction is in center
; if < 0 then lat,lon set so that -B direction is in center
; ('magf' element must exist in dat structure. See "ADD_MAGF")
; BINS: bins to use for autoscaling.
; SMOOTH: smoothing parameter. 0=none, 2 is typical
; TITLE: overrides default plot title.
; NOTITLE: if set, suppresses plot title.
; NOCOLORBAR: if set, suppresses colorbar.
; NOBORDER: if set, suppresses plot border.
;SEE ALSO: "PLOT3D_OPTIONS" to see how to set other options.
;
;;NOTE: plot3d_new has been created to eliminate a naming conflict between plot3d and a library routine that was added to IDL in version 8.1.
; The original plot3d.pro remains to preserve backwards compatibility for users using IDL 7.1 or earlier.
;
;CREATED BY: Davin Larson
;LAST MODIFICATION: @(#)plot3d.pro 1.43 02/04/17
;-
pro plot3d_new,tbindata, latitude,longitude,rotation, $
BINS=bins, $
TBINS=tbins, $
BNCENTER=bncenter, $
LABELS=labs, $
SMOOTH= smoth, $
EBINS=eb, $
SUM_EBINS=seb, $
PTLIMIT=plot_limits, $
MAP= ptmap, $
UNITS = units, $
SETLIM = setlim, $
PROC3D = proc3d, $
ZRANGE = zrange, $
ZERO = zero, $
STACK = stack, $ ;!p.multi[1:2]
ROW_MAJOR = row_major,$ ;!p.multi[4]
NOERASE = noerase, $ ;!p.multi[0]
TITLE = title, $
NOTITLE = notitle, $
NOCOLORBAR = nocolorbar, $
NOBORDER = noborder,$
log=log
@plot3d_com.pro
plot3d_options ; Set Defaults
if size(/type,tbindata) ne 8 || tbindata.valid eq 0 then begin
dprint,'Invalid data'
return
endif
if not keyword_set(units) then units='Eflux'
bindata = conv_units(tbindata,units)
str_element,bindata,'magf',value=magf
if keyword_set(magf) then begin
xyz_to_polar,magf,theta=th,phi=phi
bdir = [th,phi]
endif
str_element,bindata,'vsw',value=vsw
if n_elements(latitude) eq 0 then latitude = latitude_3d
if n_elements(longitude) eq 0 then longitude = longitude_3d
if n_elements(rotation) eq 0 then rotation = 0d
if n_elements(smoth) eq 0 then smoth = smooth_3d
if n_elements(log) eq 0 then log = logscale_3d
if keyword_set(bncenter) and n_elements(bdir) eq 2 then begin
latitude = bdir(0)
longitude = bdir(1)
if bncenter lt 0 then begin
latitude = - latitude
longitude = longitude + 180
endif
endif
n_e = bindata.nenergy
nbins = bindata.nbins
if keyword_set(bins) eq 0 then str_element,bindata,'bins',bins
if n_elements(bins) eq (nbins*n_e) then begin
if n_elements(tbins) ne nbins then tbins = total(bins,1) gt 0
bins = total(bins,2) gt 0
binsindx = where(tbins eq 0,binscnt)
if binscnt eq 0 then tbins = 0
binsindx = where(bins eq 0,binscnt)
if binscnt eq 0 then bins = 0
endif
if n_elements(bins) ne nbins then bins = replicate(1,bindata.nbins)
goodbins = where(bins,ngood)
badbins = where(bins eq 0,nbad)
if nbad ne 0 then bindata.data(*,badbins) = !values.f_nan
triang = keyword_set(triang_3d)
str_element,bindata,'map',value=ptmap
if triang eq 0 then begin
if keyword_set(ptmap) then begin
np = dimen1(ptmap)
ptmap = rotate(ptmap,7) ; this line and the next should be removed when
ptmap = rotate(ptmap,5) ; the c-code is fixed
ptmap = reverse(ptmap,2) ; assumes !order = 0
dprint,'Using supplied mapping'
endif else begin
ares = 16
ptmap = make_3dmap(bindata,64*ares,32*ares)
ptmap = shift(ptmap,32*ares,0) ; map routine expects phi: -180 to 180
endelse
endif
; handle ebins
if n_elements(eb) eq 0 then $ ; undefined ebins
start=indgen(n_e) $
else $ ; scalar or 1-D ebins
start=eb
; handle sum_ebins
if keyword_set(seb) eq 0 then $ ; SUM_EBINS keyword not set
seb=1 $
else begin
if n_elements(seb) gt 1 and (n_elements(seb) ne n_elements(eb)) then begin
message,'if SUM_EBINS is an array, it must have same dimensions'+ $
' as EBINS.'
endif
endelse
; get stop bins
stop = start + seb - 1
n_p = dimen1(start) ; number of plots
stacking= [ [1,1], [1,1], [1,2], [2,2], [2,2], [2,3], [2,3], $
[3,3], [3,3], [3,3], [3,4], [3,4], [3,4], [3,5], [3,5], $
[3,5], [4,4], [4,5], [4,5], [4,5], [4,5], [5,5], [5,5], $
[5,5], [5,5], [5,5], [5,6], [5,6], [5,6], [5,6], [5,6], $
[6,6], [6,6], [6,6], [6,6], [6,6], [6,6], [6,7], [6,7], $
[6,7], [6,7], [6,7], [6,7], [7,7], [7,7], [7,7], [7,7], $
[7,7], [7,7], [7,7], [7,8], [7,8], [7,8], [7,8], [7,8], $
[7,8], [7,8], [8,8], [8,8], [8,8], [8,8], [8,8], [8,8], $
[8,8], [8,8]]
if n_elements(stack) ne 2 then stack=stacking(*,n_p)
if not keyword_set(noerase) then noerase = 0 else noerase = 1
if not keyword_set(row_major) then row_major = 0 else row_major = 1
!p.multi=[noerase,stack(0),stack(1),0,row_major]
if !p.charsize ne 0 then charsize = !p.charsize else charsize = 1
!y.omargin = [1,4]
!x.omargin = [0,6]
if !p.multi(1) gt 2 or !p.multi(2) gt 2 then begin
; charsize = !p.charsize *2
!y.omargin = !y.omargin * 2
!x.omargin = !x.omargin * 2
endif
xposmax = 0.
yposmax = 0.
yposmin = 1.
for plot_num=0,n_p-1 do begin
e1 = start(plot_num)
en1 = total(bindata.energy(e1,*))/n_elements(bindata.energy(e1,*))
e2 = stop(plot_num)
en2 = total(bindata.energy(e2,*))/n_elements(bindata.energy(e2,*))
if e1 ne e2 then begin
ttl = string(en1)+' - '+string(en2)+' eV'
endif else begin
ttl = strcompress(string(en1,format='(F12.1)')+' eV', $
/remove_all)
endelse
ttl=strcompress(ttl)
nsteps = e2-e1+1
mapname = strlowcase(strmid(mapname_3d,0,3))
if keyword_set(setlim) then begin
map_set, latitude,longitude,0., /WHOLE_MAP, $
/noerase, $
MOLL = mapname eq 'mol' , $
CYLIN = mapname eq 'cyl' , $
ORTHO = mapname eq 'ort' , $
AITOFF= mapname eq 'ait' , $
LAMBERT=mapname eq 'lam' , $
GNOMIC =mapname eq 'gno' , $
MERCATO=mapname eq 'mer' , $
LIMIT=plot_limits , $
TITLE=ttl, $
charsize = charsize, $
XMARGIN=[1,1], YMARGIN=[1,1], noborder=noborder
return
endif
data = total(bindata.data(e1:e2,*),1) ; sum over energies
if keyword_set(zrange) then range = zrange[*] $
else begin
range = minmax(data(goodbins),/pos)
endelse
if keyword_set(triang) then begin
; data = total(bindata.data(e1:e2,*),1) ; sum over energies
; range = minmax(data(goodbins),/pos)
theta = total(bindata.theta(e1:e2,*),1)/nsteps
phi = total(bindata.phi(e1:e2,*),1)/nsteps
if keyword_set(tbins) then ind = where(tbins ne 0) $
else ind = where(finite(theta) and finite(phi))
if ind(0) eq -1 then goto, contu
data = data(ind)
theta = theta(ind)
phi = phi(ind)
triangulate,phi,theta,fvalue=data,sphere=sphere,/degrees
image = trigrid(data,sphere=sphere,[2.,2.],[-180.,-90.,178.,90.],/deg)
image = bytescale(image,log=log, range=range,zero=zero)
endif else begin
; data = total(bindata.data(e1:e2,*),1) ; sum over energies
; range = minmax(data(goodbins))
data = [!values.f_nan,data]
; data = [0.,data]
data = data(ptmap+1) ; make 2-dimensional
if( keyword_set(smoth)) then begin
dim = dimen(data) * 1*smoth ; increase the dimensions
data = rebin(data,dim(0),dim(1))
data = smooth_periodic(data,smoth)
endif
image = bytescale(data,log=log, range=range,zero=zero)
image = congrid(image,360,180) ; assumes !order eq 0
endelse
mapname = strlowcase(strmid(mapname_3d,0,3))
map_set, latitude,longitude,rotation, /advance, /WHOLE_MAP, $
MOLL = mapname eq 'mol' , $
CYLIN = mapname eq 'cyl' , $
ORTHO = mapname eq 'ort' , $
AITOFF= mapname eq 'ait' , $
LAMBERT=mapname eq 'lam' , $
GNOMIC =mapname eq 'gno' , $
MERCATO=mapname eq 'mer' , $
LIMIT=plot_limits , $
TITLE=ttl, $
charsize = charsize, $
XMARGIN=[1,1], YMARGIN=[1,1], noborder = noborder
new = map_image(image,sx,sy,xsize,ysize, $
COMPRESS=compress_3d, BILIN=0)
; new = map_patch(image,xstart=sx,ystart=sy,xsize=xsize,ysize=ysize)
tv,new,sx,sy,xsize=xsize,ysize=ysize,order=0
if n_elements(bdir) eq 2 then begin ; draw mag field
plots,bdir(1),bdir(0),psym=1
plots,bdir(1)+180,-bdir(0),psym=4
endif
if n_elements(vsw) then begin
xyz_to_polar,vsw,theta=th,phi=ph
plots,ph,th,psym=2
endif else plots,179.99,0.,psym=2 ; draw the sun.
if keyword_set(grid_3d) then map_grid,londel=grid_3d(0),latdel=grid_3d(1)
str_element,opts,'plot3d_proc',value=proc3d
str_element,bindata,'plot3d_proc',value=proc3d
if keyword_set(proc3d) then begin
proc = ''
if ndimen(proc3d) eq 0 then proc=proc3d else begin
if plot_num lt n_elements(proc3d) then proc = proc3d(plot_num)
endelse
if keyword_set(proc) then call_procedure,proc,bindata,en2,plot_num
endif
if keyword_set(labs) then begin
lab=strcompress(indgen(nbins),/rem)
colorcode = [!p.background,!p.color]
xyouts,phi,theta,lab,align=.5,COLOR=colorcode(bins)
endif
;map_grid,/label,latdel =45,londel=45,CHARSIZE=charsize
;map_continents
if xposmax lt !x.window(1) then xposmax=!x.window(1)
if yposmax lt !y.window(1) then yposmax=!y.window(1)
if yposmin gt !y.window(1) then yposmin=!y.window(0)
contu:
endfor
if keyword_set(notitle) then ttl = '' else $
if not keyword_set(title) then begin
str_element,bindata,'trange',trange
; if not keyword_set(trange) then $
trange=[bindata.time,bindata.end_time]
ttl = bindata.project_name+' '+bindata.data_name
; ttl =ttl+'!C'+ trange_str(bindata.time,bindata.end_time)
ttl =ttl+'!C'+ trange_str(trange)
endif else ttl = title
tcharsize = 1.25 * charsize
titleheight = tcharsize * float(!d.y_ch_size)/!d.y_size
xyouts, 0.5,1.-titleheight,ttl,/normal,alignment=.5,charsize=tcharsize
time_stamp
!p.multi = 0
!y.omargin=0
!x.omargin=0
space = charsize * float(!d.x_ch_size)/!d.x_size
pos =[xposmax+space,yposmin,xposmax+3*space,yposmax]
if not keyword_set(nocolorbar) then $
draw_color_scale,range = [0.,1.],pos=pos,chars=charsize
return
end