;************************************************************************ ;* ;* mono_esr.pro ;* for ESR ;* monoesr.date ;* Version 1.0 Oct 12, 1999 Copied ;* 1.1 Oct 30, 1999 990308 ;* 1.2 Oct 31, 1999 980921 ;* 1.2 Nov 1, 1999 Last edition ;* 1.3 Nov 22, 1999 lat and lon for E and wind ;* 1.4 July 24, 2000 long pulse ;* 1.5 Oct 28, 2000 aclp (AC or LP) ;* 1.6 Dec 1, 2000 add latmono_all,lonmono_all for ti files ;*--------------------------------------------------------------------------------------- ;* mono.pro for calculation of ion velocities ;* monostatic method, especially CP2 EISCAT data ;* ;* set mono.date ;* ;*--------------------------------------------------------------------------------------- ;* Version 1.0 Feb 16, 1995 Created by S. Nozawa ;* 1.0 Feb 17, 1995 Last edition ;* ;* 2.0 Nov 25, 1996 modified ;* 2.0 Dec 1, 1996 Last edition ;* 2.1 Dec 3, 1996 twice ;* 2.2 Dec 3, 1996 data at alt le 0 descarded ;* 2.3 Dec 4, 1996 n1_treatment for 890828 ;* 2.4 Dec 6, 1996 long pulse ;* 3.0 Apr 7, 1997 average_ne, Ne included ;* 4.0 July 15, 1999 GUISDAP ;* 4.1 July 17, 1999 Calculate errors as well ;* Northward, eastward and vertically downward directions are counted positive ;************************************************************************ pro print_proname,proname date='Version 1.5 on October 28, 2000' & s='-----' & s2=' ' print,s,s2,proname,s2,date,s2,s end function cal_verr,v10,v20,v30,a10,a20,a30 a1=a10*1d+0 & a2=a20*1d+0 & a3=a30*1d+0 v1=v10*1d+0 & v2=v20*1d+0 & v3=v30*1d+0 vv=(v1^2)*(a1^2) + (v2^2)*(a2^2) + (v3^2)*(a3^2) if vv lt 0 then stop v=sqrt(vv) & v=float(v) return,v end ;* mono_cal2.pro ;* 1.0 July 17, 1999 error estimation pro mono_cal2,vmono,vmonoerr,qmono_all,Emono,Amono,Vxyz,Vxyzerr,qq,fap=ion,$ count=icount if not keyword_set(ion) then ion = 0 Az1p=Amono(0) & Az2p=Amono(1) & Az3p=Amono(2) El1p=Emono(0) & El2p=Emono(1) & El3p=Emono(2) Vr1=vmono(0) & Vr2=vmono(1) & Vr3=vmono(2) Vr1err=vmonoerr(0) & Vr2err=vmonoerr(1) & Vr3err=vmonoerr(2) qmono=qmono_all if n_elements(Amono) gt 3 then begin Az4p=Amono(3) & El4p=Emono(3) & Vr4=vmono(3) & Vr4err=vmonoerr(3) endif if ion eq 1 then begin if keyword_set(icount) then if icount eq 1 then $ print,'+++ use FAP instead of Vertical position +++' Az1p=Amono(1) & Az2p=Amono(2) & Az3p=Amono(3) & Az4p=Amono(0) El1p=Emono(1) & El2p=Emono(2) & El3p=Emono(3) & El4p=Emono(0) Vr1=vmono(1) & Vr2=vmono(2) & Vr3=vmono(3) & Vr4=vmono(0) Vr1err=vmonoerr(1) & Vr2err=vmonoerr(2) Vr3err=vmonoerr(3) & Vr4err=vmonoerr(0) endif fac=!PI/180.0 ;print,format='(A16,3F7.1)','Az1, Az2, Az3 = ',Az1p,Az2p,Az3p ;print,format='(A16,3F7.1)','El1, El2, El3 = ',El1p,El2p,El3p Az1p=Az1p*fac & El1p=El1p*fac Az2p=Az2p*fac & El2p=El2p*fac Az3p=Az3p*fac & El3p=El3p*fac ; from A. Brekke, E-feilds, neutral winds and currents derived from ; Chatanika, Radar Probing of the auroral plasma, pp. 285, UNIVERSITETSFORLAGET, ; edited by A. Brekke, 1977 ; originally, (south, east, upward) K1=-cos(El1p)*cos(Az1p) & K2=cos(El1p)*sin(Az1p) & K3=sin(El1p) K4=-cos(El2p)*cos(Az2p) & K5=cos(El2p)*sin(Az2p) & K6=sin(El2p) K7=-cos(El3p)*cos(Az3p) & K8=cos(El3p)*sin(Az3p) & K9=sin(El3p) DELTA=K1*(K5*K9 - K6*K8) - K4*(K2*K9 - K3*K8) + K7*(K2*K6 - K3*K5) ;--- Calculate velocities Vx= Vr1*(K5*K9 - K6*K8) - Vr2*(K2*K9 - K3*K8) + Vr3*(K2*K6 - K3*K5) Vy=-Vr1*(K4*K9 - K6*K7) + Vr2*(K1*K9 - K3*K7) - Vr3*(K1*K6 - K3*K4) Vz= Vr1*(K4*K8 - K5*K7) - Vr2*(K1*K8 - K2*K7) + Vr3*(K1*K5 - K2*K4) Vx=Vx/delta & Vy=Vy/delta & Vz=Vz/delta ;print_parameter,k1,k2,k3,k4,k5,k6,k7,k8,k9,delta qq=total(abs(qmono(0:2))) ; originally, (south, east, upward) if qq ne 0 then Vxyz=[0.0,0.0,0.0] else Vxyz=[-Vx,Vy,-Vz] ;print,format='(4I2,3F8.1)',qq,qmono(0:2),Vxyz calculate_verr,vr1err,vr2err,vr3err,k1,k2,k3,k4,k5,k6,k7,k8,k9,delta,Vxyzerr end pro calculate_verr,vr1err,vr2err,vr3err,k1,k2,k3,k4,k5,k6,k7,k8,k9,delta,Vxyzerr d=delta a11=(K5*K9 - K6*K8)/d & a12=(K2*K9 - K3*K8)/d & a13=(K2*K6 - K3*K5)/d a21=(K4*K9 - K6*K7)/d & a22=(K1*K9 - K3*K7)/d & a23=(K1*K6 - K3*K4)/d a31=(K4*K8 - K5*K7)/d & a32=(K1*K8 - K2*K7)/d & a33=(K1*K5 - K2*K4)/d Verr2=[Vr1err,Vr2err,Vr3err] tmp=where(Verr2 lt 0,nerr) if nerr ne 0 then begin i=-32768 & Vxyzerr=[i,i,i] endif else begin vxerr=cal_verr(Vr1err,Vr2err,Vr3err,a11,a12,a13) vyerr=cal_verr(Vr1err,Vr2err,Vr3err,a21,a22,a23) vzerr=cal_verr(Vr1err,Vr2err,Vr3err,a31,a32,a33) Vxyzerr=[vxerr,vyerr,vzerr] endelse end ;970716 for p7 (Radio Science) pro print_parameter,k1,k2,k3,k4,k5,k6,k7,k8,k9,delta a11=(K5*K9 - K6*K8)/delta a12=-(K2*K9 - K3*K8)/delta a13=(K2*K6 - K3*K5)/delta a21=-(K4*K9 - K6*K7)/delta a22=(K1*K9 - K3*K7)/delta a23=-(K1*K6 - K3*K4)/delta a31=(K4*K8 - K5*K7)/delta a32=-(K1*K8 - K2*K7)/delta a33=(K1*K5 - K2*K4)/delta print,format='(3F10.3)',a11,a12,a13 print,format='(3F10.3)',a21,a22,a23 print,format='(3F10.3)',a31,a32,a33 end ; until here mono_cal2.pro function where_el,Els,El_p1,n1 dEl=0.3 & if El_p1 gt 70 then dEl=6.0 else dEl=1.1 print,format='(A12,2F5.1)','El and dE = ',El_p1,dEl E1=El_p1+dEl & E2=El_p1-dEl & tmp1=where(Els le E1 and Els ge E2,n1) return,tmp1 end pro openr_print,i,filename1 openr,i,filename1 & print,'----- Open (r) = ',filename1 end pro openw_print,i,filename1 openw,i,filename1 & print,'----- Open (w) = ',filename1 end pro mono_readdata,filename1,icount,isize1,isize2, $ dates,time1s,time2s,alts,lats,lons,quos,Vlsr,Azs,Els,Nearr,Vlsrerr,Neerr moji='string' & moji1='string' & moji2='string' & moji3='string' icount=0l & isize1=0l & isize2=0l openr_print,1,filename1 & readf,1,moji1 & readf,1,moji2 & readf,1,moji3 readu,1,icount,isize1,isize2 print,format='(I8,2I3)',icount,isize1,isize2 idata_arr=intarr(isize1,icount) & adata_arr=fltarr(isize2,icount) readu,1,idata_arr,adata_arr & close,1 dates=reform(idata_arr(0,*)) & time1s=reform(idata_arr(1,*)) time2s=reform(idata_arr(2,*)) & alts=reform(idata_arr(3,*)) lats=reform(idata_arr(4,*)) & lons=reform(idata_arr(5,*)) quos=reform(idata_arr(6,*)) & Vlsr=reform(idata_arr(8,*)) Azs=reform(adata_arr(0,*)) & Els=reform(adata_arr(1,*)) Nearr=reform(idata_arr(7,*)) ; Ne (10^(Ne*0.001)) Tiarr =reform(idata_arr(9,*)) ; Ti Tearr =reform(idata_arr(10,*)) ; Te Neerr =reform(idata_arr(11,*)) ; Ne (10^(Ne*0.001)) Vlsrerr =reform(idata_arr(12,*)) ; Vlsr Tierr =reform(idata_arr(13,*)) ; Ti Teerr =reform(idata_arr(14,*)) ; Te end pro set_az_el,cp,arr0_Az,arr0_El ; CP-2-d if cp eq 'cp2d' then begin arr0_Az=[180.0,166.5,133.3,182.6] & arr0_El=[89.9,62.9,60.4,77.5] endif ; CP-2-E if cp eq 'cp2e' then begin arr0_Az=[180.0,166.5,133.3,182.6] & arr0_El=[89.9,62.9,60.4,77.5] endif end pro set_az_el_970106,cp,arr0_Az,arr0_El arr0_Az=[177.5,166.5,133.3,182.6] & arr0_El=[90.0,57.2,54.5,77.5] end pro set_min_max,Az,dAz,Az1,Az2 Az1=Az-dAz & Az2=Az+dAz end pro set_arr_poss,k,tmp1,ipos1,ipos2,ipos3,ipos4 Case k of 0:ipos1=tmp1 1:ipos2=tmp1 2:ipos3=tmp1 3:ipos4=tmp1 endcase end pro where_az_el,arr0_Az,arr0_El,Azs,Els,alts,ipos1,ipos2,ipos3,ipos4 n=n_elements(arr0_Az)-1 for i=0,n do begin ; loop of antenna position El=arr0_El(i) & Az=arr0_Az(i) & dAz=15.0 if El gt 70 then dEl=6.0 else dEl=1.1 set_min_max,Az,dAz,Az1,Az2 & set_min_max,El,dEl,El1,El2 tmp1=where(Els ge El1 and Els le El2 and $ Azs ge Az1 and Azs le Az2 and alts gt 10) set_arr_poss,i,tmp1,ipos1,ipos2,ipos3,ipos4 endfor n1=n_elements(ipos1) & n2=n_elements(ipos2) n3=n_elements(ipos3) & n4=n_elements(ipos4) print,format='(6I7)',n1,n2,n3,n4,n1+n2+n3+n4,n_elements(Els) end pro put_date_order,year,dates days=[31,28,31,30,31,30,31,31,30,31,30,31] if year mod 4 eq 0 then days(1)=29 month=dates(0)/100 -1 dt=dates-dates(0) tmp=where(dt gt 12,n) if n gt 0 then dates(tmp)=dates(0)/100*100+days(month)+dates(tmp) mod 100 end ; last edition made on 990716: Verr and Neerr are added. pro put_arr_data,ipos1,t1,t2,alts,Azs,Els,Vlsr,quos,Nearr0,Verr,Neerr,$ t11,t21,alts1,Az1,El1,Vlsr1,quos1,Nearr1,Verr1,Neerr1 t11=t1(ipos1) & t21=t2(ipos1) & alts1=alts(ipos1) & Az1=Azs(ipos1) El1=Els(ipos1) & Vlsr1=Vlsr(ipos1) & quos1=quos(ipos1) Nearr1=Nearr0(ipos1) & Verr1=Verr(ipos1) & Neerr1=Neerr(ipos1) end pro put_data_arr,n1max,i,alt1,alt1arr n1=n_elements(alt1) ntemp=n1max put_n1_max,n1,ntemp,iflag ; iflag -1=less, 0=same, 1=larger, 2=twice Case iflag of 0:alt1arr(*,i)=alt1(*) -1:alt1arr(0:n1-1,i)=alt1(*) 1:stop ; need to make further efforts 2:alt1arr(*,i)=alt1(n1max:n1-1) endcase if iflag eq -1 then alt1arr(n1:n1max-1,i)=-9 ; no data end pro make_arrs,n1max,ic,alt1,t1,t2,A1,E1,q1,V1,Ne1,V1err,Ne1err,$ alt1arr,t1arr,t2arr,$ A1arr,E1arr ,q1arr,vlsr1arr,Ne1arr,vlsr1err,Ne1errarr,iflag n1=n_elements(alt1) & n1temp=n1max & put_n1_max,n1,n1temp,iflag if iflag eq -1 then begin print,format='(A8,I2,A8,I5)','iflag = ',iflag,' at ic =',ic & return endif put_data_arr,n1max,ic,alt1,alt1arr & put_data_arr,n1max,ic,t1,t1arr put_data_arr,n1max,ic,t2,t2arr & put_data_arr,n1max,ic,A1,A1arr put_data_arr,n1max,ic,E1,E1arr & put_data_arr,n1max,ic,q1,q1arr put_data_arr,n1max,ic,v1,vlsr1arr put_data_arr,n1max,ic,Ne1,Ne1arr put_data_arr,n1max,ic,v1err,vlsr1err put_data_arr,n1max,ic,Ne1err,Ne1errarr end pro set_arr_init,n1,isize,alt1arr,t1arr,t2arr,A1arr,E1arr,q1arr,$ vlsr1arr,Ne1,vlsr1err,Neerr1 common makebi,ibi_ver alt1arr=intarr(n1,isize) & t1arr=intarr(n1,isize) t2arr=intarr(n1,isize) & A1arr=fltarr(n1,isize) & E1arr=fltarr(n1,isize) q1arr=intarr(n1,isize) & vlsr1arr=intarr(n1,isize) & Ne1=intarr(n1,isize) vlsr1err=intarr(n1,isize) & Neerr1=intarr(n1,isize) if ibi_ver eq 3 then alt1arr=fltarr(n1,isize) end pro remake_arrs,ic1loop,alt1arr,t1arr,t2arr,A1arr,E1arr,q1arr,$ vlsr1arr,Ne1,vlsr1err,Ne1err i=ic1loop-1 & alt1arr=alt1arr(*,0:i) & t1arr=t1arr(*,0:i) t2arr=t2arr(*,0:i) & a1arr=a1arr(*,0:i) & E1arr=E1arr(*,0:i) q1arr=q1arr(*,0:i) & vlsr1arr=vlsr1arr(*,0:i) Ne1=Ne1(*,0:i) vlsr1err=vlsr1err(*,0:i) & Ne1err=Ne1err(*,0:i) end pro get_dev_n,alt1arr,is1 n1=n_elements(alt1arr(*,0))-1 & n2=n_elements(alt1arr(0,*))-1 is1=0 for i=0,n2 do begin alt=reform(alt1arr(*,i)) & a1=[alt,1000] & a2=[alt(0),alt] da12=a1-a2 & tmp=where(da12 lt 0,n) & if n ne 1 then stop if tmp(0) gt is1 then is1=tmp(0) endfor end pro get_dev_data,altorg,alt1arr,alt1arrs,alt1arrl n1=n_elements(altorg(*,0))-1 & n2=n_elements(altorg(0,*))-1 ik1=n_elements(alt1arrs(*,0)) & ik2=n_elements(alt1arrl(0,*)) & is1=0 for i=0,n2 do begin alt=reform(alt1arr(*,i)) alt0=reform(altorg(*,i)) & a1=[alt0,1000] & a2=[alt0(0),alt0] aa=where(alt0 eq -9,bb) ; if bb gt 0 then it means short of data sets da12=a1-a2 & tmp=where(da12 lt 0,n) & if n ne 1 then stop t=tmp(0) if bb gt 0 then begin alt1arrs(*,i)=alt(0:ik1-1) & alt1arrl(*,i)=alt(ik1:n1) print,i,'error data: The number of alt=0 is more than 3!' goto, J1 endif if t eq ik1 then begin alt1arrs(*,i)=alt(0:t-1) & alt1arrl(*,i)=alt(t:n1) endif else begin alt1arrs(0:t-1,i)=alt(0:t-1) & alt1arrl(0:n1-t,i)=alt(t:n1) endelse if tmp(0) gt is1 then is1=tmp(0) J1: endfor ; end of loop i end pro dev_data_sub,is1,alt,alt1arr,alt1arrs,alt1arrl n1=n_elements(alt1arr(*,0)) & n2=n_elements(alt1arr(0,*)) alt1arrs=intarr(is1,n2) & alt1arrl=intarr(n1-is1,n2) get_dev_data,alt,alt1arr,alt1arrs,alt1arrl end pro put_n1_max,n1,n1max,iflag iflag=0 ; same if n1 eq n1max then return n0=n1max*2 if n0 eq n1 then begin iflag=2 ; twice return endif if n1 lt n1max then begin iflag=-1 ; less return endif if n1 gt n1max then n1max=n1 iflag=1 ; replace end pro dev_data_subf,is1,alt,alt1arr,alt1arrs,alt1arrl n1=n_elements(alt1arr(*,0)) & n2=n_elements(alt1arr(0,*)) alt1arrs=fltarr(is1,n2) & alt1arrl=fltarr(n1-is1,n2) get_dev_data,alt,alt1arr,alt1arrs,alt1arrl end pro dev_data,is1,alt1arr,t1arr,t2arr,A1arr,E1arr,q1arr,vlsr1arr,Ne1arr,v1err,Neerr, $ alt1arrs,t1arrs,t2arrs,A1arrs,E1arrs,q1arrs,vlsr1arrs, Ne1arrs,v1errs,Neerrs,$ alt1arrl,t1arrl,t2arrl,A1arrl,E1arrl,q1arrl,vlsr1arrl, Ne1arrl,v1errl,Neerrl common makebi,ibi_ver if ibi_ver eq 3 then dev_data_subf,is1,alt1arr,alt1arr,alt1arrs,alt1arrl else $ dev_data_sub,is1,alt1arr,alt1arr,alt1arrs,alt1arrl dev_data_sub,is1,alt1arr,t1arr,t1arrs,t1arrl dev_data_sub,is1,alt1arr,t2arr,t2arrs,t2arrl dev_data_subf,is1,alt1arr,A1arr,A1arrs,A1arrl dev_data_subf,is1,alt1arr,E1arr,E1arrs,E1arrl dev_data_sub,is1,alt1arr,q1arr,q1arrs,q1arrl dev_data_sub,is1,alt1arr,vlsr1arr,vlsr1arrs,vlsr1arrl dev_data_sub,is1,alt1arr,Ne1arr,Ne1arrs,Ne1arrl dev_data_sub,is1,alt1arr,v1err,v1errs,v1errl dev_data_sub,is1,alt1arr,Neerr,Neerrs,Neerrl end pro n1_treatment,date,ic1,ic1loop,ipos1,icount,t1,t3,n1,tmp if date eq 890828 and n1 eq 101 then begin ttt=t1(ipos1(icount)) if ic1 eq 1 then print,format='(A24,3I8)',$ '890828 double data sets',ic1loop,t3,ttt icount=icount+20 & tmp=tmp(20:100) & n1=n_elements(tmp) endif end ; edited on July 16, 1999: Vlsrerr and Neerr are added pro put_make_data_arrs,date,ipos1,t1,t2,alts,Azs,Els,Vlsr,quos,Nearr,Verr,Neerr,$ alt1arrs,t1arrs,t2arrs,A1arrs,E1arrs,q1arrs,vlsr1arrs,Nearrs,Verrs,Neerrs, $ alt1arrl,t1arrl,t2arrl,A1arrl,E1arrl,q1arrl,vlsr1arrl,Nearrl,Verrl,Neerrl put_arr_data,ipos1,t1,t2,alts,Azs,Els,Vlsr,quos,Nearr,Verr,Neerr,$ t11,t21,alts1,Az1,El1,Vlsr1,quos1,Nearr1,Verr1,Neerr1 isize=1440 & n1max=0 if date eq 990701 then isize=2200 for ic1=1,2 do begin if ic1 eq 2 then $ set_arr_init,n1max,isize,alt1arr,t1arr,t2arr,A1arr,E1arr,q1arr,$ vlsr1arr,Ne1arr,vlsr1err,Ne1arrerr icount=0 & ic1loop=0 while icount lt n_elements(ipos1) do begin t3=t2(ipos1(icount)) tmp=where(t21 eq t3,n1) n1_treatment,date,ic1,ic1loop,ipos1,icount,t1,t3,n1,tmp put_arr_data,tmp,t11,t21,alts1,Az1,El1,Vlsr1,quos1,Nearr1,Verr1,Neerr1,$ t111,t211,alt1,A1,E1,V1,q1,Ne1,V1err,Ne1err Case ic1 of 1:put_n1_max,n1,n1max,iflag 2:make_arrs,n1max,ic1loop,alt1,t111,t211,A1,E1,q1,V1,Ne1,V1err,Ne1err,$ alt1arr,t1arr,t2arr,A1arr,E1arr,q1arr,vlsr1arr,Ne1arr,vlsr1err,Ne1arrerr,iflag endcase ;print,format='(4I6)',ipos1(icount),n1,n1max,ic1loop if iflag eq 2 then print,$ format='(A7,2I5,2F6.1)','Twice: ',t111(0),t211(0),A1(0),E1(0) icount=icount+n1 if iflag ne -1 then ic1loop=ic1loop+1 ; later make efforts needed. endwhile endfor ; ic1 remake_arrs,ic1loop,alt1arr,t1arr,t2arr,A1arr,E1arr,q1arr,$ vlsr1arr,Ne1arr,vlsr1err,Ne1arrerr get_dev_n,alt1arr,is1 ; reduce size dev_data,is1,alt1arr,t1arr,t2arr,A1arr,E1arr,q1arr,vlsr1arr,Ne1arr,vlsr1err,Ne1arrerr, $ alt1arrs,t1arrs,t2arrs,A1arrs,E1arrs,q1arrs,vlsr1arrs, Nearrs, Verrs,Neerrs, $ alt1arrl,t1arrl,t2arrl,A1arrl,E1arrl,q1arrl,vlsr1arrl, Nearrl, Verrl,Neerrl aa=where(alt1arrs eq -9,bb) & if bb ne 0 then stop end pro search_other_pos,t21v,t22arrs,itime2,ierr tmin=t21v(0) & tmax=tmin+5 t=reform(t22arrs(0,*)) tmp=where(t ge tmin and t le tmax,n) if n lt 1 then ierr=1 else ierr=0 itime2=tmp(0) end pro make_arr_ind_sub,itime,alt1arrs,altv altv=reform(alt1arrs(*,itime)) end pro make_arr_ind,itime,alt1arrs,t11arrs,t21arrs,A1arrs,E1arrs,q1arrs,$ vlsr1arrs,Ne1,v1err,Ne1err,$ altv,t11v,t21v,A1v,E1v,q1v,v1v,Ne1v,v1errv,Ne1errv make_arr_ind_sub,itime,alt1arrs,altv make_arr_ind_sub,itime,t11arrs,t11v make_arr_ind_sub,itime,t21arrs,t21v make_arr_ind_sub,itime,A1arrs,A1v make_arr_ind_sub,itime,E1arrs,E1v make_arr_ind_sub,itime,q1arrs,q1v make_arr_ind_sub,itime,vlsr1arrs,v1v make_arr_ind_sub,itime,Ne1,Ne1v make_arr_ind_sub,itime,v1err,v1errv make_arr_ind_sub,itime,Ne1err,Ne1errv end pro make_arr_2_sub,ih,ie,iem,ib,v1v,v1e,v1em,v1b,vion if ib ne -1 then vion=[v1v(ih),v1e(ie),v1em(iem),v1b(ib)] else $ vion=[v1v(ih),v1e(ie),v1em(iem),-999] end pro print_alt,ie,iem,ib,alt,alte,altem,altb if ib ne -1 then print,format='(A12,4I4)',$ 'Altitudes = ',alt,alte(ie),altem(iem),altb(ib) if ib eq -1 then print,format='(A12,3I4)',$ 'Altitudes = ',alt,alte(ie),altem(iem) end pro get_ialt,pulse,ih,altv,alte,ie common makebi,ibi_ver da=1 & if ibi_ver eq 3 then da=1.5 if pulse eq 'l' then da=11 alt=altv(ih) & a1=alt-da & a2=alt+da tmp=where(alte ge a1 and alte le a2,n) ;if n ge 2 then print,format='(3F6.1,A9)',alt,alte(tmp),'get_ialt' if n eq 2 and da ne 1 then begin a3=abs(alt-alte(tmp(0))) & a4=abs(alt-alte(tmp(1))) if a4 lt a3 then tmp(0)=tmp(1) endif ie=tmp(0) end pro get_ialt_main,pulse,ih,altv,alte,altem,altb,ie,iem,ib get_ialt,pulse,ih,altv,alte,ie & get_ialt,pulse,ih,altv,altem,iem get_ialt,pulse,ih,altv,altb,ib end pro put_nmax,altmono,nmax,ik if not keyword_set(ik) then ik=1 n=n_elements(altmono) if n gt nmax then begin nmax=n & if ik then print,format='(A7,I3)','nmax = ',n endif end pro make_altmono,alt,altmono if altmono(0) eq -1 then altmono=alt else altmono=[altmono,alt] end pro make_arr_2,ih,ie,iem,ib,v1v,v1e,v1em,v1b,q1v,q1e,q1em,q1b,$ E1v,E1e,E1em,E1b,A1v,A1e,A1em,A1b,Nev,Nee,Neem,Neb,$ vmono,qmono,Emono,Amono,Nemono make_arr_2_sub,ih,ie,iem,ib,v1v,v1e,v1em,v1b,vmono ; vion make_arr_2_sub,ih,ie,iem,ib,q1v,q1e,q1em,q1b,qmono ; q make_arr_2_sub,ih,ie,iem,ib,E1v,E1e,E1em,E1b,Emono ; El make_arr_2_sub,ih,ie,iem,ib,A1v,A1e,A1em,A1b,Amono ; Az make_arr_2_sub,ih,ie,iem,ib,Nev,Nee,Neem,Neb,Nemono ; Ne end pro make_arr_mono,n1,nmax,time_mono,alt_mono,q_mono,$ v_mono,Ne_mono, verr_mono,Neerr_mono ; nmax= alt, and n1t = time time_mono=intarr(n1) & alt_mono=intarr(nmax,n1) q_mono=intarr(nmax,n1) & v_mono=fltarr(3,nmax,n1) Ne_mono=intarr(nmax,n1) & verr_mono=fltarr(3,nmax,n1) Neerr_mono=intarr(nmax,n1) end pro put_data_mono,itime,ih,alt,qq,Vxyz,Vxyzerr,Ne0,Neerr,alt_mono,q_mono,v_mono,$ Ne_mono,verr_mono,Neerr_mono alt_mono(ih,itime)=alt & q_mono(ih,itime)=qq & v_mono(*,ih,itime)=Vxyz Ne_mono(ih,itime)=Ne0 verr_mono(*,ih,itime)=Vxyzerr Neerr_mono(ih,itime)=Neerr end pro get_cp2,date,cp if date lt 850000 then cp='cp2b' if date ge 850000 and date lt 861111 then cp='cp2c' if date ge 861111 and date lt 920707 then cp='cp2d' if date ge 920707 then cp='cp2e' end pro read_date,date,path,year,ext2,aclp ; ext2 is new on July 24, 2000 ; aclp is new on October 28, 2000 machine='stesun6' openr,1,'machine.dat' & readf,1,machine close,1 print,machine date=0l if machine eq 'stesun6' then begin openr,1,'monoesr.date' & readf,1,format='(I6)',date ext2='no' & if EOF(1) ne 1 then readf,1,ext2 aclp='ac' if EOF(1) ne 1 then readf,1,aclp k=strpos(ext2,';') if k ne -1 then begin temp=strmid(ext2,0,k-1) & ext2=strcompress(temp,/REMOVE_ALL) endif close,1 year='19'+strcompress(string(date/10000),/REMOVE_ALL) path='/work14/'+year+'/' path='./' & print,'data path = ',path endif else path='./' end pro average_ne,qmono,nemono,neave q=qmono(0:2) & tmp=where(q eq 0,n) if n eq 0 then begin neave=0 & return endif neave=fix(total(nemono(tmp))/n) end ; 990717 pro average_neerr,qmono,nemono,neave q=qmono(0:2) & tmp=where(q eq 0,n) if n eq 0 then begin neave=0 & return endif neave=fix(total(nemono(tmp))/n) ne0=nemono(tmp) tmp2=where(ne0 gt 0,n2) if n2 eq 0 then begin neave=-32768 return endif else neave=fix(total(ne0(tmp2))/n2) end pro cal_vel_main,date,filenameV,pulse, $ alt1arrs,t11arrs,t21arrs,A1arrs,E1arrs,q1arrs,vlsr1arrs,Ne1,$ alt2arrs,t12arrs,t22arrs,A2arrs,E2arrs,q2arrs,vlsr2arrs,Ne2,$ alt3arrs,t13arrs,t23arrs,A3arrs,E3arrs,q3arrs,vlsr3arrs,Ne3,$ alt4arrs,t14arrs,t24arrs,A4arrs,E4arrs,q4arrs,vlsr4arrs,Ne4,$ verr1,Neerr1,verr2,Neerr2,verr3,Neerr3,verr4,Neerr4 common makebi,ibi_ver ; ; Take data at the lowest altitude if there are more than two canditates ; Since 990716, however, for long-pulse and alternating code (make_bi_guisdap_3) ; the data at the closest the altitude will be taken. ; for irealloop=1,2 do begin n1=n_elements(alt1arrs(0,*)) & n2=n_elements(alt1arrs(*,0))-1 ; (60,41)=(alt,time) nmax=0 ;for altmono icount_data=0 & mono_cal_count=0 ionfap=0 for itime=0, n1-1 do begin ;vertical make_arr_ind,itime,alt1arrs,t11arrs,t21arrs,A1arrs,E1arrs,q1arrs,$ vlsr1arrs,Ne1,verr1,Neerr1,altv,t11v,t21v,A1v,E1v,q1v,v1v,Nev,v1errv,Neerrv ;east search_other_pos,t21v,t22arrs,itimee,ierr ierra=ierr & if ierr eq 0 then $ make_arr_ind,itimee,alt2arrs,t12arrs,t22arrs,A2arrs,E2arrs,q2arrs,$ vlsr2arrs,Ne2,verr2,Neerr2,$ alte,t11e,t21e,A1e,E1e,q1e,v1e,Nee,v1erre,Neerre ;eastmost search_other_pos,t21v,t23arrs,itimeem,ierr ierra=[ierra,ierr] & if ierr eq 0 then $ make_arr_ind,itimeem,alt3arrs,t13arrs,t23arrs,A3arrs,E3arrs,q3arrs,$ vlsr3arrs,Ne3,verr3,Neerr3,$ altem,t11em,t21em,A1em,E1em,q1em,v1em,Neem,v1errem,Neerrem ;field aligned search_other_pos,t21v,t24arrs,itimeb,ierr ierra=[ierra,ierr] & if ierr eq 0 then $ make_arr_ind,itimeb,alt4arrs,t14arrs,t24arrs,A4arrs,E4arrs,q4arrs,$ vlsr4arrs,Ne4,verr4,Neerr4, $ altb,t11b,t21b,A1b,E1b,q1b,v1b,Neb,v1errb,Neerrb set_nmax,pulse,altv,alte,altem,altb,nmax ; new on 990716 ; until here, data are put on right time and ready for combining together if total(ierra) ne 0 then goto, TLOOP ; no data at other antenna pos. if irealloop eq 2 then begin time_mono(icount_data)=t11v(0) print_times,itime,200,t11v(0),t11e(0),t11em(0),t11b(0) ; ^every 200^ endif altmono=-1 for iheight=0,n2 do begin ; loop of height ih=iheight & alt=altv(ih) ; search data accoring to ALTV get_ialt_main,pulse,ih,altv,alte,altem,altb,ie,iem,ib if ie eq -1 or iem eq -1 then goto, HLOOP make_arr_2,ih,ie,iem,ib,v1v,v1e,v1em,v1b,q1v,q1e,q1em,q1b,$ E1v,E1e,E1em,E1b,A1v,A1e,A1em,A1b,Nev,Nee,Neem,Neb,$ vmono,qmono,Emono,Amono,Nemono ; Newly made on July 17, 1999 for calculating errors make_arr_2,ih,ie,iem,ib,v1errv,v1erre,v1errem,v1errb,q1v,q1e,q1em,q1b,$ E1v,E1e,E1em,E1b,A1v,A1e,A1em,A1b,Neerrv,Neerre,Neerrem,Neerrb,$ vmonoerr,qmonoerr,Emonoerr,Amonoerr,Nemonoerr make_altmono,alt,altmono if irealloop eq 2 then begin mono_cal_count= mono_cal_count+1 ;mono_cal,vmono,qmono,Emono,Amono,Vxyz,qq,fap=ionfap,$ ; count=mono_cal_count mono_cal2,vmono,vmonoerr,qmono,Emono,Amono,Vxyz,Vxyzerr,qq,fap=ionfap,$ count=mono_cal_count ;print,format='(8F7.1)',Amono,Emono average_ne,qmono,nemono,neave average_neerr,qmono,Nemonoerr,Neerrave put_data_mono,icount_data,ih,alt,qq,Vxyz,Vxyzerr,neave,Neerrave,alt_mono,$ q_mono,v_mono,Ne_mono,verr_mono,Neerr_mono if ibi_ver eq 3 then altf_mono(ih,itime)=alt endif HLOOP: endfor ; iheight ;put_nmax,altmono,nmax,0 ;990716 comment out icount_data=icount_data+1 ;print,format='(2I4)',icount_data,itime TLOOP: endfor ; itime loop print,format='(5I6)',itime,t11v(0),t11e(0),t11em(0),t11b(0) if irealloop eq 1 then begin make_arr_mono,icount_data,nmax,time_mono,alt_mono,q_mono,v_mono,Ne_mono,$ verr_mono,Neerr_mono if ibi_ver eq 3 then altf_mono=fltarr(nmax,n1) endif endfor ; irealloop n1=long(n_elements(time_mono)) n2=long(n_elements(alt_mono(*,0))) openw_print,3,filenameV print,format='(A5,I5,A6,I5)','n1 = ',n1,' n2 = ',n2 writeu,3,n1,n2,time_mono,alt_mono,q_mono,v_mono,Ne_mono,verr_mono,Neerr_mono if ibi_ver eq 3 then writeu,3,altf_mono close,3 ; filenameV end pro read_head3,imake_bi_guisdap_ver iprint=0 icount=0 for i=1,8 do begin moji='moji' & readf,1,' ',moji & if iprint then print,i,moji if strmid(moji,0,1) eq 'A' then icount=icount+1 endfor if icount eq 8 then imake_bi_guisdap_ver=3 else imake_bi_guisdap_ver=2 end pro read_bi_guisdap_3,sdate0,site,cp,dat,path,iarr,farr,ibi_ver filename=site+sdate0+cp+dat ierror=0 ON_IOERROR, J1 J0: print,'----- Open (r) = ',path+filename openr,1,path+filename read_head3,ibi_ver n1=0l & n2=0l & icount=0l readu,1,icount,n1,n2 & print,icount,n1,n2 iarr=intarr(n1,icount) & farr=fltarr(n2,icount) if n2 eq 1 then farr=fltarr(icount) readu,1,iarr,farr & close,1 return J1: print,'---- File error ',!error if ierror eq 0 then begin path='/EISCAT04/cp1data/tromso/' ierror=ierror+1 goto, J0 endif print,!error stop end ;* der_guisdap.pro pro der_guisdap_bi,iarr,farr,date,time1,time2,alt,qq,Ne0,Vi0,$ Ti0,Te0,Neerr,Vierr,Tierr,Teerr,Az0,El0,ration_O_Ne,colf,range,$ time1sec,time2sec date=reform(iarr(0,*)) ; date time1=reform(iarr(1,*)) ; time1 time2=reform(iarr(2,*)) ; time2 alt=reform(iarr(3,*)) ; alt qq=reform(iarr(4,*)) ; q Ne0=reform(iarr(5,*)) ; Ne Vi0=reform(iarr(6,*)) ; Vi Ti0=reform(iarr(7,*)) ; Ti Te0=reform(iarr(8,*)) ; Te Neerr=reform(iarr(9,*)) ; Neerr Vierr=reform(iarr(10,*)) ; Vierr Tierr=reform(iarr(11,*)) ; Tierr Teerr=reform(iarr(12,*)) ; Teerr Az0=reform(iarr(13,*)) ; Az El0=reform(iarr(14,*)) ; El ration_O_Ne=reform(farr(*)) ; ration O+/Ne colf=reform(iarr(15,*)) ; colision fr. range=reform(iarr(16,*)) ; range n=n_elements(iarr(*,0)) if n eq 17 then return ; very old version before say 1998 time1sec=reform(iarr(17,*)) ; time1 sec time2sec=reform(iarr(18,*)) ; time2 sec end ;* der_guisdap.pro pro der_guisdap_bi_3,iarr,farr,date,time1,time2,alt,qq,Ne0,Vi0,$ Ti0,Te0,Neerr,Vierr,Tierr,Teerr,Az0,El0,ration_O_Ne,colf,range,$ time1sec,time2sec,altsf,lat0,lon0,Tsys,Txpw date=reform(iarr(0,*)) ; date time1=reform(iarr(1,*)) ; time1 time2=reform(iarr(2,*)) ; time2 alt=reform(iarr(3,*)) ; alt qq=reform(iarr(4,*)) ; q Ne0=reform(iarr(5,*)) ; Ne Vi0=reform(iarr(6,*)) ; Vi Ti0=reform(iarr(7,*)) ; Ti Te0=reform(iarr(8,*)) ; Te Neerr=reform(iarr(9,*)) ; Neerr Vierr=reform(iarr(10,*)) ; Vierr Tierr=reform(iarr(11,*)) ; Tierr Teerr=reform(iarr(12,*)) ; Teerr Az0=reform(iarr(13,*)) ; Az El0=reform(iarr(14,*)) ; El ration_O_Ne=reform(farr(0,*)) ; ration O+/Ne colf=reform(iarr(15,*)) ; colision fr. range=reform(iarr(16,*)) ; range n=n_elements(iarr(*,0)) time1sec=reform(iarr(17,*)) ; time1 sec time2sec=reform(iarr(18,*)) ; time2 sec lat0=reform(iarr(19,*)) ; latitude lon0=reform(iarr(20,*)) ; longitude altsf=reform(farr(1,*)) ; ration O+/Ne n0=n_elements(farr(*,0)) if n0 ge 3 then Tsys=reform(farr(2,*)) else Tsys=0 ; Tsys n1=n_elements(iarr(*,0)) if n1 ge 22 then Txpw=reform(iarr(21,*)) else Txpw=-1 ; Tx Power end pro where_az_el_970106,sdate0,arr0_Az_2,arr0_El_2,Azs,Els,alts,ipos2,ipos3 if sdate0 ne '970106' and sdate0 ne '980921' then return where_az_el,arr0_Az_2,arr0_El_2,Azs,Els,alts,ipos1_2,ipos2_2,ipos3_2,ipos4_2 ipos2=[ipos2,ipos2_2] & ipos3=[ipos3,ipos3_2] ipos2=ipos2(sort(ipos2)) & ipos3=ipos3(sort(ipos3)) end ; made on 990716 pro set_nmax,pulse,altv,alte,altem,altb,nmax maxalt1=max(alte) & maxalt2=max(altem) & maxalt3=max(altb) ;print,maxalt1,maxalt2,maxalt3 maxalt4=min([maxalt1,maxalt2,maxalt3]) da=1 & if pulse eq 'l' then da=11 ; see get_ialt.pro tmp=where(altv le maxalt4+da,nmax0) if nmax0 gt nmax then begin nmax=nmax0 & print,format='(A7,I3)','nmax = ',nmax endif end pro print_times,itime,i,t11v,t11e,t11em,t11b if (itime mod i) ne 0 then return print,format='(5I6)',itime,t11v,t11e,t11em,t11b end pro plot_check1,sdate0,Els2 setps,12 !p.multi=[0,1,15] i=0l while i lt n_elements(Els2)-100 do begin plot,Els2(i:i+100),yrange=[60,90],title=sdate0 i=i+100 endwhile end function time_to_min,t t2=(t/100)*60 + (t mod 100) & return,t2 end pro time_to_min_all,t1arrs,t2arrs,t3arrs,t4arrs common timeback,t1,t2,t3,t4 t1=t1arrs & t2=t2arrs & t3=t3arrs & t4=t4arrs t1arrs=time_to_min(t1arrs) t2arrs=time_to_min(t2arrs) t3arrs=time_to_min(t3arrs) t4arrs=time_to_min(t4arrs) end pro print_check2,i,i2,i3,i4,ndatacount common timeback,t1arrsorg,t2arrsorg,t3arrsorg,t4arrsorg k=-1 t1=t1arrsorg mod 2400 t2=t2arrsorg mod 2400 t3=t3arrsorg mod 2400 t4=t4arrsorg mod 2400 if i2 eq -1 and i3 eq -1 and i4 eq -1 then begin print,format='(9I6)',ndatacount,t3arrsorg(i),k,k,k,t3(i),k,k,k return endif if i2 ne -1 and i3 eq -1 and i4 eq -1 then begin print,format='(9I6)',ndatacount,t3arrsorg(i),t2arrsorg(i2),k,k,t3(i),t2(i2),k,k return endif if i2 ne -1 and i3 eq -1 and i4 ne -1 then begin print,format='(9I6)',ndatacount,t3arrsorg(i),t2arrsorg(i2),k,t1arrsorg(i4),t3(i),$ t2(i2),k,t1(i4) return endif if i2 eq -1 and i3 ne -1 and i4 eq -1 then begin print,format='(9I6)',ndatacount,t3arrsorg(i),k,t4arrsorg(i3),k,t3(i),k,t4(i3),k return endif if i2 eq -1 and i3 ne -1 and i4 ne -1 then begin print,format='(9I6)',ndatacount,t3arrsorg(i),k,t4arrsorg(i3),t4arrsorg(i3),t3(i),$ k,t4(i3),t1(i4) return endif if i4 ne -1 then print,format='(9I6)',ndatacount,t3arrsorg(i),t2arrsorg(i2),t4arrsorg(i3),$ t1arrsorg(i4),t3(i),t2(i2),t4(i3),t1(i4) if i4 eq -1 then print,format='(9I6)',ndatacount,t3arrsorg(i), t2arrsorg(i2),t4arrsorg(i3),$ k,t3(i),t2(i2),t4(i3),k end function check_alt,a0,alt,tmp1 n=n_elements(tmp1) dalt=abs(alt(tmp1)-a0) i=sort(dalt) t=tmp1(i(0)) ;print,a0,alt(t),alt(tmp1) return,t end function search_alt,alt01,n1,altst3,dalt a1=altst3-dalt & a2=altst3+dalt tmp1=where(alt01 gt a1 and alt01 lt a2,n1) if n1 ge 2 then tmp1=check_alt(altst3,alt01,tmp1) return,tmp1 end pro search_alt_main,altst3,dalt,alt01,alt02,alt03,tmp1,tmp2,tmp3,n1,n2,n3 tmp1=search_alt(alt01,n1,altst3,dalt) ; bit west tmp2=search_alt(alt02,n2,altst3,dalt) ; West tmp3=search_alt(alt03,n3,altst3,dalt) ; FA end pro set_data,t,alts,q,V,Verr,alt1,q1,v1,verr1 alt1=alts(*,t) q1=q(*,t) v1=v(*,t) verr1=verr(*,t) end pro set_data2,t,Ne0,Ti0,Te0,Ne1,Ti1,Te1 Ne1=Ne0(*,t) & Ti1=Ti0(*,t) & Te1=Te0(*,t) end pro set_par,alt3arrs,qq3arrs,vlsr3arrs,vlsrerr3arrs,i,alt00,q0,v0,v0err ; vertical alt00=reform(alt3arrs(*,i)) q0=reform(qq3arrs(*,i)) v0=reform(vlsr3arrs(*,i)) v0err=reform(vlsrerr3arrs(*,i)) end pro set_par2,d1,d2,d3,i,r1,r2,r3 r1=reform(d1(*,i)) & r2=reform(d2(*,i)) & r3=reform(d3(*,i)) end pro make_output_arrays,nend,time_all,alt_all,vxyz_all,vxyzerr_all,vmono_all,vmonoerr_all,$ ; qmono_all,temp_all,temp_av,temp_sigma,temperr_all,temp_n qmono_all,temp_all,temp_av,temp_sigma,temperr_all,temp_n,latmono,lonmono,Nemono,Neerrmono common ext,ext2,numalt ; newly added on July 24, 2000 ; the second line of monoesr.date i=numalt time_all=intarr(nend) & alt_all=fltarr(4,i,nend) vxyz_all=fltarr(3,i,nend) & vxyzerr_all=fltarr(3,i,nend) vmono_all=intarr(4,i,nend) & vmonoerr_all=intarr(4,i,nend) qmono_all=intarr(4,i,nend) & qmono_all(*,*,*)=-999 ; if -999 remains, it means no data. temp_all=intarr(4,i,nend) & temp_av=fltarr(i,nend) temp_sigma=fltarr(i,nend) & temperr_all=intarr(4,i,nend) temp_n=intarr(i,nend) latmono=intarr(4,i,nend) & lonmono=intarr(4,i,nend) Nemono=intarr(4,i,nend) & Neerrmono=intarr(4,i,nend) end pro average_temp,qmono,tmono,tmonoerr,tmonoave,tmonosigma,tmonon tmp=where(qmono eq 0 and tmono ge 100 and tmonoerr gt 0 and tmonoerr lt 500,n) if n eq 0 then begin tmonon=0 & tmonoave=-1.0 & tmonosigma=-1.0 & return endif tmonon=n tmonoave=(total(tmono(tmp)))/n if n eq 1 then tmonosigma=0 t=(tmono(tmp)-tmonoave)^2 if n ne 1 then tmonosigma=sqrt(total(t)/(n-1)) end ;******************************************************************** ;----- Main, main, MAIN ----- ;******************************************************************** proname='mono_esr.pro' & print_proname,proname common guisdap,igup,dat common makebi,ibi_ver common timeback,t1arrsorg,t2arrsorg,t3arrsorg,t4arrsorg common ext,ext2,numalt ; newly added on July 24, 2000 ; the second line of monoesr.date ; First edit "monoesr.date" to specify date read_date,date,path,year,ext2,aclp & cp='cp2l' ; ext2 is new for no of iteration as well as Ti=Te-altitude settings year=fix(date/10000l) & sdate0=strcompress(date,/REMOVE_ALL) outfile='vmono'+sdate0+'esr.dat' if ext2 ne 'no' then outfile='vmono'+sdate0+'esr'+ext2+'.dat' readesr_cp2, sdate0,aclp,t1s,t2s,t1sec,t2sec,Azs2,Els2, $ alts,altsorg,lats,lons,quos,Nearr,Vlsr,Ti0t,Te0t,$ Neerr,Vlsrerr,Tierrt,Teerrt,ranget,ration_O_Net,$ tmp1,tmp2,tmp3,tmp4,n1,n2,n3,n4 ; tmp1, tmp2, tmp3, tmp4 ; FA , bit west, Vertical, west ;alt1arrs=alts(*,tmp1) ; 4:FA ;alt2arrs=alts(*,tmp2) ; 2:bit west ;alt3arrs=alts(*,tmp3) ; 1:Vertical ;alt4arrs=alts(*,tmp4) ; 3:West set_data,tmp1,alts,quos,Vlsr,Vlsrerr,alt1arrs,qq1arrs,vlsr1arrs,vlsrerr1arrs ;fa set_data,tmp2,alts,quos,Vlsr,Vlsrerr,alt2arrs,qq2arrs,vlsr2arrs,vlsrerr2arrs ;bt west set_data,tmp3,alts,quos,Vlsr,Vlsrerr,alt3arrs,qq3arrs,vlsr3arrs,vlsrerr3arrs ;v set_data,tmp4,alts,quos,Vlsr,Vlsrerr,alt4arrs,qq4arrs,vlsr4arrs,vlsrerr4arrs ;west set_data2,tmp1,Nearr,Ti0t,Te0t,Ne1arr,Ti1arr,Te1arr set_data2,tmp2,Nearr,Ti0t,Te0t,Ne2arr,Ti2arr,Te2arr set_data2,tmp3,Nearr,Ti0t,Te0t,Ne3arr,Ti3arr,Te3arr set_data2,tmp4,Nearr,Ti0t,Te0t,Ne4arr,Ti4arr,Te4arr set_data2,tmp1,Neerr,Tierrt,Teerrt,Ne1err,Ti1err,Te1err set_data2,tmp2,Neerr,Tierrt,Teerrt,Ne2err,Ti2err,Te2err set_data2,tmp3,Neerr,Tierrt,Teerrt,Ne3err,Ti3err,Te3err set_data2,tmp4,Neerr,Tierrt,Teerrt,Ne4err,Ti4err,Te4err set_data2,tmp1,lats,lons,ranget,lat1arr,lon1arr,rang1arr set_data2,tmp2,lats,lons,ranget,lat2arr,lon2arr,rang2arr set_data2,tmp3,lats,lons,ranget,lat3arr,lon3arr,rang3arr set_data2,tmp4,lats,lons,ranget,lat4arr,lon4arr,rang4arr t1arrs=t1s(tmp1) & t2arrs=t1s(tmp2) & t3arrs=t1s(tmp3) & t4arrs=t1s(tmp4) time_to_min_all,t1arrs,t2arrs,t3arrs,t4arrs Azarr1=Azs2(tmp1) & Azarr2=Azs2(tmp2) & Azarr3=Azs2(tmp3) & Azarr4=Azs2(tmp4) Elarr1=Els2(tmp1) & Elarr2=Els2(tmp2) & Elarr3=Els2(tmp3) & Elarr4=Els2(tmp4) ;plot_check1,sdate0,Els2 case sdate0 of '980921':cycle=8 ; until 1635 on 980921 cycle=4; since 2212 on 980922 cycle=12 '990308':cycle=8 else :cycle=10 end ; Vertical position is used as standard. t=t3arrs ; vertical nend=n_elements(Elarr3) make_output_arrays,nend,time_all,alt_all,vxyz_all,vxyzerr_all,vmono_all,vmonoerr_all,$ qmono_all,temp_all,temp_av,temp_sigma,temperr_all,temp_n,latmono_all,lonmono_all,$ Nemono_all,Neerrmono_all altst=alt3arrs ; standard altitute is at vertical position ndatacount=0l icount=0l ;stop ; check here latmono, lonmono for iloop=0,nend-1 do begin i=iloop t0=t(i) time_all(i)=t0 tend=t(i)+cycle tmp1=where(t2arrs gt t0 and t2arrs lt tend,n1) ; bit west tmp2=where(t4arrs gt t0 and t4arrs lt tend,n2) ; west tmp3=where(t1arrs gt t0 and t1arrs lt tend,n3) ; FA i1=i & i2=tmp1(0) & i3=tmp2(0) & i4=tmp3(0) print_check2,i1,i2,i3,i4,ndatacount if n1 ne 1 or n2 ne 1 then goto, JUMP1 ndatacount=ndatacount+1 set_par,alt3arrs,qq3arrs,vlsr3arrs,vlsrerr3arrs,i,alt00,q0,v0,v0err ; vertical set_par,alt2arrs,qq2arrs,vlsr2arrs,vlsrerr2arrs,i2,alt01,q1,v1,v1err ; bit west set_par,alt4arrs,qq4arrs,vlsr4arrs,vlsrerr4arrs,i3,alt02,q2,v2,v2err ; west set_par2,Ne3arr,Ti3arr,Te3arr,i,Ne00,Ti00,Te00 set_par2,Ne2arr,Ti2arr,Te2arr,i2,Ne01,Ti01,Te01 set_par2,Ne4arr,Ti4arr,Te4arr,i3,Ne02,Ti02,Te02 set_par2,Ne3err,Ti3err,Te3err,i,Ne00err,Ti00err,Te00err set_par2,Ne2err,Ti2err,Te2err,i2,Ne01err,Ti01err,Te01err set_par2,Ne4err,Ti4err,Te4err,i3,Ne02err,Ti02err,Te02err set_par2,lat3arr,lon3arr,rang3arr,i,lat00,lon00,rang00 ; Vertical set_par2,lat2arr,lon2arr,rang2arr,i2,lat01,lon01,rang01 ; bit west set_par2,lat4arr,lon4arr,rang4arr,i3,lat02,lon02,rang02 ; west if i4 ne -1 then begin Amono=[Azarr3(i1),Azarr2(i2),Azarr4(i3),Azarr1(i4)] Emono=[Elarr3(i1),Elarr2(i2),Elarr4(i3),Elarr1(i4)] set_par,alt1arrs,qq1arrs,vlsr1arrs,vlsrerr1arrs,i4,alt03,q3,v3,v3err ; FA set_par2,Ne1arr,Ti1arr,Te1arr,i4,Ne03,Ti03,Te03 set_par2,Ne1err,Ti1err,Te1err,i4,Ne03err,Ti03err,Te03err set_par2,lat1arr,lon1arr,rang1arr,i4,lat03,lon03,rang03 ; FA endif else begin print,'No data at Field-aligned position.' Amono=[Azarr3(i1),Azarr2(i2),Azarr4(i3),180.6] Emono=[Elarr3(i1),Elarr2(i2),Elarr4(i3),81.6] q3(*)=-1 ; assuming q3 existing (parameters come from the previous one) ;lat03=0 & lon03=0 ; asume previous one exist, if not then error would occur Ne03(*)=0 & Ne03err(*)=0 endelse for ialt=0,n_elements(alt00)-1 do begin altst3=alt00(ialt) dalt=25.0 search_alt_main,altst3,dalt,alt01,alt02,alt03,tmp1,tmp2,tmp3,n1,n2,n3 ;print,format='(4F7.1)',altst3,alt01(tmp1(0)),alt02(tmp2(0)),alt03(tmp3(0)) if n1 eq 0 or n2 eq 0 or n3 eq 0 then goto, JUMP2 i2=tmp1(0) & i3=tmp2(0) & i4=tmp3(0) vmono=[v0(ialt),v1(i2),v2(i3),v3(i4)] vmonoerr=[v0err(ialt),v1err(i2),v2err(i3),v3err(i4)] qmono=[q0(ialt),q1(i2),q2(i3),q3(i4)] mono_cal2,vmono,vmonoerr,qmono,Emono,Amono,Vxyz,Vxyzerr,qq,fap=0 tmono=[Ti00(ialt),Ti01(i2),Ti02(i3),Ti03(i4)] tmonoerr=[Ti00err(ialt),Ti01err(i2),Ti02err(i3),Ti03err(i4)] temp_all(*,ialt,i)=tmono temperr_all(*,ialt,i)=tmonoerr average_temp,qmono,tmono,tmonoerr,tmonoave,tmonosigma,tmonon temp_av(ialt,i)=tmonoave temp_sigma(ialt,i)=tmonosigma temp_n(ialt,i)=tmonon alt_all(*,ialt,i)=[altst3,alt01(i2),alt02(i3),alt03(i4)] ;v, bw,w,fa latmono_all(*,ialt,i)=[lat00(ialt),lat01(i2),lat02(i3),lat03(i4)] ;v, bw,w,fa lonmono_all(*,ialt,i)=[lon00(ialt),lon01(i2),lon02(i3),lon03(i4)] ;v, bw,w,fa Nemono_all(*,ialt,i)=[Ne00(ialt),Ne01(i2),Ne02(i3),Ne03(i4)] ;v, bw,w,fa Neerrmono_all(*,ialt,i)=[Ne00err(ialt),Ne01err(i2),Ne02err(i3),Ne03err(i4)] ;v, bw,w,fa vxyz_all(*,ialt,i)=Vxyz vxyzerr_all(*,ialt,i)=Vxyzerr vmono_all(*,ialt,i)=vmono vmonoerr_all(*,ialt,i)=vmonoerr qmono_all(*,ialt,i)=qmono JUMP2: ;altloop endfor ; ialt icount=icount+1l if icount eq 1 then altstandard=alt00 JUMP1: endfor ; time loop openw,1,outfile & print,'----- Open (w) = ',outfile i1=long(n_elements(vxyz_all(0,*,0))) ; altitute i2=long(n_elements(vxyz_all(0,0,*))) ; time writeu,1,i1,i2,time_all,vxyz_all,vxyzerr_all,altstandard,alt_all,vmono_all,vmonoerr_all,$ qmono_all,latmono_all,lonmono_all,Nemono_all,Neerrmono_all close,1 outfile2='ti'+sdate0+cp+'_esr.dat' if ext2 ne 'no' then outfile2='ti'+sdate0+cp+'_esr'+ext2+'.dat' openw,1,outfile2 & print,'----- Open (w) = ',outfile2 i1=long(n_elements(temp_all(0,*,0))) ; altitute i2=long(n_elements(temp_all(0,0,*))) ; time writeu,1,i1,i2,time_all,altstandard,alt_all,temp_av,temp_sigma,temp_n,temp_all,$ temperr_all,latmono_all,lonmono_all,qmono_all ; add latmono_all,lonmono_all on 001201 close,1 end