PRO DF_FLUX_PLOT,yyyydddb,hhb,yyyyddde,hhe,dtmin=dtmin,fmin=fmin,fmax=fmax,$
vmin=vmin,vmax=vmax,snrmin=snrmin,betamin=betamin,no_rfi=no_rfi,$
view=view,latmin=latmin,titre=titre,output_path=output_path,$
smin=smin,smax=smax,$
dB=dB,occurrence=occurrence,position=position,$
pts_limbe=pts_limbe,err_bar=err_bar,horizon=horizon,mfl=mfl,$
spdyn_flux=spdyn_flux,spdyn_polar=spdyn_polar,high_res=high_res,$
add_sls=add_sls,mark_spdyn=mark_spdyn
if ~keyword_set(dtmin) then dtmin=(hhe-hhb)*60.
if ~keyword_set(vmin) then vmin=0. else vmin=abs(float(vmin))
if ~keyword_set(vmax) then vmax=1.1 else vmax=abs(float(vmax))
if vmax lt vmin then begin
v1=vmax & v2=vmin
vmin=v1 & vmax=v2
endif
if ~keyword_set(fmin) then fmin=100. else fmin=float(fmin)
if ~keyword_set(fmax) then fmax=1000. else fmax=float(fmax)
if ~keyword_set(dtmin) then dtmin=60.d0 else dtmin=double(dtmin)
if ~keyword_set(snrmin) then snrmin=20.
if ~keyword_set(betamin) then betamin=20.
if ~keyword_set(titre) then titre='SKR_map'
if ~keyword_set(output_path) then output_path=''
if ~keyword_set(view) then view='90'
if ~keyword_set(latmin) then latmin=50.
if ~keyword_set(smin) then smin=-230.
if ~keyword_set(smax) then smax=-215.
smin1=smin & smax1=smax
smin2=smin & smax2=smax
debut = long(aj_t97(yyyydddb))
fin = long(aj_t97(yyyyddde))
reso = 2.
nerr = 8.
source='sq'
read_data_src_list,src_list
src_list = (src_list(where(src_list.df_name eq source)))(0)
ua = 150.e6/src_list.body_rad_value
set_plot,'ps'
device,filename=output_path+titre+'.ps',/color,/landscape,bits=8
!p.font=0
if keyword_set(spdyn_flux) or keyword_set(spdyn_polar) then begin
if keyword_set(high_res) then begin
high_res_plot,yyyydddb,0,yyyyddde,24,ssi,vvi,tt,ff,fminmax=[3,1500],add_sls=add_sls,$
vmin=0.,vmax=2.,snrmin=0.,betamin=betamin,/no_rfi,/no_plot,smin=1e-23
endif else begin
high_res_plot,yyyydddb,hhb,yyyyddde,hhe,ssi,vvi,tt,ff,fminmax=[3,1500],add_sls=add_sls,$
vmin=vmin,vmax=vmax,snrmin=snrmin,betamin=betamin,no_rfi=no_rfi,/no_plot,smin=1e-23
endelse
endif
for jour=debut,fin do begin
if fin*1440+hhe*60-debut*1440-hhb*60 le dtmin then begin
jour_ajb=yyyydddb & h1=hhb
jour_aje=yyyyddde & h2=hhe
endif else begin
jour_ajb=t97_aj(jour) & jour_aje=jour_ajb
h1=0.d0 & h2=24.d0
if jour eq debut then h1=double(hhb)
if jour eq fin then h2=double(hhe)
endelse
hh1=long(h1) & hh2=long(h2+0.999)
ptrDataFileList = make_file_list(jour_ajb,hh1, jour_aje,hh2, level='n2')
if (ptr_valid(ptrDataFileList) eq 0) then begin
n_n2DataFile=0
endif else begin
n2DataFileList=*ptrDataFileList
n_n2DataFile = n_elements(n2DataFileList)
ptr_free,ptrDataFileList
endelse
ptrDataFileList = make_file_list(jour_ajb,hh1, jour_aje,hh2, level='n3e', $
source='sq', antset='d')
if (ptr_valid(ptrDataFileList) eq 0) then begin
n_n3DataFile=0
print,"# Warning ! No N3e data file selected..."
endif else begin
n3eDataFileList=*ptrDataFileList
n_n3eDataFile = n_elements(n3eDataFileList)
ptr_free,ptrDataFileList
endelse
n2=[{data_n2}]
n3e=[{data_n3e}]
index=[{data_index}]
loc=[{data_loc}]
loc_error=[{data_loc}]
for i=0l,n_n3eDataFile-1l do begin
while strmid(n3eDataFileList(i),9,/rev) ne strmid(n2DataFileList(i),9,/rev) do begin
if i eq 0 then n2DataFileList=n2DataFileList(1:*) else $
n2DataFileList=[n2DataFileList(0:i-1),n2DataFileList(i+1:*)]
n_n2DataFile = n_n2DataFile-1
endwhile
read_data_binary,n2DataFileList(i),n2t,level='n2',valid
if valid eq 1 then begin
read_data_binary,n3eDataFileList(i),n3et,level='n3e'
n3e = [n3e,n3et]
n2 = [n2,n2t(n3et.num)]
period_path = strmid(n3eDataFileList(i),0,strlen(n3eDataFileList(i))-21)
str_aj=strmid(n3eDataFileList(i),9,/rev)
dataFileIn = period_path+'index/INDEX_'+str_aj
read_data_binary,dataFileIn,indext,level='index'
index = [index,indext]
if n_elements(n3et) eq 1 then begin
loct=[{data_loc}]
loc_errort=replicate([{data_loc}],nerr)
endif else begin
dataFileIn = period_path+'loc/loc_SPV_'+str_aj
read_data_binary,dataFileIn,loct,level='loc'
dataFileIn = period_path+'loc/loc_err_SPV_n'+string(format='(i2.2)',nerr)+'_'+str_aj
read_data_binary,dataFileIn,loc_errort,level='loc'
endelse
loc = [loc,loct]
loc_error=[loc_error,loc_errort]
endif
endfor
if n_elements(n2) eq 1 then return
n3e=n3e[1:*] & n2=n2[1:*] & index=index[1:*] & loc=loc[1:*]
if n_elements(loc) ne n_elements(n3e) then stop,'*** Warning *** Number of loc and n3e events incorrect'
loc_error=loc_error[1:*]
loc_error=reform(loc_error,nerr,n_elements(n3e))
wt=where(n2.t97 ge (long(n2(0).t97)+h1/24.d0) and n2.t97 lt (long(n2(n_elements(n2)-1l).t97)+h2/24.d0))
if wt(0) ne -1 then begin
n2=n2(wt)
n3e=n3e(wt)
loc=loc(wt)
index=index(wt)
loc_error=loc_error(*,wt)
endif else stop,'*** Intervalle horaire vide de donnees! ***'
f_hf=[(findgen(15)*100+393.750),(findgen(15)*100+406.250),525,1143.75,1225,1231.25,1743.75,$
1756.25,643.750,650.000,656.250,1343.75,1425,1718.75,(findgen(11)*100+400)]
f_parasites=bytarr(n_elements(n2))
for i=0,n_elements(f_hf)-1l do f_parasites=f_parasites+(n2.f eq f_hf(i))
f_bf=[21.4779,25.9205,101.302,200.265,205.027]
for i=0,n_elements(f_bf)-1l do f_parasites=f_parasites+(n2.f gt (f_bf(i)-0.1) and n2.f lt (f_bf(i)+0.1))
wf=where(f_parasites eq 0)
if keyword_set(no_rfi) and wf(0) ne -1 then begin
n2=n2(wf)
n3e=n3e(wf)
loc=loc(wf)
index=index(wf)
loc_error=loc_error(*,wf)
endif
foot = loc.foot
temp = where(index.ind eq 0)
if temp(0) ne -1 then n3e(temp).s=0.
foot(1,*) = 90.-foot(1,*)
foot(2,*) = (foot(2,*)-180.)*24./360.
temp=where(foot(2,*) lt 0) & if temp(0) ne -1 then foot(2,temp)=foot(2,temp)+24.
loc_error=loc_error.foot
loc_error(1,*,*)=90.-loc_error(1,*,*)
loc_error(2,*,*)=(loc_error(2,*,*)-180.)*24./360.
for i=0,n_elements(loc_error(0,*,0))-1l do begin
temp=where(loc_error(2,i,*) lt 0)
if temp(0) ne -1 then loc_error(2,i,temp)=loc_error(2,i,temp)+24.
endfor
sn=10*alog10(n3e.sn>1e-2)
rs=index.rr
tl=index.tl
lat = reform(asin(index.zz/index.rr)/!dtor)
xyz=fltarr(3,n_elements(n2.t97))
xyz(0,*)=-rs*cos(lat*!dtor)*cos(tl/24.*2*!pi)
xyz(1,*)=-rs*cos(lat*!dtor)*sin(tl/24.*2*!pi)
xyz(2,*)=rs*sin(lat*!dtor)
dr=reform(sqrt((xyz(0,*)-loc.xyz(0))^2+(xyz(1,*)-loc.xyz(1))^2+(xyz(2,*)-loc.xyz(2))^2))
load_data_ephem,jour_ajb,h1,jour_aje,h2,veph,qeph,veph='vsat',qeph='qssq',t97=n2.t97
qcenter = q_vmake(veph.r)
qc0 = q_rot(qeph.q,qcenter)
rr0 = sqrt(total(qc0(1:3,*)^2.,1))
qr1 = q_make(reform(-atan(qc0(2,*),qc0(1,*))),[0,0,1])
qr2 = q_make(reform(!pi/2-acos(qc0(3,*)/rr0)),[0,1,0])
qrot0 = q_prod(q_prod(qr2,qr1),qeph.q)
qcenter2= q_rot(qrot0,qcenter)
veph0 = veph
veph.r = qcenter2(1:3,*)
read_antenna_set,ant,'calDec04_H12',/rad
dthdip = angular_distance(veph0.r,make_vect_sph(1,ant.dip.al,ant.dip.be))
dthz = angular_distance(veph0.r,make_vect_sph(1,ant.z.al,ant.z.be))
dthxp = angular_distance(veph0.r,make_vect_sph(1,ant.xp.al,ant.xp.be))
dthxm = angular_distance(veph0.r,make_vect_sph(1,ant.xm.al,ant.xm.be))
dthdipz = fltarr(n_elements(veph0.r(0)))
wwdz = where(n2.ant eq 3)
if wwdz(0) ne -1 then begin
vv = crossP1(make_vect_sph(1,ant.dip.al,ant.dip.be),make_vect_sph(1,ant.z.al,ant.z.be),/norm)
dthdipz(wwdz) = !pi/2.-angular_distance(veph0(wwdz).r,vv)
endif
wwxpz = where(n2.ant mod 10 eq 1)
if wwxpz(0) ne -1 then begin
vv = crossP1(make_vect_sph(1,ant.xp.al,ant.xp.be),make_vect_sph(1,ant.z.al,ant.z.be),/norm)
dthdipz(wwxpz) = !pi/2.-angular_distance(veph0(wwxpz).r,vv)
endif
wwxmz = where(n2.ant mod 10 eq 2)
if wwxmz(0) ne -1 then begin
vv = crossP1(make_vect_sph(1,ant.xm.al,ant.xm.be),make_vect_sph(1,ant.z.al,ant.z.be),/norm)
dthdipz(wwxmz) = !pi/2.-angular_distance(veph0(wwxmz).r,vv)
endif
if fin*1440+hhe*60-debut*1440-hhb*60 le dtmin then nn=1 else nn=round((h2-h1)*60./dtmin)
for k=0,nn-1l do begin
if keyword_set(spdyn_flux) or keyword_set(spdyn_polar) then begin
!p.multi=[3,1,2]
if (max(tt) - min(tt)) le 1. and hhb lt hhe then begin
ti=(min(tt)-long(min(tt)))*24. & tf=(max(tt)-long(min(tt)))*24.
titrex = 'DOY '+strmid(t97_aj(min(tt)),7,4)+'-'+strmid(t97_aj(min(tt)),11,3)+' (h)'
endif else begin
ti=(min(tt)-aj_t97(long(strmid(t97_aj(min(tt)),7,4)+'000'))) & tf=(max(tt)-aj_t97(long(strmid(t97_aj(min(tt)),7,4)+'000')))
titrex = 'DOY '+strmid(t97_aj(min(tt)),7,4)
endelse
wtt=where(tt ge (jour+(k*dtmin+h1*60.)/1440.) and tt lt (jour+((k+1)*dtmin+h1*60.)/1440.),ntt)
wff=where(ff ge fmin and ff le fmax,nff)
sstemp=ssi
vvtemp=vvi
if keyword_set(mark_spdyn) and wtt(0) ne -1 and wff(0) ne -1 then begin
i1=wtt(0) & i2=wtt(ntt-1l)
j1=wff(0) & j2=wff(nff-1l)
sstemp((i1-1l)>0:i1,j1:j2)=max(ssi) & sstemp(i2:(i2+1l)<((size(ssi))(1)-1l),j1:j2)=max(ssi)
sstemp(i1:i2,j1:j1+1l)=max(ssi) & sstemp(i1:i2,j2-1l:j2)=max(ssi)
vvtemp((i1-1l)>0:i1,j1:j2)=1.1 & vvtemp(i2:(i2+1)<((size(vvi))(1)-1l),j1:j2)=1.1
vvtemp(i1:i2,j1:j1+1l)=1.1 & vvtemp(i1:i2,j2-1l:j2)=1.1
endif
if keyword_set(spdyn_flux) and ~keyword_set(spdyn_polar) then begin
spdynps,sstemp,ti,tf,3,1500,titrex,'Frequency (kHz)','Flux density',0,0,0,0,0.98,/log,$
posxmin=0.05,posxmax=0.95,posymin=0.62,posymax=1
endif else if keyword_set(spdyn_polar) and ~keyword_set(spdyn_flux) then begin
spdynps,vvtemp,ti,tf,3,1500,titrex,'Frequency (kHz)','Circular polarization',0,0,0,0,0.98,0,'.',/log,$
posxmin=0.05,posxmax=0.95,posymin=0.62,posymax=1
endif else begin
!p.multi=[4,1,3]
spdynps,sstemp,ti,tf,3,1500,'','','',0,0,0,0,0.98,/log,$
posxmin=0.05,posxmax=0.95,posymin=0.805,posymax=1.
spdynps,vvtemp,ti,tf,3,1500,'','','',0,0,0,0,0.98,0,'.',/log,$
posxmin=0.05,posxmax=0.95,posymin=0.6,posymax=0.795
xyouts,0.9815,0.62,'RH',charsize=.7,/normal
xyouts,0.9815,0.761,'LH',charsize=.7,/normal,col=255
xyouts,0.5,0.55,titrex,charsize=1.,/normal,align=0.5
xyouts,0.,0.8,'Frequency (kHz)',charsize=1.,/normal,orient=90,align=0.5
xyouts,0.06,0.815,'Flux density',charsize=1.,/normal
xyouts,0.06,0.61,'Circular polarization',charsize=1.,/normal,col=255
endelse
endif else !p.multi=[0,2,1]
ww=where(n2.t97 ge (jour+(k*dtmin+h1*60.d0)/1440.) and n2.t97 lt (jour+((k+1)*dtmin+h1*60.d0)/1440.),count)
if count gt 1 then begin
wwi=where(n2.f ge fmin and n2.f le fmax and abs(n3e.v) ge vmin and abs(n3e.v) le vmax and loc.dist eq 0 and foot(0,*) ne 0. and $
sn(0,*) ge snrmin and sn(1,*) ge snrmin and abs(dthdipz) ge (betamin*!dtor) and $
n2.t97 ge (jour+(k*dtmin+h1*60.d0)/1440.) and n2.t97 lt (jour+((k+1)*dtmin+h1*60.d0)/1440.),counti)
wwj=where(n2.f ge fmin and n2.f le fmax and abs(n3e.v) ge vmin and abs(n3e.v) le vmax and loc.dist gt 0 and foot(0,*) ne 0. and $
sn(0,*) ge snrmin and sn(1,*) ge snrmin and abs(dthdipz) ge (betamin*!dtor) and $
n2.t97 ge (jour+(k*dtmin+h1*60.d0)/1440.) and n2.t97 lt (jour+((k+1)*dtmin+h1*60.d0)/1440.),countj)
wwk=where(n2.f ge fmin and n2.f le fmax and abs(n3e.v) ge vmin and abs(n3e.v) le vmax and loc.dist ge 0 and foot(0,*) ne 0. and $
sn(0,*) ge snrmin and sn(1,*) ge snrmin and abs(dthdipz) ge (betamin*!dtor) and $
n2.t97 ge (jour+(k*dtmin+h1*60.d0)/1440.) and n2.t97 lt (jour+((k+1)*dtmin+h1*60.d0)/1440.),countk)
if keyword_set(spdyn_flux) or keyword_set(spdyn_polar) then begin
!p.multi(0)=2
!p.multi(2)=2
endif
if view eq '45' then lim=mean(rs(ww))/2. else $
if view eq '90' then lim=mean(rs(ww)) else $
if view eq '120' then lim=mean(rs(ww))*tan(60.*!dtor) else $
stop,'Choose field of view 45/90/120'
xmin=-lim & xmax=lim
ymin=-lim & ymax=lim
posxmin=0.0 & posymin=posxmin
posxmax=0.5 & posymax=posxmax
plot,[xmin,xmax],[ymin,ymax],xr=[xmin,xmax],yr=[ymin,ymax],xstyle=13,ystyle=13,$
xtit='',ytit='',tit='',charsize=1.,$
ticklen=0.2,/nodata,/iso,pos=[posxmin,posymin,posxmax,posymax]
xyouts,0.185,0.51,'View from Cassini (FOV = '+view+'!eo!n)',charsize=1.3,col=0,/normal,align=0.5
if countk gt 1 then begin
field_view=fltarr(1024,1024)
field_viewdt=field_view
fill_tab_error,field_view,index(wwk).ydf,index(wwk).zdf,$
n3e(wwk).s*n2(wwk).dt*n2(wwk).df*(dr(wwk)/ua)^2/1064.04,reso,nerr,lim,dB=dB
fill_tab_error,field_viewdt,index(wwk).ydf,index(wwk).zdf,n2(wwk).dt/1e3,reso,nerr,lim,/sum
if ~keyword_set(dB) then field_view = field_view/total(n2(wwk).dt*n2(wwk).df)
if keyword_Set(dB) then begin
smin1=min(field_view)
smax1=max(field_view)
tab=bytscl(field_view)
endif else begin
ws=where(field_view gt 0)
if ~keyword_set(set_smin) and ws(0) ne -1 then smin1=10*alog10(min(field_view(ws))) else smin1=smin
if ~keyword_set(set_smax) and ws(0) ne -1 then smax1=10*alog10(max(field_view(ws))) else smax1=smax
tab=bytscl(10*alog10(field_view>1e-40),min=smin1,max=smax1)
endelse
if keyword_set(occurrence) then begin
smin1=min(field_viewdt)
smax1=max(field_viewdt)
tab=bytscl(field_viewdt,min=smin1,max=smax1)
endif
tab=255-tab
if keyword_set(position) then begin
wwirh=[where(n3e(wwi).v lt 0.,compl=wwilh)]
wwjrh=[where(n3e(wwj).v lt 0.,compl=wwjlh)]
loadct,1
if wwilh(0) ne -1 then begin
oplot,[index(wwi(wwilh)).ydf],[index(wwi(wwilh)).zdf],psym=1,symsize=0.5,col=140
if keyword_set(pts_limbe) and countj gt 0 then $
oplot,[index(wwj(wwjlh)).ydf],[index(wwj(wwjlh)).zdf],psym=4,symsize=0.5,col=140
endif
loadct,3
if wwirh(0) ne -1 then begin
oplot,[index(wwi(wwirh)).ydf],[index(wwi(wwirh)).zdf],psym=1,symsize=0.5,col=140
if keyword_set(pts_limbe) and countj gt 0 then $
oplot,[index(wwj(wwjrh)).ydf],[index(wwj(wwjrh)).zdf],psym=4,symsize=0.5,col=140
endif
loadct,0
endif else $
tv,tab,!x.window(0),!y.window(0),xsize=!x.window(1)-!x.window(0),ysize=!y.window(1)-!y.window(0),/normal
endif
if keyword_set(mfl) and counti gt 0 then begin
tl_mfl=[mean(tl(ww))-3,mean(tl(ww))+3]
tl_mfl=[5,7,9] & lat_mfl=[-75,-75,-75]
for mfl=0,n_elements(tl_mfl)-1l do begin
LOAD_MFL,mean(lat(ww)),mean(tl(ww)),lat_mfl(mfl),tl_mfl(mfl),x,y
oplot,x,y,line=2,thick=2,col=150
endfor
endif
axis, xaxis=0, xra=[xmin,xmax], xstyle=1, xticklen=0.015, $
xtitle='R!iS!n', charsize = 1.1
axis, xaxis=1, xra=[xmin,xmax], xstyle=1, xticklen=0.015, $
xtickname=replicate(' ',9)
axis, yaxis=0, yra=[ymin,ymax], ystyle=1, yticklen=0.015, $
ytitle='R!iS!n', charsize = 1.1
axis, yaxis=1, yra=[ymin,ymax], ystyle=1, yticklen=0.015, $
ytickname=replicate(' ',9)
rdeg=tan(reso*!dtor)*lim
if keyword_set(position) then $
oplot,0.9*lim+rdeg*cos(findgen(30)/30.*2.1*!pi),0.9*lim+rdeg*sin(findgen(30)/30.*2.1*!pi),thick=2
loadct,3
grille_sat,alpha=0,beta=mean(lat(ww)),1.,0.098,!pi-mean(tl(ww))/24.*2*!pi,$
/plot,/equ,/mer_cml,/par,/meri
loadct,0
if keyword_set(spdyn_flux) or keyword_set(spdyn_polar) then !p.multi(0)=1
rmax=sqrt(1./(cos(latmin*!dtor)^2.+(sin(latmin*!dtor)/(1-0.098))^2.))*cos(latmin*!dtor)
xmin=-rmax & xmax=rmax
ymin=-rmax & ymax=rmax
posxmin=0.525 & posymin=0.025
posxmax=0.975 & posymax=0.475
plot,[xmin,xmax],[ymin,ymax],xr=[xmin,xmax],yr=[ymin,ymax],xstyle=13,ystyle=13,$
xtit='',ytit='',tit='',charsize=1.,$
ticklen=0.2,/nodata,/iso,pos=[posxmin,posymin,posxmax,posymax]
xyouts,0.69,0.51,'Magnetic polar projection',charsize=1.3,col=0,/normal,align=0.5
if counti gt 0 then begin
proj=fltarr(1024,1024)
proj_dt=proj
y=foot(0,wwi)*cos(foot(1,wwi)*!dtor)*cos(foot(2,wwi)/24.*2*!pi)
x=-foot(0,wwi)*cos(foot(1,wwi)*!dtor)*sin(foot(2,wwi)/24.*2*!pi)
err=fltarr(2,nerr,counti)
err(0,*,*)=-loc_error(0,*,wwi)*cos(loc_error(1,*,wwi)*!dtor)*sin(loc_error(2,*,wwi)/24.*2*!pi)
err(1,*,*)=loc_error(0,*,wwi)*cos(loc_error(1,*,wwi)*!dtor)*cos(loc_error(2,*,wwi)/24.*2*!pi)
fill_tab_error,proj,x,y,n3e(wwi).s*n2(wwi).dt*n2(wwi).df*(dr(wwi)/ua)^2/1064.04,err,nerr,1.,dB=dB
fill_tab_error,proj_dt,x,y,n2(wwi).dt/1e3,err,nerr,1.,/sum
if ~keyword_set(dB) then proj = proj/(total(n2(wwi).dt*n2(wwi).df))
rmax_pix=round(rmax*512)
for a=0l,1023l do begin
for b=0l,1023l do begin
if sqrt((a-512l)^2+(b-512l)^2) gt rmax_pix then proj(a,b)=0.
if sqrt((a-512l)^2+(b-512l)^2) gt rmax_pix then proj_dt(a,b)=0.
endfor
endfor
tab=proj(*,512-rmax_pix:511+rmax_pix) & tab=tab(512-rmax_pix:511+rmax_pix,*)
tabdt=proj_dt(*,512-rmax_pix:511+rmax_pix) & tabdt=tabdt(512-rmax_pix:511+rmax_pix,*)
if keyword_set(dB) then begin
smin2=min(tab)
smax2=max(tab)
tab=bytscl(tab,smin2,smax2)
endif else begin
ws=where(tab gt 0)
if ~keyword_set(set_smin) and ws(0) ne -1 then smin2=10*alog10(min(tab(ws))) else smin2=smin-20.
if ~keyword_set(set_smax) and ws(0) ne -1 then smax2=10*alog10(max(tab(ws))) else smax2=smax-20.
tab=bytscl(10*alog10(tab>1e-40),min=smin2,max=smax2)
endelse
if keyword_set(occurrence) then begin
smin2=min(tabdt)
smax2=max(tabdt)
tab=bytscl(tabdt,min=smin2,max=smax2)
endif
tab=255-tab
if keyword_set(position) then begin
ypol=reform(foot(0,*)*cos(foot(1,*)*!dtor)*cos(foot(2,*)/24.*2*!pi))
xpol=reform(-foot(0,*)*cos(foot(1,*)*!dtor)*sin(foot(2,*)/24.*2*!pi))
if wwilh(0) ne -1 then begin
loadct,1
oplot,[xpol(wwi(wwilh))],[ypol(wwi(wwilh))],col=140,psym=1
if keyword_set(pts_limbe) and countj gt 0 then oplot,[xpol(wwj(wwjlh))],[ypol(wwj(wwjlh))],col=140,psym=4
endif
if wwirh(0) ne -1 then begin
loadct,3
oplot,[xpol(wwi(wwirh))],[ypol(wwi(wwirh))],col=140,psym=1
if keyword_set(pts_limbe) and countj gt 0 then oplot,[xpol(wwj(wwjrh))],[ypol(wwj(wwjrh))],col=140,psym=4
endif
if keyword_set(err_bar) then begin
loadct,1
if wwilh(0) ne -1 then begin
for tmp=0,n_elements(wwilh)-1l do begin
y2deg=loc_error(0,indgen(nerr+1) mod nerr,wwi(wwilh(tmp)))*cos(loc_error(1,indgen(nerr+1) mod nerr,wwi(wwilh(tmp)))*!dtor)*cos(loc_error(2,indgen(nerr+1) mod nerr,wwi(wwilh(tmp)))/24.*2*!pi)
x2deg=-loc_error(0,indgen(nerr+1) mod nerr,wwi(wwilh(tmp)))*cos(loc_error(1,indgen(nerr+1) mod nerr,wwi(wwilh(tmp)))*!dtor)*sin(loc_error(2,indgen(nerr+1) mod nerr,wwi(wwilh(tmp)))/24.*2*!pi)
w2deg=where(x2deg ne 0. and y2deg ne 0.)
if w2deg(0) gt 2 then begin
x2deg=[x2deg(w2deg),x2deg(w2deg(0))]
y2deg=[y2deg(w2deg),y2deg(w2deg(0))]
endif
oplot,x2deg,y2deg,col=150,thick=2
oplot,[xpol(wwi(wwilh(tmp)))],[ypol(wwi(wwilh(tmp)))],psym=1,col=140
endfor
endif
loadct,3
if wwirh(0) ne -1 then begin
for tmp=0,n_elements(wwirh)-1l do begin
y2deg=loc_error(0,indgen(nerr+1) mod nerr,wwi(wwirh(tmp)))*cos(loc_error(1,indgen(nerr+1) mod nerr,wwi(wwirh(tmp)))*!dtor)*cos(loc_error(2,indgen(nerr+1) mod nerr,wwi(wwirh(tmp)))/24.*2*!pi)
x2deg=-loc_error(0,indgen(nerr+1) mod nerr,wwi(wwirh(tmp)))*cos(loc_error(1,indgen(nerr+1) mod nerr,wwi(wwirh(tmp)))*!dtor)*sin(loc_error(2,indgen(nerr+1) mod nerr,wwi(wwirh(tmp)))/24.*2*!pi)
w2deg=where(x2deg ne 0. and y2deg ne 0.)
if w2deg(0) gt 2 then begin
x2deg=[x2deg(w2deg),x2deg(w2deg(0))]
y2deg=[y2deg(w2deg),y2deg(w2deg(0))]
endif
oplot,x2deg,y2deg,col=150,thick=2
oplot,[xpol(wwi(wwirh(tmp)))],[ypol(wwi(wwirh(tmp)))],psym=1,col=140
endfor
endif
endif
loadct,0
endif else tv,tab,!x.window(0),!y.window(0),xsize=!x.window(1)-!x.window(0),ysize=!y.window(1)-!y.window(0),/normal
endif
for i=0,23 do oplot,[0,rmax],([1,1]*float(i)/12)*!pi,/polar,line=1,col=100,thick=2
for i=latmin,90,10 do oplot,fltarr(101)+sqrt(1./(cos(i*!dtor)^2.+(sin(i*!dtor)/(1-0.098))^2.))*cos(i/!radeg),findgen(101)/100*2.*!pi,/polar,line=1,col=100,thick=2
if latmin eq 0. then begin
for i=latmin,latmin do oplot,fltarr(101)+sqrt(1./(cos(i*!dtor)^2.+(sin(i*!dtor)/(1-0.098))^2.))*cos(i/!radeg),findgen(101)/100*2.*!pi,/polar,col=100,thick=2
endif
loadct,2 & oplot,[0,-rmax*sin(tl(ww(count/2-1l))/24.*2*!pi)],[0,rmax*cos(tl(ww(count/2-1l))/24.*2*!pi)],thick=8,line=3,col=15
arrow,-0.95*rmax*sin(tl(ww(count/2-1l))/24.*2*!pi),0.95*rmax*cos(tl(ww(count/2-1l))/24.*2*!pi),$
-rmax*sin(tl(ww(count/2-1l))/24.*2*!pi),rmax*cos(tl(ww(count/2-1l))/24.*2*!pi),/solid,hthick=2,thick=2,col=15,/data
loadct,0
xyouts,0.69,0.482,'24:00',charsize=1.,/normal,align=0.5
loadct,3 & xyouts,0.69,0.003,'12:00',charsize=1.,col=150,/normal,align=0.5 & loadct,0
xyouts,0.52,0.25,'06:00',charsize=1.,/normal,align=0.5,orient=90
xyouts,0.862,0.25,'18:00',charsize=1.,/normal,align=0.5,orient=-90
if latmin eq 50. then begin
xyouts,0.815,0.105,'50!eo!n',charsize=0.8,col=0,/normal
xyouts,0.75,0.172,'70!eo!n',charsize=0.8,col=0,/normal
endif
if keyword_set(horizon) then begin
loadct,2
if ~ keyword_set(ntheta) then ntheta=180
if ~ keyword_set(nphi) then nphi=360
test=horizon_radio([mean(xyz(0,ww)),mean(xyz(1,ww)),mean(xyz(2,ww))],fmin,mfl_model='SPV',horizon_foot=horizon_foot,ntheta=ntheta,nphi=nphi)
wmfl=where(horizon_foot(1,*) gt 90.,complement=wmflcomp)
if lat(ww(count/2-1l)) gt 0. then wmfl=wmflcomp
x_horiz=-horizon_foot(0,wmfl)*cos((90.-horizon_foot(1,wmfl))*!dtor)*sin((horizon_foot(2,wmfl)+180.)*!dtor)
y_horiz=horizon_foot(0,wmfl)*cos((90.-horizon_foot(1,wmfl))*!dtor)*cos((horizon_foot(2,wmfl)+180.)*!dtor)
test_horiz=where(sqrt(x_horiz^2+y_horiz^2) le rmax)
if test_horiz(0) ne -1 then begin
x_horiz=x_horiz(test_horiz) & y_horiz=y_horiz(test_horiz)
endif
w_horiz = where(x_horiz ne 0 and y_horiz ne 0)
if w_horiz(0) ne -1 then begin
x_horiz = x_horiz(w_horiz)
y_horiz = y_horiz(w_horiz)
lon_ref = ((tl(ww(count/2-1l))/24.*2*!pi+0.5*!pi) mod (2.*!pi))
x1_horiz=x_horiz*cos(-lon_ref)-y_horiz*sin(-lon_ref)
y1_horiz=x_horiz*sin(-lon_ref)+y_horiz*cos(-lon_ref)
x0_horiz=double(max(x1_horiz)+min(x1_horiz))/2.
y0_horiz=double(max(y1_horiz)+min(y1_horiz))/2.
test_horiz=atan((y1_horiz-y0_horiz),(x1_horiz-x0_horiz))
wxy=where(test_horiz lt 0) & if wxy(0) ne -1 then test_horiz(wxy)+=2.*!pi
s_horiz=sort(test_horiz)
oplot,x_horiz(s_horiz),y_horiz(s_horiz),thick=6,col=10
endif
if fmax ne fmin then begin
test=horizon_radio([mean(xyz(0,ww)),mean(xyz(1,ww)),mean(xyz(2,ww))],fmax,mfl_model='SPV',horizon_foot=horizon_foot,ntheta=ntheta,nphi=nphi)
wmfl=where(horizon_foot(1,*) gt 90.,complement=wmflcomp)
if lat(ww(count/2-1l)) gt 0. then wmfl=wmflcomp
x_horiz=-horizon_foot(0,wmfl)*cos((90.-horizon_foot(1,wmfl))*!dtor)*sin((horizon_foot(2,wmfl)+180.)*!dtor)
y_horiz=horizon_foot(0,wmfl)*cos((90.-horizon_foot(1,wmfl))*!dtor)*cos((horizon_foot(2,wmfl)+180.)*!dtor)
test_horiz=where(sqrt(x_horiz^2+y_horiz^2) le rmax)
if test_horiz(0) ne -1 then begin
x_horiz=x_horiz(test_horiz) & y_horiz=y_horiz(test_horiz)
endif
w_horiz = where(x_horiz ne 0 and y_horiz ne 0)
if w_horiz(0) ne -1 then begin
x_horiz = x_horiz(w_horiz)
y_horiz = y_horiz(w_horiz)
lon_ref = ((tl(ww(count/2-1l))/24.*2*!pi+0.5*!pi) mod (2.*!pi))
x1_horiz=x_horiz*cos(-lon_ref)-y_horiz*sin(-lon_ref)
y1_horiz=x_horiz*sin(-lon_ref)+y_horiz*cos(-lon_ref)
x0_horiz=double(max(x1_horiz)+min(x1_horiz))/2.
y0_horiz=double(max(y1_horiz)+min(y1_horiz))/2.
test_horiz=atan((y1_horiz-y0_horiz),(x1_horiz-x0_horiz))
wxy=where(test_horiz lt 0) & if wxy(0) ne -1 then test_horiz(wxy)+=2.*!pi
s_horiz=sort(test_horiz)
oplot,x_horiz(s_horiz),y_horiz(s_horiz),thick=6,col=20
endif
endif
loadct,0
endif
if ~keyword_set(position) then begin
posxmin=0.38 & posymin=0
posxmax=0.4 & posymax=0.25
bar_titre=''
if keyword_set(dB) then bar_tick = [' ','!9S!3 dB!Inorm!N',' '] else bar_tick = [' ','dB',' ']
plot, [0,1], [0,1], /nodata, /noerase, xst=1, yst=12, $
pos = [posxmin,posymin,posxmax,posymax], $
/norm, xticklen=0, xticks=2, xtickv=[0,0.5,1], $
xtickname=bar_tick, xtitle=bar_titre, charsize=0.8
axis, /yaxis, yra=[smin1,smax1], ticklen=-0.02, /ysty, charsize=0.8
bcb = 255-bytscl(replicate(1,10)#bindgen(256))
tv,bcb,!x.window(0),!y.window(0),xsize = !x.window(1)-!x.window(0), $
ysize = !y.window(1) - !y.window(0), /normal
oplot, [1,0,0,1], [1,1,0,0]
posxmin=0.9 & posymin=0
posxmax=0.92 & posymax=0.25
bar_titre=''
if keyword_set(dB) then bar_tick = [' ','!9S!3 dB!Inorm!N',' '] else bar_tick = [' ','dB',' ']
plot, [0,1], [0,1], /nodata, /noerase, xst=1, yst=12, $
pos = [posxmin,posymin,posxmax,posymax], $
/norm, xticklen=0, xticks=2, xtickv=[0,0.5,1], $
xtickname=bar_tick, xtitle=bar_titre, charsize=0.8
axis, /yaxis, yra=[smin2,smax2], ticklen=-0.02, /ysty, charsize=0.8
bcb = 255-bytscl(replicate(1,10)#bindgen(256))
tv,bcb,!x.window(0),!y.window(0),xsize = !x.window(1)-!x.window(0), $
ysize = !y.window(1) - !y.window(0), /normal
oplot, [1,0,0,1], [1,1,0,0]
endif
if keyword_set(position) then begin
loadct,3
posxmin=0.9 & posymin=0.1
posxmax=0.92 & posymax=0.25
plot, [0,1], [0,1], /nodata, /noerase, xst=1, yst=12, $
pos = [posxmin,posymin,posxmax,posymax], $
/norm, xticklen=0, xticks=2, xtickv=[0,0.5,1], $
xtickname=[' ','f (kHz)',' '], xtitle='RH', charsize=0.8
axis, /yaxis, yra=[fmin,fmax], ticklen=-0.02, /ysty, charsize=0.8
bcb = bytscl(replicate(1,10)#(bindgen(256)/255.*150+70))
tv,bcb,!x.window(0),!y.window(0),xsize = !x.window(1)-!x.window(0), $
ysize = !y.window(1) - !y.window(0), /normal
oplot, [1,0,0,1], [1,1,0,0]
loadct,1
posxmin=0.97 & posymin=0.1
posxmax=0.99 & posymax=0.25
plot, [0,1], [0,1], /nodata, /noerase, xst=1, yst=12, $
pos = [posxmin,posymin,posxmax,posymax], $
/norm, xticklen=0, xticks=2, xtickv=[0,0.5,1], $
xtickname=[' ','f (kHz)',' '], xtitle='LH', charsize=0.8
axis, /yaxis, yra=[fmin,fmax], ticklen=-0.02, /ysty, charsize=0.8
bcb = bytscl(replicate(1,10)#(bindgen(256)/255.*150+70))
tv,bcb,!x.window(0),!y.window(0),xsize = !x.window(1)-!x.window(0), $
ysize = !y.window(1) - !y.window(0), /normal
oplot, [1,0,0,1], [1,1,0,0]
endif
xyouts,0.38,0.48,'Ephemeris:',charsize=1.1,col=0,/normal
xyouts,0.38,0.478,'__________',charsize=1.,col=0,/normal
aaaa=strmid(string(t97_aj(n2(ww(count/2-1l)).t97)),7,4)
jjj=strmid(string(t97_aj(n2(ww(count/2-1l)).t97)),11,3)
hh=strtrim(string(long(h1+k*dtmin/60.)),1) & if (h1+k*dtmin/60.) lt 10. then hh='0'+hh
mm=strtrim(string(long(h1*60.+k*dtmin-hh*60.)),1) & if float(mm) lt 10. then mm='0'+mm
xyouts,0.375,0.44,'Day',charsize=1.1,col=0,/normal
xyouts,0.425,0.44,': '+aaaa+'-'+jjj+' ',charsize=1.1,col=0,/normal
xyouts,0.375,0.41,'Time',charsize=1.1,col=0,/normal
xyouts,0.425,0.41,': '+hh+':'+mm,charsize=1.1,col=0,/normal
hhtl=strmid(strtrim(string(tl(ww(count/2-1l))),1),0,2) & if hhtl lt 10. then hhtl='0'+strmid(hhtl,0,1)
mmtl=strtrim(string(long(60*(tl(ww(count/2-1l))-long(tl(ww(count/2-1l)))))),1) & if mmtl lt 10. then mmtl='0'+strmid(mmtl,0,1)
xyouts,0.375,0.37,'r!iS/C!n (R!iS!n)',charsize=1.1,col=0,/normal
xyouts,0.435,0.37,'= '+strmid(strtrim(string(rs(ww(count/2-1l))),1),0,4),charsize=1.1,col=0,/normal
xyouts,0.375,0.34,'!9l!3!iS/C!n (!eo!n)',charsize=1.1,col=0,/normal
xyouts,0.435,0.34,'= '+strmid(strtrim(string(lat(ww(count/2-1l))),1),0,5),charsize=1.1,col=0,/normal
loadct,2 & xyouts,0.375,0.31,'TL!iS/C!n',charsize=1.,col=15,/normal
xyouts,0.435,0.31,'= '+hhtl+':'+mmtl,charsize=1.1,col=15,/normal & loadct,0
xyouts,0.9,0.48,'Data selection:',charsize=1.1,/normal
xyouts,0.9,0.478,'_____________',charsize=1.,/normal
xyouts,0.9,0.44,'!9d!3t (min)',charsize=1.,/normal
xyouts,0.97,0.44,'= '+strmid(dtmin,7,2),charsize=1.,/normal
if keyword_set(horizon) then begin
loadct,2
xyouts,0.9,0.41,'f!imin!n (kHz)',charsize=1.,col=10,/normal
xyouts,0.97,0.41,'= '+strtrim(strmid(long(fmin),8,4),1),charsize=1.,col=10,/normal
xyouts,0.9,0.38,'f!imax!n (kHz)',charsize=1.,col=20,/normal
xyouts,0.97,0.38,'= '+strtrim(strmid(long(fmax),8,4),1),charsize=1.,col=20,/normal
loadct,0
endif else begin
xyouts,0.9,0.41,'f!imin!n (kHz)',charsize=1.,/normal
xyouts,0.97,0.41,'= '+strtrim(strmid(long(fmin),8,4),1),charsize=1.,col=10,/normal
xyouts,0.9,0.38,'f!imax!n (kHz)',charsize=1.,/normal
xyouts,0.97,0.38,'= '+strtrim(strmid(long(fmax),8,4),1),charsize=1.,col=20,/normal
endelse
xyouts,0.9,0.35,'V',charsize=1.,/normal
xyouts,0.97,0.35,'= ['+strmid(strtrim(string(vmin),1),0,4)+','+strmid(strtrim(string(vmax),1),0,4)+']',charsize=1.,/normal
xyouts,0.9,0.32,'SN!iX,Z!n (dB)',charsize=1.,/normal
!p.font=-1 & xyouts,0.98,0.32,'!9'+string("142b)+'!3',charsize=1.,/normal & !p.font=0
xyouts,0.995,0.32,strmid(strtrim(string(snrmin),1),0,2),charsize=1.,/normal
xyouts,0.9,0.29,'!9b!3 (!eo!n)',charsize=1.,/normal
!p.font=-1 & xyouts,0.98,0.29,'!9'+string("142b)+'!3',charsize=1.,/normal & !p.font=0
xyouts,0.995,0.29,strmid(strtrim(string(betamin),1),0,2),charsize=1.,/normal
endif else begin
if ~keyword_set(xmin) then xmin=-1.
if ~keyword_set(xmax) then xmax=+1.
ymin=xmin & ymax=xmax
if keyword_set(spdyn_flux) or keyword_set(spdyn_polar) then begin
!p.multi(0)=2
!p.multi(2)=2
endif
plot,[0],[0],xr=[xmin,xmax],yr=[ymin,ymax],xstyle=13,ystyle=13,$
xtit='R!iS!n',ytit='R!iS!n',tit='',charsize=1.1,ticklen=0.015,/iso,pos=[0,0,0.5,0.5]
xyouts,0.185,0.51,'View from Cassini (FOV = '+view+'!eo!n)',charsize=1.3,col=0,/normal,align=0.5
xyouts,0.19,0.25,'No data',charsize=1.3,col=0,/normal,align=0.5
if keyword_set(spdyn_flux) or keyword_set(spdyn_polar) then !p.multi(0)=1
rmax=sqrt(1./(cos(latmin*!dtor)^2.+(sin(latmin*!dtor)/(1-0.098))^2.))*cos(latmin*!dtor)
xmin=-rmax & xmax=rmax
ymin=-rmax & ymax=rmax
plot,[xmin,xmax],[ymin,ymax],xr=[xmin,xmax],yr=[ymin,ymax],xstyle=13,ystyle=13,$
xtit='',ytit='',tit='',charsize=1.,ticklen=0.2,/nodata,/iso,pos=[0.525,0.025,0.975,0.475]
xyouts,0.69,0.51,'Magnetic polar projection',charsize=1.3,col=0,/normal,align=0.5
for i=0,23,1 do oplot,[0,rmax],([1,1]*float(i)/12)*!pi,/polar,line=1,col=100,thick=2
for i=latmin,90,10 do oplot,fltarr(101)+sqrt(1./(cos(i*!dtor)^2.+(sin(i*!dtor)/(1-0.098))^2.))*cos(i/!radeg),findgen(101)/100*2.*!pi,/polar,line=1,col=100,thick=2
if latmin eq 0. then begin
for i=latmin,latmin do oplot,fltarr(101)+sqrt(1./(cos(i*!dtor)^2.+(sin(i*!dtor)/(1-0.098))^2.))*cos(i/!radeg),findgen(101)/100*2.*!pi,/polar,col=100,thick=2
endif
xyouts,0.69,0.482,'24:00',charsize=1.,col=150,/normal,align=0.5
loadct,3 & xyouts,0.69,0.003,'12:00',charsize=1.,col=150,/normal,align=0.5 & loadct,0
xyouts,0.52,0.25,'06:00',charsize=1.,col=150,/normal,align=0.5,orient=90
xyouts,0.862,0.25,'18:00',charsize=1.,col=150,/normal,align=0.5,orient=-90
if latmin eq 50. then begin
xyouts,0.815,0.105,'50!eo!n',charsize=0.8,col=0,/normal
xyouts,0.75,0.172,'70!eo!n',charsize=0.8,col=0,/normal
endif
xyouts,0.38,0.48,'Ephemeris:',charsize=1.1,col=0,/normal
xyouts,0.38,0.478,'__________',charsize=1.,col=0,/normal
aaaa=strmid(string(jour_aj),7,4)
jjj=strmid(string(jour_aj),11,3)
hh=strtrim(string(long(h1+k*dtmin/60.)),1) & if (h1+k*dtmin/60.) lt 10. then hh='0'+hh
mm=strtrim(string(long(h1*60.+k*dtmin-hh*60.)),1) & if (h1*60.+k*dtmin-hh*60.) lt 10. then mm='0'+mm
xyouts,0.375,0.44,'Day',charsize=1.1,col=0,/normal
xyouts,0.425,0.44,': '+aaaa+'-'+jjj+' ',charsize=1.1,col=0,/normal
xyouts,0.375,0.41,'Time',charsize=1.1,col=0,/normal
xyouts,0.425,0.41,': '+hh+':'+mm,charsize=1.1,col=0,/normal
xyouts,0.9,0.48,'Data selection:',charsize=1.1,col=150,/normal
xyouts,0.9,0.478,'_____________',charsize=1.,col=150,/normal
xyouts,0.9,0.44,'No data',charsize=1.,col=150,/normal
endelse
endfor
endfor
device,/close & set_plot,'x'
spawn,'ps2pdf '+titre+'.ps '+titre+'.pdf'
spawn,'rm -f '+titre+'.ps'
end