[cig-commits] r4481 - geodyn/3D/MAG/trunk/idl
polson at geodynamics.org
polson at geodynamics.org
Wed Sep 6 11:44:10 PDT 2006
Author: polson
Date: 2006-09-06 11:44:09 -0700 (Wed, 06 Sep 2006)
New Revision: 4481
Added:
geodyn/3D/MAG/trunk/idl/magts.pro
Log:
new idl for mag l-files
Added: geodyn/3D/MAG/trunk/idl/magts.pro
===================================================================
--- geodyn/3D/MAG/trunk/idl/magts.pro 2006-09-06 18:39:06 UTC (rev 4480)
+++ geodyn/3D/MAG/trunk/idl/magts.pro 2006-09-06 18:44:09 UTC (rev 4481)
@@ -0,0 +1,541 @@
+;PROCEDURE MAGTS takes data from l-files generated by MAG
+;and creates time series plots and statistics.
+;This version reads an l-file consisting of 17 time series,
+;the first record being dimensionless time. Energies and
+;rms magnetic field and velocity are scaled as in MAG; tilt is dipole
+;vector colatitude; pole longitude is dipole vector longitude
+
+
+;SET UP X-WINDOW FOR COLOR DISPLAY
+SET_PLOT,'X'
+;SET DEVICE FOR COLOR on Mac
+DEVICE, RETAIN=2, DECOMPOSED=0
+
+;CONSTANTS
+pi=3.14159
+
+;ARRAY DIMENSIONS
+bigarr = FLTARR(17,6000)
+arr = FLTARR(17)
+
+;SCALE FACTORS (default)
+tscale=1.0
+smsize=1.0
+
+!P.MULTI=0
+
+;READ TIMESERIES L-FILE
+fname = ' ' ;filename of timeseries data
+PRINT,' '
+PRINT,'Enter l-file name:'
+READ, fname
+
+;READ TIMESERIES DATA FROM l-file FOR n TIMESTEPS
+OPENR, 1, fname
+n = 0
+WHILE (NOT EOF(1)) DO BEGIN
+ n = n + 1
+ READF, 1, arr
+ bigarr(*, n-1) = arr
+ENDWHILE
+CLOSE, 1
+
+PRINT, ' '
+PRINT, 'NUMBER OF STEPS=',n
+PRINT, ' '
+
+;DEFINE TIME SERIES ARRAYS
+tseries = FLTARR(17, n)
+tseries = bigarr(*, 0:n-1)
+
+time = tseries(0,*)
+ke = tseries(1,*)
+kepol = tseries(2,*)
+me = tseries(3,*)
+mepol = tseries(4,*)
+keaxtor = tseries(5,*)
+keaxpol = tseries(6,*)
+meaxipol = tseries(7,*)
+meaxitor = tseries(8,*)
+nutop = tseries(9,*)
+nubot = tseries(10,*)
+bmean = tseries(11,*)
+dipole = tseries(12,*)
+dipaxi = tseries(13,*)
+tilt = tseries(14,*)
+diplong = tseries(15,*)
+vmean = tseries(16,*)
+
+;UUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUU
+;INTERPOLATE TO UNIFORM SERIES
+dtime=time(n-1)-time(0)
+;create uniform time series
+timeu=time(0)+(dtime*findgen(n)/n)
+;create new y-data
+resu=interpol(ke,time,timeu) & ke=resu
+resu=interpol(kepol,time,timeu) & kepol=resu
+resu=interpol(me,time,timeu) & me=resu
+resu=interpol(mepol,time,timeu) & mepol=resu
+resu=interpol(keaxtor,time,timeu) & keaxtor=resu
+resu=interpol(meaxipol,time,timeu) & meaxipol=resu
+resu=interpol(meaxitor,time,timeu) & meaxitor=resu
+resu=interpol(nutop,time,timeu) & nutop=resu
+resu=interpol(nubot,time,timeu) & nubot=resu
+resu=interpol(bmean,time,timeu) & bmean=resu
+resu=interpol(vmean,time,timeu) & vmean=resu
+resu=interpol(dipole,time,timeu) & dipole=resu
+resu=interpol(dipaxi,time,timeu) & dipaxi=resu
+resu=interpol(tilt,time,timeu) & tilt=resu
+resu=interpol(diplong,time,timeu) & diplong=resu
+
+;replace old ts
+res=resu & time=timeu
+
+
+;UUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUUU
+
+
+;SSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSS
+;COMPUTE STATISTICS
+;find extreme values
+kemax=max(ke) & kemin=min(ke)
+memax=max(me) & memin=min(me)
+mepmax=max(mepol) & mepmin=min(mepol)
+apmax=max(meaxipol) & apmin=min(meaxipol)
+atmax=max(meaxitor) & atmin=min(meaxitor)
+dmax=max(dipole) & dmin=min(dipole)
+admax=max(dipaxi) & admin=min(dipaxi)
+nutmax=max(nutop) & nutmin=min(nutop)
+nubmax=max(nubot) & nubmin=min(nubot)
+vmax=max(vmean) & vmin=min(vmean)
+bmax=max(bmean) & bmin=min(bmean)
+tilmax=max(tilt) & tilmin=min(tilt)
+lonmax=max(diplong) & lonmin=min(diplong)
+
+PRINT, ' '
+PRINT, 'Parameter Average Stnd Dev Max Min'
+res1=moment(ke,sdev=sd1) & print, 'KE: ' ,res1(0),sd1,kemax,kemin
+res2=moment(kepol,sdev=sd2) & print, 'KE(pol): ',res2(0),sd2
+res3=moment(me,sdev=sd3) & print, 'ME: ', res3(0),sd3,memax,memin
+res4=moment(mepol,sdev=sd4) & print, 'ME(pol): ', res4(0),sd4
+res5=moment(keaxtor,sdev=sd5) & print, 'KE(axitor): ', res5(0),sd5
+res6=moment(keaxpol,sdev=sd6) & print, 'KE(axipol): ', res6(0),sd6
+res8=moment(meaxipol,sdev=sd8) & print, 'ME(axipol): ', res8(0),sd8,apmax,apmin
+res9=moment(meaxitor,sdev=sd9) & print, 'ME(axitor): ', res9(0),sd9,atmax,atmin
+res10=moment(nutop,sdev=sd10) & print, 'Nu(top): ', res10(0),sd10,nutmax,nutmin
+res11=moment(nubot,sdev=sd11) & print, 'Nu(bot): ', res11(0),sd11,nubmax,nubmin
+res12=moment(bmean,sdev=sd12) & print, 'B(ave): ', res12(0),sd12,bmax,bmin
+res13=moment(dipole,sdev=sd13) & print, 'B(dip): ', res13(0),sd13,dmax,dmin
+res14=moment(dipaxi,sdev=sd14) & print, 'B(axi): ', res14(0),sd14,admax,admin
+res15=moment(tilt,sdev=sd15) & print, 'Tilt: ', res15(0),sd15,tilmax,tilmin
+res16=moment(diplong,sdev=sd16) & print, 'Lon(ave): , res16(0),sd16,lonmax,lonmin
+res17=moment(vmean,sdev=sd17) & print, 'V(ave): ', res17(0),sd17,vmax,vmin
+
+PRINT, ' '
+PRINT,'Correlations'
+PRINT, 'KE vs ME: ', correlate(ke,me)
+PRINT, 'KE vs Nu(top): ', correlate(ke,nutop)
+PRINT, 'KE vs Nu(bot): ', correlate(ke,nubot)
+PRINT, 'ME vs Nu(top): ', correlate(me,nutop)
+PRINT, 'ME vs Nu(bot): ', correlate(me,nubot)
+PRINT, 'dipole vs Nu(top): ', correlate(dipole,nutop)
+PRINT, 'dipole vs Nu(bot): ', correlate(dipole,nubot)
+PRINT, 'dipole vs vmean: ', correlate(dipole,vmean)
+PRINT, 'dipole vs tilt: ', correlate(dipole,tilt)
+;SSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSSS
+
+
+;#######################################################################
+; Here begins plotting material:
+;#######################################################################
+; Define window size and plotting frames
+ XDIM=18.0
+ YDIM=24.0
+ IGIFPS=0
+
+ LCOLOR = 0
+ LOADCT, LCOLOR ;loads B&W color table
+
+;** WINDOW SIZE IN PIXELS (ALSO SIZE OF GIF FILE)
+ XWINDOW=500
+ YWINDOW=575
+ SCWINDOW=1.25
+ XWIND=XWINDOW*SCWINDOW
+ YWIND=YWINDOW*SCWINDOW
+ CSZ=1.00 & CSB=1.50 & CSS=0.720 ; CHARACTER SIZES
+ SZF=1.2*SCWINDOW ; MAKES GIF AND PS CHARACTERS SAME SIZE
+
+;*** Normalized coordinates for 2x2 plots
+ XY1=[ 1.5/XDIM,14./YDIM, 8.5/XDIM,22./YDIM] & XT1=5./XDIM & YT1=22.2/YDIM
+ XP1 = 6.5/XDIM & YP1 = 20.5/YDIM
+ XY2=[ 10.5/XDIM,14./YDIM,17.5/XDIM,22./YDIM] & XT2=14./XDIM & YT2=22.2/YDIM
+ XP2 = 15.5/XDIM & YP2 = 20.5/YDIM
+ XY3=[ 1.5/XDIM, 3.5/YDIM, 8.5/XDIM,11.5/YDIM] & XT3=5./XDIM & YT3=11.7/YDIM
+ XY4=[ 10.5/XDIM, 3.5/YDIM,17.5/XDIM,11.5/YDIM] & XT4=14./XDIM & YT4=11.7/YDIM
+ XCT=8/XDIM & YCT=16/YDIM
+
+;*************************************************************************
+;************** THE PLOTTING MENU STARTS HERE **************************
+;*************************************************************************
+LABEL0: IF IGIFPS EQ 1 THEN PRINT,"USE --GOTO, LABELOUT-- STATEMENT HERE."
+ PRINT,' ' & PRINT,' ' & PRINT,' ' & PRINT,' ' & PRINT,' '
+ PRINT, "ENTER OPTION: "
+ PRINT, "EXIT PROCEDURE = -1"
+ PRINT, "ENERGIES I = 1; ENERGIES II = 2"
+ PRINT, "FIELD PLOTS = 3; STATISTICS = 4"
+ PRINT, "RESCALE = 5; SPECTRA = 6; POLE PLOT =7"
+ PRINT, "CHANGE WINDOW SIZE = 11"
+ PRINT, "ANNOTATIONS = 12"
+ PRINT, "PRINT POSTSCRIPT=21"
+ READ,IOPTION
+
+ IF (IOPTION GE 1 AND IOPTION LE 2) THEN BEGIN
+ WINDOW,0,xsize=XWIND,ysize=YWIND,title='GRAPHIC WINDOW'
+ !P.CHARTHICK=1.5
+ !P.FONT=0
+ ENDIF
+
+ CASE IOPTION OF
+ -1: GOTO, LABEL98
+ 0: GOTO, LABEL0
+ 1: GOTO, LABEL1
+ 2: GOTO, LABEL2
+ 3: GOTO, LABEL3
+ 4: GOTO, LABEL4
+ 5: GOTO, LABEL5
+ 6: GOTO, LABEL6
+ 7: GOTO, LABEL7
+ 11: GOTO, LABEL11
+ 12: GOTO, LABEL0
+ 21: GOTO, LABEL21
+ 22: GOTO, LABEL0
+ ELSE: GOTO, LABEL0
+ ENDCASE
+
+;####################################################################
+
+;11111111111111111111111111111111111111111111111111111111111111111111111
+LABEL1: & IPAGE=1
+ERASE
+!P.MULTI = [0,1,4,0]
+TL=' Scale='+string(tscale)+' Smooth='+string(smsize)
+
+PLOT, time, ke, LINESTYLE = 0, TITLE='KE', YRANGE=[kemin,kemax]
+OPLOT, time, replicate(res1(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res1(0)-sd1,n),LINESTYLE = 1
+OPLOT, time, replicate(res1(0)+sd1,n),LINESTYLE = 1
+
+PLOT, time, me, LINESTYLE = 0, TITLE='ME',YRANGE=[memin,memax]
+OPLOT, time, replicate(res3(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res3(0)-sd3,n),LINESTYLE = 1
+OPLOT, time, replicate(res3(0)+sd3,n),LINESTYLE = 1
+
+PLOT, time, nutop, LINESTYLE = 0, TITLE='Nu(top)', YRANGE=[nutmin,nutmax]
+OPLOT, time, replicate(res10(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res10(0)-sd10,n),LINESTYLE = 1
+OPLOT, time, replicate(res10(0)+sd10,n),LINESTYLE = 1
+
+PLOT, time, nubot, LINESTYLE = 0, TITLE='Nu(bot)', YRANGE=[nubmin,nubmax],$
+ XTITLE='Time'+TL
+OPLOT, time, replicate(res11(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res11(0)-sd11,n),LINESTYLE = 1
+OPLOT, time, replicate(res11(0)+sd11,n),LINESTYLE = 1
+
+
+GOTO, LABEL99
+;111111111111111111111111111111111111111111111111111111111111111111111111
+
+
+;222222222222222222222222222222222222222222222222222222222222222222222222
+LABEL2: IPAGE=2
+ERASE
+!P.MULTI = [0,1,4,0]
+TL=' Scale='+string(tscale)+' Smooth='+string(smsize)
+
+PLOT, time, mepol, LINESTYLE = 0, TITLE='ME(pol)',$
+ YRANGE=[mepmin,mepmax]
+OPLOT, time, replicate(res4(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res4(0)-sd4,n),LINESTYLE = 1
+OPLOT, time, replicate(res4(0)+sd4,n),LINESTYLE = 1
+
+PLOT, time, meaxipol, LINESTYLE = 0, TITLE='MEaxi(pol)',$
+ YRANGE=[apmin,apmax]
+OPLOT, time, replicate(res8(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res8(0)-sd8,n),LINESTYLE = 1
+OPLOT, time, replicate(res8(0)+sd8,n),LINESTYLE = 1
+
+PLOT, time, meaxitor, LINESTYLE = 0, TITLE='MEaxi(tor)',$
+ YRANGE=[atmin,atmax]
+OPLOT, time, replicate(res9(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res9(0)-sd9,n),LINESTYLE = 1
+OPLOT, time, replicate(res9(0)+sd9,n),LINESTYLE = 1
+
+PLOT, time, dipole, LINESTYLE = 0, TITLE='Dipole', XTITLE='Time'+TL,$
+ YRANGE=[dmin,dmax]
+OPLOT, time, replicate(res13(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res13(0)-sd13,n),LINESTYLE = 1
+OPLOT, time, replicate(res13(0)+sd13,n),LINESTYLE = 1
+
+GOTO, LABEL99
+;22222222222222222222222222222222222222222222222222222222222222222222222
+
+
+;33333333333333333333333333333333333333333333333333333333333333333333333
+LABEL3: IPAGE=3
+
+ptilt=tilt
+PRINT,'PLOT REVERSED TILT? (YES=1):'
+READ,revtilt
+if (revtilt eq 1) then ptilt=180-tilt
+pres=moment(ptilt,sdev=sdp)
+
+ERASE
+!P.MULTI = [0,1,4,0]
+TL=' Scale='+string(tscale)+' Smooth='+string(smsize)
+
+PLOT, time, vmean, LINESTYLE = 0, TITLE='Vrms(vol)',YRANGE=[vmin,vmax]
+OPLOT, time, replicate(res17(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res17(0)-sd17,n),LINESTYLE = 1
+OPLOT, time, replicate(res17(0)+sd17,n),LINESTYLE = 1
+
+PLOT, time, bmean, LINESTYLE = 0, TITLE='Brms(vol)',YRANGE=[bmin,bmax]
+OPLOT, time, replicate(res12(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res12(0)-sd12,n),LINESTYLE = 1
+OPLOT, time, replicate(res12(0)+sd12,n),LINESTYLE = 1
+
+PLOT, time, dipole, LINESTYLE = 0, TITLE='Axial Dipole(rms)',YRANGE=[dmin,dmax]
+OPLOT, time, replicate(res14(0),n),LINESTYLE = 2
+OPLOT, time, replicate(res14(0)-sd14,n),LINESTYLE = 1
+OPLOT, time, replicate(res14(0)+sd14,n),LINESTYLE = 1
+
+PLOT, time, ptilt, LINESTYLE = 0, TITLE='Tilt',XTITLE='Time'+TL
+OPLOT, time, replicate(pres(0),n),LINESTYLE = 2
+OPLOT, time, replicate(pres(0)-sdp,n),LINESTYLE = 1
+OPLOT, time, replicate(pres(0)+sdp,n),LINESTYLE = 1
+
+GOTO, LABEL99
+;33333333333333333333333333333333333333333333333333333333333333333333333
+
+
+;4444444444444444444444444444444444444444444444444444444444444444444
+;COMPUTE STATISTICS
+LABEL4: IPAGE=4
+
+kemax=max(ke) & kemin=min(ke)
+memax=max(me) & memin=min(me)
+apmax=max(meaxipol) & apmin=min(meaxipol)
+atmax=max(meaxitor) & apmin=min(meaxitor)
+dmax=max(dipole) & dmin=min(dipole)
+admax=max(dipaxi) & admin=min(dipaxi)
+nutmax=max(nutop) & nutmin=min(nutop)
+nubmax=max(nubot) & nubmin=min(nubot)
+vmax=max(vmean) & vmin=min(vmean)
+bmax=max(bmean) & bmin=min(bmean)
+tilmax=max(tilt) & tilmin=min(tilt)
+lonmax=max(diplong) & lonmin=min(diplong)
+
+PRINT, ' '
+PRINT, 'Parameter Average Stnd Dev Max Min'
+res1=moment(ke,sdev=sd) & print, 'KE: ' ,res1(0),sd1,kemax,kemin
+res2=moment(kepol,sdev=sd) & print, 'KE(pol): ',res2(0),sd2
+res3=moment(me,sdev=sd) & print, 'ME: ', res3(0),sd3,memax,memin
+res4=moment(mepol,sdev=sd) & print, 'ME(pol): ', res4(0),sd4
+res5=moment(keaxtor,sdev=sd) & print, 'KE(axitor): ', res5(0),sd5
+res6=moment(keaxpol,sdev=sd) & print, 'KE(axipol): ', res6(0),sd6
+res8=moment(meaxipol,sdev=sd) & print, 'ME(axipol): ', res8(0),sd8,apmax,apmin
+res9=moment(meaxitor,sdev=sd) & print, 'ME(axitor): ', res9(0),sd9,atmax,atmin
+res10=moment(nutop,sdev=sd) & print, 'Nu(top): ', res10(0),sd10,nutmax,nutmin
+res11=moment(nubot,sdev=sd) & print, 'Nu(bot): ', res11(0),sd11,nubmax,nubmin
+res12=moment(bmean,sdev=sd) & print, 'B(ave): ', res12(0),sd12,bmax,bmin
+res13=moment(dipole,sdev=sd) & print, 'B(dip): ', res13(0),sd13,dmax,dmin
+res14=moment(dipaxi,sdev=sd) & print, 'B(axi): ', res14(0),sd14,admax,admin
+res15=moment(tilt,sdev=sd) & print, 'Tilt: ', res15(0),sd15,tilmax,tilmin
+res16=moment(diplong,sdev=sd16) & print, 'Lon(ave): ', res16(0),sd16,lonmax,lonmin
+res17=moment(vmean,sdev=sd17) & print, 'V(ave): ', res17(0),sd17,vmax,vmin
+
+PRINT, ' '
+PRINT,'Correlations'
+PRINT, 'KE vs ME: ', correlate(ke,me)
+PRINT, 'KE vs Nu(top): ', correlate(ke,nutop)
+PRINT, 'KE vs Nu(bot): ', correlate(ke,nubot)
+PRINT, 'ME vs Nu(top): ', correlate(me,nutop)
+PRINT, 'ME vs Nu(bot): ', correlate(me,nubot)
+PRINT, 'dipole vs Nu(top): ', correlate(dipole,nutop)
+PRINT, 'dipole vs Nu(bot): ', correlate(dipole,nubot)
+PRINT, 'dipole vs vmean: ', correlate(dipole,vmean)
+PRINT, 'dipole vs tilt: ', correlate(dipole,tilt)
+
+;POLE STATISTICS
+rad=pi*abs(cos(tilt))/180.
+phi=pi*diplong/180.
+xpole=rad*cos(phi) & ypole=rad*sin(phi)
+resx=moment(xpole,sdev=sdx)
+resy=moment(ypole,sdev=sdy)
+aver=sqrt(resx(0)^2 + resy(0)^2)
+aveincl=180*atan(aver,1)/pi
+avelong=180*atan(resy(0),resx(0))/pi
+rmstilt=res15(0)
+PRINT, ' '
+PRINT, 'POLE STASTICS'
+PRINT, ' Xave Xsd Yave Ysd Rave'
+PRINT, resx(0),sdx,resy(0),sdy,aver
+PRINT,'Dipole Long(ave) Tilt(rms,deg) Incl(ave,deg)'
+PRINT, avelong,rmstilt,aveincl
+PRINT,' '
+
+GOTO, LABEL99
+;4444444444444444444444444444444444444444444444444444444444444444444444
+
+
+;5555555555555555555555555555555555555555555555555555555555555555555
+LABEL5: IPAGE=5
+PRINT, 'CHANGE TIME SCALE AND SMOOTH'
+PRINT, ' '
+PRINT, 'TIMESCALE FACTOR=' , tscale
+PRINT, 'ENTER TIME SCALE FACTOR (def=1):' & read,tscale
+ time=tscale*time
+PRINT, ' '
+ tzro=time(0)
+PRINT, 'SHIFT INITIAL TIME TO 0? (yes=1):' & read,tshift
+ if tshift eq 1 then time=time-tzro
+PRINT, 'ENTER SMOOTHING FACTOR (odd; def=1):'
+READ,smsize
+
+ke=smooth(ke,smsize,/edge_truncate)
+kepol=smooth(kepol,smsize,/edge_truncate)
+me=smooth(me,smsize,/edge_truncate)
+mepol=smooth(mepol,smsize,/edge_truncate)
+keaxtor=smooth(keaxtor,smsize,/edge_truncate)
+keaxpol=smooth(keaxpol,smsize,/edge_truncate)
+meaxipol=smooth(meaxipol,smsize,/edge_truncate)
+meaxitor=smooth(meaxitor,smsize,/edge_truncate)
+nutop=smooth(nutop,smsize,/edge_truncate)
+nubot=smooth(nubot,smsize,/edge_truncate)
+bmean=smooth(bmean,smsize,/edge_truncate)
+dipole=smooth(dipole,smsize,/edge_truncate)
+dipaxi=smooth(dipaxi,smsize,/edge_truncate)
+tilt=smooth(tilt,smsize,/edge_truncate)
+vmean=smooth(vmean,smsize,/edge_truncate)
+
+PRINT, ' '
+PRINT,'NEW TIMESCALE FACTOR=',tscale
+PRINT,'NEW SMOOTHING FACTOR=',smsize
+
+GOTO, LABEL99
+;5555555555555555555555555555555555555555555555555555555555555555555
+
+;666666666666666666666666666666666666666666666666666666666666666666
+;SPECTRA
+LABEL6: IPAGE=6
+PRINT, 'Enter number of frequencies:'
+Read, fmax
+
+;REMOVE MEANS
+bzm=bmean-res12(0)
+dzm=dipole-res13(0)
+vzm=vmean-res16(0)
+nzm=(nutop-res10(0)+nubot-res11(0))/2
+
+;POWER SPECTRA
+bp=0.5*abs(fft(bzm,-1))^2
+dp=0.5*abs(fft(dzm,-1))^2
+vp=0.5*abs(fft(vzm,-1))^2
+np=0.5*abs(fft(nzm,-1))^2
+
+!P.MULTI = [0,2,2,0]
+;!P.POSITION = XY1
+plot, /ylog,/xlog,bp,title='AVE FIELD',ytitle='Power',$
+ xtitle='Frequency',xrange=[1,fmax]
+
+;!P.POSITION = XY2
+plot, /ylog,/xlog,dp,title='DIPOLE FIELD',ytitle='Power',$
+ xtitle='Frequency',xrange=[1,fmax]
+
+;!P.POSITION = XY3
+plot, /ylog,/xlog,vp,title='AVE VELOCITY',ytitle='Power',$
+ xtitle='Frequency',xrange=[1,fmax]
+
+;!P.POSITION = XY4
+plot, /ylog,/xlog,np,title='AVE NU',ytitle='Power',$
+ xtitle='Frequency',xrange=[1,fmax]
+
+GOTO, LABEL99
+;666666666666666666666666666666666666666666666666666666666666666666
+
+
+;77777777777777777777777777777777777777777777777777777777777777777
+LABEL7: IPAGE=7
+!P.MULTI = [0,2,2,0]
+
+LOADCT, 39 ;loads color table 39
+
+;CONVERT FROM TILT TO LATITUDE
+diplat=90-TILT
+
+;NORTH POLAR
+MAP_SET,90,0,0,/ORTHOGRAPHIC,/ADVANCE,/NOBORDER,/NOERASE,/ISOTROPIC,$
+ TITLE='DIPOLE AXIS - NORTH', LIMIT=[70,-180,90,180]
+MAP_GRID,LATDEL=5,LONDEL=45,LONLAB=72,LATLAB=45,GLINETHICK=1,LABEL=2
+PLOTS,diplong,diplat,PSYM=3,SYMSIZE=1,COLOR=250
+
+;SOUTH POLAR
+MAP_SET,-90,0,0,/ORTHOGRAPHIC,/ADVANCE,/NOBORDER,/NOERASE,/ISOTROPIC,$
+ TITLE='DIPOLE AXIS - SOUTH', LIMIT=[-90,-180,-70,180]
+MAP_GRID,LATDEL=5,LONDEL=45,LONLAB=-72,LATLAB=45,GLINETHICK=1,LABEL=2
+PLOTS,diplong,diplat,PSYM=3,SYMSIZE=1,COLOR=250
+
+;PRIMARY MERIDIAN
+MAP_SET,0,0,0,/ORTHOGRAPHIC,/ADVANCE,/NOBORDER,/NOERASE,/ISOTROPIC,TITLE='0 E'
+MAP_GRID,LATDEL=30,LONDEL=30,GLINETHICK=1,/HORIZON
+PLOTS,diplong,diplat,PSYM=3,SYMSIZE=1,COLOR=250
+
+;SECOND MERIDIAN
+MAP_SET,0,180,0,/ORTHOGRAPHIC,/ADVANCE,/NOBORDER,/NOERASE,/ISOTROPIC,TITLE='180 E'
+MAP_GRID,LATDEL=30,LONDEL=30,GLINETHICK=1,/HORIZON
+PLOTS,diplong,diplat,PSYM=3,SYMSIZE=1,COLOR=250
+
+
+GOTO, LABEL99
+;77777777777777777777777777777777777777777777777777777777777777777
+
+;11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11
+LABEL11: PRINT, FORMAT='("WINDOW SCALE FACTOR? CURRENT=",F7.3)',SCWINDOW
+ READ, SCWINDOW
+ XWIND=XWINDOW*SCWINDOW & YWIND=YWINDOW*SCWINDOW
+ CSZ=SCWINDOW*0.8 & CSB=1.12*SCWINDOW & CSS=0.60*SCWINDOW
+ GOTO, LABEL99
+;11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11 11
+
+;21 21 21 21 21 21 21 21 21 21 21 21 21 21 21 21 21
+LABEL21: &
+OUTFILE= ' '
+PRINT,'ENTER FULL .PS FILE NAME:'
+READ, OUTFILE
+SET_PLOT,'PS'
+DEVICE,FILENAME=OUTFILE,/COLOR
+DEVICE,XSIZE=XDIM,YSIZE=YDIM,XOFFSET=1.8,YOFFSET=2.0
+
+ CASE IPAGE OF
+ 1: GOTO, LABEL1
+ 2: GOTO, LABEL2
+ 3: GOTO, LABEL3
+ 4: GOTO, LABEL4
+ 5: GOTO, LABEL5
+ 6: GOTO, LABEL6
+ 7: GOTO, LABEL7
+ ELSE: GOTO, LABEL0
+ ENDCASE
+
+;21 21 21 21 21 21 21 21 21 21 21 21 21 21 21 21 21
+
+;99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99
+LABEL99: IF IOPTION EQ 21 THEN DEVICE,/CLOSE
+ !P.MULTI = 0
+ SET_PLOT, 'X'
+ LOADCT,LCOLOR
+ GOTO,LABEL0
+;99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99 99
+
+LABEL98: PRINT,'EXIT MAGTS'
+
+END
More information about the cig-commits
mailing list