* load <baseline> <ext> <cp1> <cp2> ...
 macro load 
   baseline=[1] ; shift
   ext=[1] ; shift
   while [#] .ne. 0 do
     cp=[1] ; shift
     ve/read en,nue[cp],numu[cp],nutau[cp] vec/numu-[baseline]-[cp]-[ext]
     ve/read en,anue[cp],anumu[cp],anutau[cp] vec/antinumu-[baseline]-[cp]-[ext]
     ve/read en,nnue[cp],nnumu[cp],nnutau[cp] vec/nomatter-numu-[baseline]-[cp]-[ext]
     ve/read en,nanue[cp],nanumu[cp],nanutau[cp] vec/nomatter-antinumu-[baseline]-[cp]-[ext]
     ve/copy nue[cp] v1		| work around lameness in sigma
     ve/copy anue[cp] v2
     sigma a[cp]=(v1-v2)/(v1+v2)
     ve/copy nnue[cp] v1	| work around lameness in sigma
     ve/copy nanue[cp] v2
     sigma na[cp]=(v1-v2)/(v1+v2)
   endwhile
   sigma en=en/1.0e9
 return
 
 macro numunue cp nue
   if [nue] = nue then
     title 'P([n]?[m]!"5#[n]?e!) (matter) CP phase='//[cp]
   else
     title 'P([n]?[m]!"5#[n]?e!) (anti-matter) CP phase='//[cp]
   endif
   null 0.5 10 0 .1
   set mscf 1

   set mtype 20
   set pmci 1
   ve/pl [nue]0%en ! s
   
   set mtype 24
   set pmci 2
   ve/pl [nue][cp]%en ! s
 return
 
 macro hist id vec
   s=$vdim([vec],1)
   1d [id] [vec] [s] 0.5 10
   put/cont [id] [vec]
 return 
 
 macro plot-prob max nu1 nu2 baseline dm2sol dm2atm s22t12 s22t23 s22t13
   ystart=[max]*0.95		| Y start value for multi lables
   ystep=[max]*0.05		| Y step value for multi labels
   xstart1=4.5
   xstart2=7.0
   exec cpa#hist 1 n[nu2]0
   exec cpa#hist 2 n[nu2]45
   exec cpa#hist 3 n[nu2]90
   exec cpa#hist 4 n[nu2]135
   exec cpa#hist 5 n[nu2]180
   exec cpa#hist 6 n[nu2]225
   exec cpa#hist 7 n[nu2]270
   exec cpa#hist 8 n[nu2]315

   exec neutrino#fancyname [nu1]
   name1=[@]
   exec neutrino#fancyname [nu2]
   name2=[@]

   title 'P('//[name1]//'!"5#'//[name2]//')'
   null 0.5 10 0 [max]
   exec multi#plot-dmod 0 1 1 2 2 3 3 4 4 5 5 6 6 7 7 8 8
   exec multi#line lin [xstart1] [ystart] [ystep] 1 1 '[d]=0' 2 2 '[d]=45' 3 3 '[d]=90' 4 4 '[d]=135' 5 5 '[d]=180' 6 6 '[d]=225' 7 7 '[d]=270' 8 8 '[d]=315'
   L='L= '//[baseline]//' km'
   m1='[D]m^2!?sol!= '//[dm2sol]//' eV^2!'
   m2='[D]m^2!?atm!= '//[dm2atm]//' eV^2!'
   s1='sin^2!(2[Q]?12!)= '//[s22t12]
   s2='sin^2!(2[Q]?23!)= '//[s22t23]
   s3='sin^2!(2[Q]?13!)= '//[s22t13]
   exec multi#itx lin [xstart2] [ystart] [ystep] 1 [L] 1 [m1] 1 [m2] 1 [s1] 1 [s2] 1 [s3]
   atitle 'E (GeV)' 'Probability'
   
 return

 macro plot-asym baseline dm2sol dm2atm s22t12 s22t23 s22t13
   start=0.0
   step=.1

   exec cpa#hist 1 a0
   exec cpa#hist 2 a45
   exec cpa#hist 3 a90
   exec cpa#hist 4 a180

   title 'CP Asymmetry'
   null 0.5 10 -1 1
   exec multi#plot-dmod 0 1 1 2 2 6 3 4 4
   exec multi#line lin 4.5 [start] [step] 1 1 '[d]=0' 2 2 '[d]=45' 6 3 '[d]=90' 4 4 '[d]=180'
   L='L= '//[baseline]//' km'
   m1='[D]m^2!?sol!= '//[dm2sol]//' eV^2!'
   m2='[D]m^2!?atm!= '//[dm2atm]//' eV^2!'
   s1='sin^2!(2[Q]?12!)= '//[s22t12]
   s2='sin^2!(2[Q]?23!)= '//[s22t23]
   s3='sin^2!(2[Q]?13!)= '//[s22t13]
   exec multi#itx lin 7.0 [start] [step] 1 [L] 1 [m1] 1 [m2] 1 [s1] 1 [s2] 1 [s3]
   atitle 'E (GeV)' 'Asymmetry'
   
   hi/del 1
   hi/del 2
   hi/del 3
   hi/del 4
 return
 
 macro plot-asym-diff baseline dm2sol dm2atm s22t12 s22t23 s22t13
   start=1.7
   step=0.1

   sigma a45d=a45-a0
   sigma a90d=a90-a0
   sigma a180d=a180-a0
   exec cpa#hist 2 a45d
   exec cpa#hist 3 a90d
   exec cpa#hist 4 a180d
   
   title 'Asymmetry differences w/ and w/out CP phase'
   null 0.5 10 -0.1 2
   exec multi#plot-dmod 1 2 2 6 3 4 4
   exec multi#line lin 4.5 [start] [step] 2 2 '[d]=45' 6 3 '[d]=90' 4 4 '[d]=180'
   L='L= '//[baseline]//' km'
   m1='[D]m^2!?sol!= '//[dm2sol]//' eV^2!'
   m2='[D]m^2!?atm!= '//[dm2atm]//' eV^2!'
   s1='sin^2!(2[Q]?12!)= '//[s22t12]
   s2='sin^2!(2[Q]?23!)= '//[s22t23]
   s3='sin^2!(2[Q]?13!)= '//[s22t13]
   exec multi#itx lin 7.0 [start] [step] 1 [L] 1 [m1] 1 [m2] 1 [s1] 1 [s2] 1 [s3]
   atitle 'E (GeV)' 'A([d])-A(0)'
   
   hi/del 2
   hi/del 3
   hi/del 4

 return
 
 macro do-baseline baseline dm2sol dm2atm s22t12 s22t23 s22t13
   hi/del *
   ve/del *
   basename=[s22t12]//'-'//[s22t23]//'-'//[s22t13]//'-'//[dm2sol]//'-'//[dm2atm]

   exec cpa#load [baseline] [basename]

   exec plotit plots/numu-nue-[baseline]-[basename].eps -113 cpa#plot-prob 0.2 numu nue [baseline] [dm2sol] [dm2atm] [s22t12] [s22t23] [s22t13]
   exec plotit plots/anumu-anue-[baseline]-[basename].eps -113 cpa#plot-prob 0.2 anumu anue [baseline] [dm2sol] [dm2atm] [s22t12] [s22t23] [s22t13]
   exec plotit plots/numu-numu-[baseline]-[basename].eps -113 cpa#plot-prob 1.0 numu numu [baseline] [dm2sol] [dm2atm] [s22t12] [s22t23] [s22t13]
   exec plotit plots/anumu-anumu-[baseline]-[basename].eps -113 cpa#plot-prob 1.0 anumu anumu [baseline] [dm2sol] [dm2atm] [s22t12] [s22t23] [s22t13]
   exec plotit plots/asym-[baseline]-[basename].eps -113 cpa#plot-asym [baseline] [dm2sol] [dm2atm] [s22t12] [s22t23] [s22t13]
   exec plotit plots/asym-diff-[baseline]-[basename].eps -113 cpa#plot-asym-diff [baseline] [dm2sol] [dm2atm] [s22t12] [s22t23] [s22t13]

 return
 
 macro test baseline dm2sol dm2atm s22t12 s22t23 s22t13
   hi/del *
   ve/del *
   basename=[s22t12]//'-'//[s22t23]//'-'//[s22t13]//'-'//[dm2sol]//'-'//[dm2atm]
   message [basename]
   exec cpa#load [baseline] [basename]
   message 'done loading'
   exec plotit none -113 cpa#plot-prob 0.2 numu nue [baseline] [dm2sol] [dm2atm] [s22t12] [s22t23] [s22t13]
 return
 

 macro all
   exec init 'seminar'
   exec do-baseline 2540 5.0e-5 3.0e-3 0.8 1.0 0.1
   exec do-baseline 2920 5.0e-5 3.0e-3 0.8 1.0 0.1
   exec do-baseline 1920 5.0e-5 3.0e-3 0.8 1.0 0.1
   exec do-baseline 1290 5.0e-5 3.0e-3 0.8 1.0 0.1
   exec do-baseline 1770 5.0e-5 3.0e-3 0.8 1.0 0.1
   exec do-baseline  735 5.0e-5 3.0e-3 0.8 1.0 0.1
   
   exec do-baseline 2540 5.0e-5 3.0e-3 0.8 1.0 0.04
   exec do-baseline 2540 5.0e-5 3.0e-3 0.8 1.0 0.2

 return