/Users/ros/RSAP/v6mac > /sw/bin/gfortran -c samover6.f
samover6.f: In function `samover':
samover6.f:625: error: could not split insn
(insn 5754 5753 5755 (set (reg:SI 0 r0 [2328])
(const_int 4288247288 [0xff9975f8])) 313 {*movsi_internal1} (nil)
(nil))
samover6.f:625: internal compiler error: in final_scan_insn, at final.c:2515
Please submit a full bug report,
with preprocessed source if appropriate.
source file follows:
c==========================================================================
c
subroutine samover (luin, device, ipass, numodf_files)
c
c Overlay numodf_files plots a.la. samstack 9 Jan 2004
c
use rsap_constants
use rsap_BLK2
use rsap_BLK3
use rsap_BLK4
use rsap_BLK5
use rsap_BLK6
use rsap_BLK7
use rsap_BLK8
use rsap_BLK12
use rsap_yaxislabels
use rsap_union
use rsap_viewport
c
implicit real*8 (a-h,o-z)
implicit integer*8 (i-n)
c
integer*8 numodf, nrecdo
integer*4 k ! 11 Mar 2005
c
integer*4 redpsc, bluepsc, greenpsc, red, blue, green
data redpsc,bluepsc,greenpsc, red,blue,green /255,255,255, 0,0,0/
c
real*8 yplot(nmax), timep(nmax)
character*16 ihead, ylabel, ylab(80)
character*8 ndum(10), device
character*80 ndum80
c
character*16 ihead2, ylabel2, ylab2(80)
character*8 ndum2(10)
c
cccccccc character*24 xaxlab, yaxlab ! 27 Jan 2004
character*48 xaxlab, yaxlab
c
character*80 version
data number/20/
character*24 fdate; external fdate
c
character*2 yaxlabrun, nrunid
c
character*3 typeplotupper(numplotypes)
character*3 typeplotlower(numplotypes)
character*3 typeplot, blank
data typeplotupper/'D ','DE ','DT ','DB ','DET',
& 'DEB','RT ','RB ','BOD','DTB','ETB','TMD','BMD','PD','INT',
& 'TOT','BOB','TMT','BMB','TPT','BPB','TOD','BOT','DOD','DMD',
& 'E ','EOD'/
data typeplotlower/'d ','de ','dt ','db ','det',
& 'deb','rt ','rb ','bod','dtb','etb','tmd','bmd','pd','int',
& 'tot','bob','tmt','bmb','tpt','bpb','tod','bot','dod','dmd',
& 'e ','eod'/
data blank/' '/
c------------------------------------------------
c
c knumrun(i) = 1,2 says plot from odfname, odf2name
c knumrun(i) = 12 says plot from both odf files
c
c kover > 1: read a group of nxp*nyp*kover cards
c kover = 1: read a group of nxp*nyp*kstack cards
400 i = 0
420 i = i + 1
cccccccc print "($,
print "(
&' Enter run, type, emin,emax, ymin,ymax (i2,a3,4f10) > ')"
if (kcmf.ne.0) then
read (lucmf,fmt='(a80)', err=450, end=500) ndum80
else
read ( luin,fmt='(a80)', err=450, end=500) ndum80
endif
if (kdebug.gt.1) print *, ' ndum80 =', ndum80
read (unit=ndum80, fmt='(a2)', err=450) nrunid
if (nrunid.eq.'q' .or. nrunid.eq.'Q' .or.
& nrunid.eq.blank) go to 500
c
if (nrunid.eq.'1 ' .or. nrunid.eq.'2 ' .or.
& nrunid.eq.'3 ' .or. nrunid.eq.'4 ' .or.
& nrunid.eq.'5 ' .or. nrunid.eq.'6 ' .or.
& nrunid.eq.'7 ' .or. nrunid.eq.'8 ' .or.
& nrunid.eq.'12') then
ktransrun(i) = 0
read (unit=ndum80, fmt='(i2,a3, 4f10.4)', err=450)
& knumrun(i), typeplot, xmn(i),xmx(i), ymn(i),ymx(i)
else
if (nrunid.ne.'a ' .and. nrunid.ne.'b ' .and.
& nrunid.ne.'c ' .and. nrunid.ne.'d ' .and.
& nrunid.ne.'e ' .and. nrunid.ne.'f ' .and.
& nrunid.ne.'g ' .and. nrunid.ne.'h ' .and.
& nrunid.ne.'ab') go to 450
ktransrun(i) = 1
if (nrunid.eq.'a ') knumrun(i) = 1
if (nrunid.eq.'b ') knumrun(i) = 2
if (nrunid.eq.'c ') knumrun(i) = 3
if (nrunid.eq.'d ') knumrun(i) = 4
if (nrunid.eq.'e ') knumrun(i) = 5
if (nrunid.eq.'f ') knumrun(i) = 6
if (nrunid.eq.'g ') knumrun(i) = 7
if (nrunid.eq.'h ') knumrun(i) = 8
if (nrunid.eq.'ab') knumrun(i) = 12
read (unit=ndum80, fmt='(2x, a3, 4f10.4)', err=450)
& typeplot, xmn(i),xmx(i), ymn(i),ymx(i)
endif
if (kdebug.gt.0) print
& "(' ------- i ktransrun knumrun typeplot xmn xmx ymn ymx'
& /8x,2i3,1x,i3,1x,a3, 4f10.4)",
& i,ktransrun(i),knumrun(i),typeplot, xmn(i),xmx(i), ymn(i),ymx(i)
go to 460
450 print *,' ***** format error. try again '
i = i - 1
go to 420
460 continue
c
if (knumrun(i).eq.99 .or.
& typeplot.eq.blank) then
go to 500
endif
c
if (knumrun(i).eq.0 ) knumrun(i) = 1
if (knumrun(i).eq.10) knumrun(i) = 1
if (knumrun(i).eq.20) knumrun(i) = 2
if (knumrun(i).lt.0) then
i = i - 1
go to 420
endif
if (knumrun(i).ne.1 .and. knumrun(i).ne.2 .and.
& knumrun(i).ne.3 .and. knumrun(i).ne.4 .and.
& knumrun(i).ne.5 .and. knumrun(i).ne.6 .and.
& knumrun(i).ne.7 .and. knumrun(i).ne.8 .and.
& knumrun(i).ne.12) then
print *,' ***** ignored illegal run number ', knumrun(i)
i = i - 1
go to 420
endif
if (knumrun(i).ne.12 .and.
& knumrun(i).gt.numodf_files) then
print *, ' ****** samover error *****'
print *,' numodf_files =',numodf_files
i = i - 1
go to 420
endif
numplot(i) = 1
do k=1,numplotypes
if (typeplot.eq.typeplotupper(k)) numplot(i) = k
if (typeplot.eq.typeplotlower(k)) numplot(i) = k
enddo
print "(i3, 2x, a3, 2x, i3, 2x, 4f10.4)",
& knumrun(i),typeplot, numplot(i),xmn(i),xmx(i), ymn(i),ymx(i)
c
c Have we read nxp*nyp*kstack cards ?
c
if (kover.eq.1 .and. i.lt.nxplots*nyplots*kstack) go to 420
c
c Have we read nxplots*nyplots*kover cards ?
if (kover.gt.1 .and. i.lt.nxplots*nyplots*kover) go to 420
i = i + 1
500 nplot = i - 1
c
if (nplot.le.0) then
if (kdebug > 0) print *, ' samover 500, nplot .le. 0'
cccccccc if (kcmf.ne.0) kcmf = 0
return
endif
c------------------------------------------------------------
if (kover.eq.1 .and. kstack.gt.1) then
ixlimit = 1
do i=2,nplot
xmn(i) = xmn(ixlimit)
xmx(i) = xmx(ixlimit)
if (nxplots.lt.2) cycle
if (i.eq.kstack .or. i.eq.2*kstack .or. i.eq.3*kstack) then
ixlimit = i + 1
print "('x limits set to plot',i3,' limits for plots',i3,' -',i4)"
& , i+1-kstack, i+2-kstack, i
endif
if (ixlimit.gt.nplot) ixlimit = nplot
enddo
if (nxplots.lt.2) print *,
&' x axis limits set to plot 1 limits for plots 2 -', nplot
endif
c------------------------------------------------------------
c------------------------------------------------------------
c ! start of device_pass loop
do k=1,kdevice_pass
select case (k)
case (2)
device = 'psc'
if (kdebug > 0) print *,
& ' ---------- plot output to file rsap.psc'
call plsfnam ('rsap.psc')
case (3)
device = 'ps'
if (kdebug > 0) print *,
& ' ---------- plot output to file rsap.ps'
call plsfnam ('rsap.ps')
case (4)
device = 'plmeta'
if (kdebug > 0) print *,
& ' ---------- plot output to file rsap.meta'
call plsfnam ('rsap.meta')
end select
c initialize PLPLOT
if (device.eq.'psc') call plscolbg (redpsc,bluepsc,greenpsc)
if (device.ne.'psc') call plscolbg (red,blue,green)
if (kfont.eq.2) call plfontld (kfont)
call plstart (device, nxplots, nyplots)
if (kfont.eq.2) call plfont (kfont)
c------------------------------------------------------------
c set x & y axis limits
if (kover > 1) call setaxlim
c----------------------------------------------- Start of loop on frames
name2(1) = 'SAMMY2.'
name2(2) = 'testxxxx'
c
xaxlab = 'E (keV)'
if (kplev.gt.0) xaxlab = 'E (eV)'
if (kuserxaxis.gt.0) xaxlab = userxaxislab
koveraxset = 0
do 888 i=1,nplot
kerrbar = 0
kcurves = 2
kdedet = 0
kaverage = kavgodf(i)
cccccccc if (kover .eq. 1) then
cccccccc xmin=xmn(i) ; xmax=xmx(i)
cccccccc ymin=ymn(i) ; ymax=ymx(i)
cccccccc endif
koveraxset = koveraxset + 1
if (kover .eq. 1 .or. koveraxset > kover) then
koveraxset = 0
xmin=xmn(i) ; xmax=xmx(i)
ymin=ymn(i) ; ymax=ymx(i)
endif
numrun = knumrun(i)
numodf = numrun
if (numodf.lt.1 .or. numodf.gt.nodfmax) numodf = 1
ktrans = ktransrun(i)
if (numrun.eq.12) then
nrecdo= nodfenergies(1)
else
nrecdo = nodfenergies(numrun)
endif
c knumrun(i) = 1,2 says plot from odfname, odf2name
c knumrun(i) = 12 or 21 says plot from both odf files
c.........
cccccccc if (numplot(i).lt.20) then
if (numplot(i).lt.28) then
cc numrun = 1 or 12 says load big into yplot
c numrun = 2 says load big2 into yplot
call fillup (timep,yplot,2)
if (ktrans.le.0) then
yaxlab = yaxislabel(numplot(i))
if (kaverage.gt.1)
& yaxlab = yaxislabelavg(numplot(i))
else
yaxlab = transyaxislabel(numplot(i))
if (kaverage.gt.1)
& yaxlab = transyaxislabelavg(numplot(i))
endif
write (unit=yaxlabrun, fmt='(i2)') numrun
kkk = kgetlen (yaxlab)
if (kkk.gt.22) kkk = 22
if (numplot(i) < 16 .or. numplot(i) .eq. 23)
& yaxlab(kkk+1:kkk+2) = yaxlabrun(1:2)
if (kuseryaxis.gt.0) yaxlab = useryaxislab(numodf)
if (kdebug > 0) print "('kkk, yaxlab:', i4,2x, a24)",
& kkk, yaxlab
endif
c
numcase = numplot(i)
c
if (numcase > 40 .and. numcase < 51) kcurves = 1
select case (numcase)
c********************************************************* Data
case (1)
kcurves = 1
if (numrun.eq.12) then
call fillup2(timep,temp1,2)
call fillup (timep,yplot,2)
endif
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************************* Data + Err
case (2)
kcurves = 1
kerrbar = 1
if (numrun.eq.12) then
call fillup2(timep,temp1,2)
call fillup (timep,yplot,2)
endif
call fillup (timep, yerr,3)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************************* Data + Th1
case (3)
kdedet = 1
if (numrun.eq.12) call fillup2(timep,temp1,4)
call fillup (timep, temp,4)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************************* Data +Bayes
case (4)
if (numrun.eq.12) call fillup2(timep,temp1,5)
call fillup (timep, temp,5)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************************** Data + Err + Th
case (5)
kdedet = 1
kerrbar = 1
if (numrun.eq.12) call fillup2(timep,temp1, 4)
call fillup (timep, yerr, 3)
call fillup (timep, temp, 4)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************************** Data + Err + Bayes
case (6)
kerrbar = 1
if (numrun.eq.12) call fillup2(timep,temp1,5)
call fillup (timep, yerr,3)
call fillup (timep, temp,5)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c******************************************* Residuals ((TH-Exp)/Err1)
case (7)
kcurves = 3
call fillup (timep, yerr, 3)
call fillup (timep, temp, 4)
do j=1,nrecdo
if (yerr(j) .gt. 0.0) then
yplot(j) = (temp(j)-yplot(j))/yerr(j)
else
yplot(j) = 0.0
endif
yerr(j) = 0.0
enddo
do j=1,nrecdo
temp(j) = yplot(j)
yplot(j) = 0.0
enddo
call pdover (timep,yplot,nrecdo, xaxlab, yaxlab)
c**************************************** Residuals ((Bayes-Exp)/Err1)
case (8)
kcurves = 3
call fillup (timep, yerr, 3)
call fillup (timep, temp, 5)
do j=1,nrecdo
if (yerr(j) .gt. 0.0) then
yplot(j) = (temp(j)-yplot(j))/yerr(j)
else
yplot(j) = 0.0
endif
yerr(j) = 0.0
enddo
cccccccc do j=1,nrecdo
cccccccc temp(j) = yplot(j)
cccccccc yplot(j) = 0.0
cccccccc enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************************ Bayes / Exp "BOD"
c************************************************ Theory / Exp "TOD"
case (9, 22)
kcurves = 1
jthbay = 5
if (numplot(i).eq.22) jthbay = 4
if (numrun.eq.12) then
call fillup2(timep,yplot,2) !load Data2
call fillup2(timep,temp1,jthbay) !load Bayes2
c ! or Theory2
do j=1,nrecdo
if (yplot(j) .gt. 0.0) then
temp1(j) =100.*(temp1(j)/yplot(j))
else
temp1(j) = 0.0
endif
enddo
call fillup (timep,yplot,2) ! load Data1
endif
c
cccccccc kcurves = 2
call fillup (timep, temp, jthbay) ! load Bayes1 or Theory1
do j=1,nrecdo
if (yplot(j) .gt. 0.0) then
yplot(j) =temp(j)/yplot(j) - one
else
yplot(j) = 0.0
endif
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************************** Data1 + Th + Bayes
case (10)
kcurves = 3
call fillup (timep, temp, 5)
call fillup (timep, temp1,4)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************** Data + Err + Th + Bayes
case (11)
kcurves = 3
kerrbar = 1
call fillup (timep, yerr, 3)
call fillup (timep, temp, 5)
call fillup (timep, temp1,4)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************************ (TH-Exp) "TMD"
case (12)
kcurves = 1
call fillup (timep, temp, 4)
do j=1,nrecdo
if (kdebug.gt.0) print *, yplot(j), temp(j)
yplot(j) = temp(j)-yplot(j)
if (kdebug.gt.0) print *, yplot(j)
enddo
call pdover (timep,yplot,nrecdo, xaxlab, yaxlab)
c************************************************ (Bayes-Exp) "BMD"
case (13)
kcurves = 1
call fillup (timep, temp,5)
do j=1,nrecdo
if (kdebug.gt.0) print *, yplot(j), temp(j)
yplot(j) = temp(j)-yplot(j)
if (kdebug.gt.0) print *, yplot(j)
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************************ (Bayes/Theory) "BOT"
case (23)
kcurves = 1
call fillup (timep, temp,5)
call fillup (timep, temp1,4)
do j=1,nrecdo
if (kdebug.gt.0) print *, temp1(j), temp(j)
yplot(j) = 0.0
if (temp1(j) > 0.0) yplot(j) = temp(j)/temp1(j)
if (kdebug.gt.0) print *, yplot(j)
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c*********************************************** (dsig/d(sqrt(E)) "PD"
case (14)
kcurves = 2
defac = xmin
c so pdo will autoscale x axis limits
xmin = 0.0
print *, 'dsig/d(sqrt(E). defac =', defac
do j=1,nrecdo
dsig = defac * (big(j,4,numodf+1) - big(j,4,numodf))
yplot(j) = dsig
if (kdebug.gt.0)
& print *, big(j,4,numodf+1), big(j,4,numodf), yplot(j)
enddo
c load SAMMY deriv. into temp
c write rsap.deriv
call pdsam (defac, timep, temp ,4)
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************** integrate (sig/E) "INT"
case (15)
kcurves = 2
call int_sig_over_e (numodf, nrecdo, yplot, timep)
c
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************** TH2 / TH1 "TOT", "TPT"
case (16, 20)
kcurves = 1
if (numrun.ne.12) then
print *, ' **** Error, numrun .NE. 12'
go to 888
endif
print *, 'nrecdo = ', nrecdo
if (kdebug.gt.1) write (33,*),
&' E TH1 TH2 Ratio Ratio-1'
call fillup (timep, temp1,4) ! load TH1
call fillup2(timep, yplot,4) ! load TH2
do j=1,nrecdo
yratio = zero
if (temp1(j) > zero) yratio = yplot(j)/temp1(j)
if (kdebug > 1) write (33, '(1p3e15.6, 0pf12.7, 1pe12.3)')
& timep(j), yplot(j),temp1(j), yratio, yratio-1.d0
yplot(j) = yratio
c TH2/TH1 - 1
if (numplot(i).eq.20) yplot(j) = yratio - one
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************* Bayes2 / Bayes1 "BOB", "BPB"
case (17, 21)
kcurves = 1
if (numrun.ne.12) then
print *, ' **** Error, numrun .NE. 12'
go to 888
endif
call fillup (timep,temp1,5) ! load Bayes1
call fillup2(timep,yplot,5) ! load Bayes2
do j=1,nrecdo
if (temp1(j) > zero) then
yplot(j) = yplot(j)/temp1(j)
c Bayes2/Bayes1 - 1
if (numplot(i).eq.21) yplot(j) = yplot(j) - one
else
yplot(j) = zero
endif
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c************************************************* TH2 - TH1 "TMT"
case (18)
kcurves = 1
if (numrun.ne.12) then
print *, ' **** Error, numrun .NE. 12'
go to 888
endif
call fillup (timep, temp1,4) ! load TH1
call fillup2(timep, yplot,4) ! load TH2
c
if (kdebug > 0) then
call delopen (9, 'rsap.tmt')
call system ('pwd > rsap.tmt')
read (9, '(a24)') idummy
write (9, '(a24)') fdate()
write (9, '(4a8)') (name(j), j=1,4)
write (9,*) 'j ANORM ODF File'
do j=1,2
write (9, '(i2, f8.4, 1x, a72)')
& j, anormodf(j), odfnames(j)
enddo
write (9,*) 'NRM =', anormodf(3)
write (9,*),
&' E TH1 TH2 dT=TH2-TH1 NRM*dT'
endif
c
do j=1,nrecdo
ydiff = yplot(j) - temp1(j)
if (kdebug > 0) write (9, '(f12.5, 1p2e15.6, 1p2e13.4)')
& timep(j), temp1(j), yplot(j), ydiff, anormodf(3)*ydiff
yplot(j) = ydiff
enddo
c
if (kdebug > 0) then
close (9)
print * , ' Th1, Th2, differences written to rsap.tmt'
endif
c
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c*********************************************** Bayes2 - Bayes1 "BMB"
case (19)
kcurves = 1
if (numrun.ne.12) then
print *, ' **** Error, numrun .NE. 12'
go to 888
endif
call fillup (timep,temp1,5) ! load Bayes1
call fillup2(timep,yplot,5) ! load Bayes2
do j=1,nrecdo
yplot(j) = yplot(j) - temp1(j)
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************** Data2 / Data1 "DOD"
case (24)
kcurves = 1
if (numrun.ne.12) then
print *, ' **** Error, numrun .NE. 12'
go to 888
endif
print *, 'nrecdo = ', nrecdo
call fillup (timep, temp1,2) ! load Data1
call fillup2(timep, yplot,2) ! load Data2
c
if (kdebug.gt.1) write (33,*),
&' E Data1 Data2 Ratio Ratio-1'
do j=1,nrecdo
yratio = zero
if (temp1(j) > zero) yratio = yplot(j)/temp1(j)
if (kdebug > 1) write (33, '(1p3e15.6, 0pf12.7, 1pe12.3)')
& timep(j), yplot(j),temp1(j), yratio, yratio-1.d0
yplot(j) = yratio
c Data2/Data1 - 1
cccccccc if (numplot(i).eq.26) yplot(j) = yratio - one
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c*********************************************** Data2 - Data1 "DMD"
case (25)
kcurves = 1
if (numrun.ne.12) then
print *, ' **** Error, numrun .NE. 12'
go to 888
endif
call fillup (timep,temp1,2) ! load Data1
call fillup2(timep,yplot,2) ! load Data2
do j=1,nrecdo
yplot(j) = yplot(j) - temp1(j)
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c*********************************************** Uncertainty "E"
case (26)
kcurves = 1
call fillup (timep,yplot,3) ! load Uncertainty
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c**************************************** Uncertainty / Data "EOD"
case (27)
kcurves = 1
call fillup (timep,temp1,2) ! load Data1
call fillup (timep,yplot,3) ! load Uncertainty
do j=1,nrecdo
if (dabs(temp1(j)) .le. zero) then
yplot(j) = zero
cycle
else
yplot(j) = yplot(j) / temp1(j)
endif
enddo
call pdover (timep,yplot,nrecdo,xaxlab, yaxlab)
c********************************************************* s1/2
case (41)
call fillup (timep, yplot,6)
call pdover (timep,yplot,nrecdo,xaxlab, 's1/2')
c********************************************************* p1/2
case (42)
call fillup (timep, yplot,7)
call pdover (timep,yplot,nrecdo,xaxlab, 'p1/2')
c********************************************************* p3/2
case (43)
call fillup (timep, yplot,8)
call pdover (timep,yplot,nrecdo,xaxlab, 'p3/2')
c********************************************************* d3/2
case (44)
call fillup (timep, yplot,9)
call pdover (timep,yplot,nrecdo,xaxlab, 'd3/2')
c********************************************************* d5/2
case (45)
call fillup (timep, yplot,10)
call pdover (timep,yplot,nrecdo,xaxlab, 'd5/2')
c********************************************************* f5/2
case (46)
call fillup (timep, yplot,11)
call pdover (timep,yplot,nrecdo,xaxlab, 'f5/2')
c********************************************************* f7/2
case (47)
call fillup (timep, yplot,12)
call pdover (timep,yplot,nrecdo,xaxlab, 'f7/2')
c********************************************************* g7/2
case (48)
call fillup (timep, yplot,13)
call pdover (timep,yplot,nrecdo,xaxlab, 'g7/2')
end select
c*********************************************************
888 continue
c---------------------------------------------------------------------
c Fig caption for nyplots > 1
if (kcaption > 0 .and. nxplots > 1) then
call caption_write(device)
endif
call plend()
if (kdebug.gt.0) print *, ' samover call plend():'
ncall = 0 ! 27 Nov 2000
if (k == 3) device = 'x'
enddo ! end of device_pass loop
if (nplot < nxplots*nyplots*kstack) return
c
go to 400 ! back to read another group of kover cards
end
c*****************************************************************
subroutine setaxlim
c
use rsap_constants
use rsap_BLK8
use rsap_BLK12
c
implicit real*8 (a-h,o-z)
implicit integer*8 (i-n)
c
xmin=xmn(1) ; xmax=xmx(1)
ymin=ymn(1) ; ymax=ymx(1)
do k=2,kover
if (xmn(k).ne.xmx(k)) then
if (xmn(k).lt.xmin) xmin = xmn(k)
if (xmx(k).gt.xmax) xmax = xmx(k)
endif
if (ymin.eq.ymax) then
if (ymn(k) < ymin) ymin = ymn(k)
if (ymx(k) > ymax) ymax = ymx(k)
endif
enddo
if(kdebug > 0) then
print *,'setaxlim: xmin xmax ymin ymax'
print "('setaxlim:', 4f16.5)", xmin,xmax, ymin,ymax
endif
c
if (xmin.eq.xmax) then
call xhilover
if(kdebug > 0)
& print "('setaxlim:', 4f16.5)", xmin,xmax, ymin,ymax
endif
if (ymin.eq.ymax) then
c Assume all overlay plots are trans or total
ktrans = ktransrun(1)
call yhilover (numplot)
if(kdebug > 0)
& print "('setaxlim:', 4f16.5)", xmin,xmax, ymin,ymax
endif
c-----
if (kyaxlog.le.0) then
ydiff = ymax - ymin
if (dabs(ydiff).lt.1.d-5) ydiff = 1.d-5
ymin = ymin - 0.02 * abs(ydiff)
ymax = ymax + 0.02 * abs(ydiff)
endif
if(kdebug.gt.0) print *,'setaxlim: ymin, ymax =', ymin, ymax
c
if (kyaxlog.gt.0) then
if (ymax.gt.0.0) then
ymax = dlog10(ymax)
else
ymax = 1.0
endif
if (ymin.gt.0.0) then
ymin = dlog10(ymin)
else
ymin = -3.
endif
endif
c------------------------------------------------------------------------
if (kxaxlog.gt.0) then
if (xmax.gt.0.0) then
xmax = dlog10(xmax)
else
xmax = 1.0
endif
if (xmin.gt.0.0) then
xmin = dlog10(xmin)
else
xmin = -3.
endif
endif
c------------------------------------------------------------------------
if (kdebug.gt.0) then
print *, 'setaxlim: kxaxlog xmin xmax'
print '(i15, 2f16.5)', kxaxlog, xmin, xmax
print *, 'setaxlim: kyaxlog ymin ymax'
print '(i15, 2f16.5)', kyaxlog, ymin, ymax
endif
return
end
c*****************************************************************
subroutine xhilover
c
use rsap_constants
use rsap_BLK2
use rsap_BLK12
c
implicit real*8 (a-h,o-z)
implicit integer*8 (i-n)
c
energyscale = one
cccccccc if (kplev.gt.0) energyscale = 1000.d0
if (kplev.eq.1) energyscale = 1000.d0
c
if (kover.lt.1 .or. kover.gt.nodfmax) then
print *, ' xhilover: kover OUT OF RANGE 1 -', nodfmax
return
endif
xmin = 1.d9
xmax = zero
do k=1,kover
nrecord = nodfenergies(k)
do i=1,nrecord
if (big(i,1, k) .lt. xmin) xmin = big(i,1, k)
if (big(i,1, k) .gt. xmax) xmax = big(i,1, k)
enddo
enddo
xmin = energyscale * xmin
xmax = energyscale * xmax
return
end
c*****************************************************************
subroutine yhilover (numplot)
c
use rsap_constants
use rsap_BLK2
use rsap_BLK12
c
implicit real*8 (a-h,o-z)
implicit integer*8 (i-n)
c
integer*8 numplot(1)
c
energyscale = one
cccccccc if (kplev.gt.0) energyscale = 1000.d0
if (kplev.eq.1) energyscale = 1000.d0
c
if (kover.lt.1 .or. kover.gt.nodfmax) then
print *, ' yhilover: kover OUT OF RANGE 1 -', nodfmax
return
endif
c total, TRANS for ktrans = 0,1
numfill = 2
if (kdebug > 0) print *, 'yhilover: ktrans =', ktrans
if (ktrans.gt.0) numfill = 2 + 4
c
ymin = 1.d+9
ymax = zero
j = numfill
jth = numfill + 2
jb = numfill + 3
c
if (kdebug > 0) print *, 'yhilover: xmin, xmax =', xmin, xmax
if (kdebug.gt.0) print *,
&'yhilover: k numplot i1 i2 ylo yhi ymin',
&' ymax anormodf'
do k=1,kover ! start loop on kover
ylo = 1.d+9
yhi = zero
nrecord = nodfenergies(k)
if (kdebug > 0) then
print *, 'yhilover: k, nrecord =', k,nrecord
print *,'yhilover: k, big(1,1, k) =', k, big(1,1, k)
print *,'yhilover: k, big(nrecord,1, k) =', k,big(nrecord,1, k)
endif
if (big(1,1, k) .gt. xmax .or.
& big(nrecord,1, k) .lt. xmin) cycle
if (kdebug > 1) print *, 'yhilover: k, nrecord =', k,nrecord
do i=1,nrecord
i1 = i
if (big(i,1, k) .gt. xmin) exit
enddo
if (kdebug > 1) print *, 'yhilover: k, i1 =', k,i1
do i=i1,nrecord
i2 = i
if (big(i,1, k) .gt. xmax) then
i2=i-1
exit
endif
enddo
if (kdebug > 1) print *, 'yhilover: k, i1,i2 =', k,i1,i2
if (i2.lt.i1) i2=i1
if (kdebug > 1) print *, 'yhilover: k, i1,i2 =', k,i1,i2
do i=i1,i2
if (big(i,j, k) .lt. ylo) ylo = big(i,j, k)
if (big(i,j, k) .gt. yhi) yhi = big(i,j, k)
enddo
if (kdebug > 0) print *, 'yhilover: ylo, yhi =', ylo, yhi
if (numplot(k). lt.3) go to 60
c Theory and Bayes
j1 = jth
j2 = jb
c Theory
if (numplot(k).eq.3 .or. numplot(k).eq.5) j2=jth
c Bayes
if (numplot(k).eq.4 .or. numplot(k).eq.6) j1=jb
if (kdebug.gt.1) print *, ' yhilover: j1, j2 =', j1, j2
c
do jj=j1,j2
do i=i1,i2
if (big(i,jj,k) .lt. ylo) ylo = big(i,jj,k)
if (big(i,jj,k) .gt. yhi) yhi = big(i,jj,k)
enddo
enddo
60 continue
ylo = anormodf(k) * ylo
yhi = anormodf(k) * yhi
if (ylo.lt.ymin) ymin = ylo
if (yhi.gt.ymax) ymax = yhi
if (kdebug.gt.0) print "(' yhilover:', 4i4, 4f12.5, f12.3)",
& k, numplot(k),i1,i2, ylo,yhi, ymin,ymax, anormodf(k)
enddo
c ! end loop on kover
return
end
c==========================================================================
c
subroutine pdover (x,y,npts, ix,iy)
c REWRITTEN to use PLPLOT
c
c* this subroutine draws the graph of y vs x c
c* x,y = arrays of x,y values
c* npts = number of points passed to pdo
c ix, iy = label for x, y axis
c
use rsap_constants
use rsap_pdo_constants
use rsap_BLK2
use rsap_BLK3
use rsap_BLK4
use rsap_BLK5
use rsap_BLK6
use rsap_BLK7
use rsap_BLK12
use rsap_BLK17
use rsap_color
use rsap_BLK20
use rsap_sgv
use rsap_viewport
c
implicit real*8 (a-h,o-z)
implicit integer*8 (i-n)
integer*4 indx(ngrps*nres)
real*4 ersave(ngrps*nres), erorder(ngrps*nres)
c
integer*4 kpts, kpts2, kup
integer*4 ncallann, karrow
c
real*4 bias, peak(nemax), pkarea(nemax), pkwidth(nemax)
real*4 gammaneut(nemax)
c
character*16 ylabel
character*18 kekevtext
character*4 kavgtext
character*(*) ix,iy
c
character*32 khead1, khead2
character*8 ihead, ihead2, blank
equivalence (khead2 (1:8), ihead2(1:8))
c
real*8 x(1),y(1)
real*4 xmin4,xmax4, ymin4,ymax4, x4(nmax),y4(nmax)
real*4 yeb4lo(nmax), yeb4hi(nmax)
cccccccc real*4 xvmin, xvmax, yvmin,yvmax, xtick, ytick, zero4, chsiz
real*4 xtick, ytick, zero4, chsiz
real*4 yone, y4temp(nmax), y4temp1(nmax), y4err(nmax)
real*4 xj1, xj2, yj1, yj2, samerout, ysam, ysamlo, reswidth
real*4 x4his(nmax), y4his(nmax)
real*4 chisqn
c
character*24 fdate; external fdate
c
data nslbw, zero4, yone /0, 0.0, 1.0/
c
integer*4 klength, klength2, numpeaks
integer*4 kcol4b, kcol4th, kcol4dat, kptscyanline, kcol4box
character*13 iytest
c 13 Nov 2000
integer*4 :: kheadlength = 32
call getheadstring (kheadlength, khead1)
c
if (numrun.ge.1 .and. numrun.le.nodfmax) then
odfname = odfnames(numrun)
kpoints = ksymodf(numrun)
endif
klength = kgetlen (odfname)
klength2 = kgetlen (odf2name)
if (kdebug.gt.0) then
print *,' pdover:----------------------------------------------'
if (kdebug.gt.1) print "(' pdover: klength, klength2 = ', 2i4)",
& klength, klength2
print "(' pdover: klength, odfname, klength2, odf2name :'
& /i4,a72 /i4,a72)", klength, odfname, klength2, odf2name
print "(' pdover: numrun npts xmin xmax ymin ymax:'
& /4x, 2i6, 4f12.3)", numrun, npts, xmin,xmax, ymin,ymax
print "(' pdover: x(1), x(npts), kpoints:'/f12.3,f9.3,i4)",
& x(1), x(npts), kpoints
endif
c
ihead2 = name2(1)
c
ncall = ncall + 1
ncallann = ncall
c-----
if (kover.eq.1) call setpdoaxlim (npts, x, y)
c------------------------------------------------------------------------
c kpts = number of points plotted (xmin < x(n) < xmax)
c
c kplo = first data point plotted xmin < x(kplo)
c kphi = last data point plotted xmax <
x(kphi)
c
kplo = 1
kphi = npts
xmintest = xmin
xmaxtest = xmax
if (kxaxlog.gt.0) then
xmintest = 10.d0**xmin
xmaxtest = 10.d0**xmax
endif
if (kxaxistime > 0) then
do n=1,npts
if (x(n).lt.xmintest) then
kphi = n - 1
exit
endif
if (x(n).gt.xmaxtest) kplo = n + 1
enddo
else
do n=1,npts
if (x(n).gt.xmaxtest) then
kphi = n -1
exit
endif
if (x(n).lt.xmintest) kplo = n + 1
enddo
endif
if (kplo < 1) kplo = 1
if (kphi > npts) kphi = npts
kpts = kphi - kplo + 1
if (kdebug.gt.0) then
print *,' pdover: kpts kplo kphi xlo xhi xmin xmax ymin ymax'
print "(' pdover:',3i5,1x, 6(f10.3,1x))",
& kpts,kplo,kphi, x(kplo),x(kphi), xmin,xmax, ymin,ymax
endif
c------------------------------------------------------------------------
c------------------------------------------------------------------------
if (kover.eq.1 .and. kbias.gt.0) then
c
c Find peaks for (xmin < x(n) < xmax)
bias = float(kbias)
k1find = kplo
k2find = kphi
c Compute kfwhm from fwhm (eV)
c
kfwhm = fwhm * float(kpts-1) / (x(kphi)-x(kplo))
if (kfwhm.lt.3) kfwhm = 3
if (kdebug > 0)
& print "('PKFIND: FWHM =', i5, ' data points')", kfwhm
c
call pkfind (k1find, k2find, kfwhm, bias, chisqn,
& peak, pkarea, pkwidth, gammaneut, numpeaks)
c save chisq for plot annotation
if (chisqn > 0.0)
& write (unit=kchisqn, fmt="('#gx#fr#u2#d/N =', f7.2)") chisqn
ichisqn = 1
c
if (kdebug > 0) print
& "('PDO: k1find k2find fwhm kfwhm kbias numpeaks'
& /'PDO: ',2i6, f12.5, 3i8)",
& k1find,k2find, fwhm, kfwhm, kbias, numpeaks
c
if (numpeaks > 0) call writepeaksamin (
& peak, pkarea, pkwidth, gammaneut, numpeaks)
c
endif
c------------------------------------------------------------------------
c------------------------------------------------------------------------
c draw graph using PLPLOT
c Set up the viewport and window using PLENV. The axes are
c scaled separately (just = 0), and we just draw a labelled
c box (axis = 0).
c
kcolax = jcolodfaxes(numrun) ! 9 Jan 2004
kcoltext = jcolodftext(numrun) ! 9 Jan 2004
c
call plcol(kcolax) ; call pllsty(1)
call plwid(klinewid)
c character size
chsiz = 0.5
cccccccc if (nxplots*nyplots.eq.1 .and. kstack.lt.3) chsiz = 0.7
if (nxplots*nyplots.le.3 .and. kstack.lt.3) chsiz = 0.65
if (nxplots*nyplots.le.2 .and. kstack.lt.3) chsiz = 0.7
chsiz = chsiz * axlabelscale
call plschr (zero4, chsiz)
xmin4 = xmin ; xmax4 = xmax
ymin4 = ymin ; ymax4 = ymax
c---------------------------------------------------
if (kstack.gt.1) then
if (ncall.eq.1) then
if (kdebug > 0) print
& "(' before pladv - nstack, kstack, ncall =', 3i4)",
& nstack, kstack, ncall
if (numrun.ne.92) call pladv(0)
nstack = 0
call set_viewport
endif
nstack = nstack + 1
nslbw = nslbw + 1
if (nstack.gt.kstack) then
nstack = 1
if (kdebug > 0) print
& "(' before pladv - nstack, kstack, ncall =', 3i4)",
& nstack, kstack, ncall
if (numrun.ne.92) call pladv(0)
call set_viewport
endif
c 3 Dec 2004
if (kstack == 2 .and. kdebug == 1) then
if (nstack == 1) then
deltay = 0.5 * deltay
yvmin = yvmin + deltay
endif
if (nstack == 2) then
yvmin = yvmin - 2.*deltay
deltay = 3. * deltay
endif
endif
c 3 Dec 2004
c 6 Dec 2004
if (yaxscale < 1.0 .and. yaxscale > 0.0) then
if (nstack == 1) then
yvmin = yvmin + (1.-yaxscale)*deltay
deltay = yaxscale * deltay
endif
if (nstack == 2) then
deltay = deltay * (2.-yaxscale)/yaxscale
yvmin = yvmin - deltay*(1.-yaxscale)/(1. -0.5*yaxscale)
endif
endif
c 6 Dec 2004
yvmin = yvmin + deltay
yvmax = yvmin + 0.94 * deltay
if (kdebug.gt.0)
& print "(' *** ncall,kstack,nstack, yvmin,yvmax =',3i5, 2f10.5)",
& ncall, kstack, nstack, yvmin, yvmax
call plvpor (xvmin, xvmax, yvmin,yvmax)
call plwind (xmin4, xmax4, ymin4, ymax4)
if (kyaxlog.gt.0) then
if (kxaxlog.le.0) then
if (nstack.eq.1) call plbox('bcnst',0,0,'bcnstvl',0,0)
if (nstack.gt.1) call plbox('bcst', 0,0,'bcnstvl',0,0)
else
if (nstack.eq.1) call plbox('bcnstl',0,0,'bcnstvl',0,0)
if (nstack.gt.1) call plbox('bcstl', 0,0,'bcnstvl',0,0)
endif
else
if (kxaxlog.le.0) then
if (nstack.eq.1) call plbox('bcnst',0,0,'bcnstv',0,0)
if (nstack.gt.1) call plbox('bcst', 0,0,'bcnstv',0,0)
else
if (nstack.eq.1) call plbox('bcnstl',0,0,'bcnstv',0,0)
if (nstack.gt.1) call plbox('bcstl', 0,0,'bcnstv',0,0)
endif
endif
if (nstack.eq.1) call pllab (ix,iy," ")
if (nstack.gt.1) call pllab (" ", iy, " ")
c-----
if (nstack.eq.kstack) then
if (ncall.eq.kstack ) then
call plcol(kcolax)
call pllab (" ",iy," ")
call plcol(kcoltext)
chsizsav = chsiz
chsiz = 0.7 * chsiz
call plschr (zero4, chsiz)
call plmtex ('t',dispm, posm, justm, khead1)
call plcol(kcoltext)
if (kdatewrite.gt.0)
& call plmtex ('t',dispm2, posm2, justm2, fdate())
chsiz = chsizsav
call plcol(kcolax)
endif
endif
c-----
else
c Here kstack = 1
if (kdebug > 0) print
& "(' Here kstack = 1 - nstack, kstack, ncall =', 3i4)",
& nstack, kstack, ncall
if (kover.eq.1 .and. numrun.ne.92) call pladv(0)
if (kover.gt.1 .and. ncall.eq.1) call pladv(0)
call set_viewport
call plvpor (xvmin, xvmax, yvmin,yvmax)
call plwind (xmin4, xmax4, ymin4, ymax4)
if (kyaxlog.le.0) then
if (kxaxlog.le.0) then
call plbox('bcnst',0,0,'bcnstv',0,0)
else
call plbox('bcnstl',0,0,'bcnstv',0,0)
endif
else
if (kxaxlog.le.0) then
call plbox('bcnst',0,0,'bcnstvl',0,0)
else
call plbox('bcnstl',0,0,'bcnstvl',0,0)
endif
endif
cccccccc call plenv (xmin4, xmax4, ymin4, ymax4, 0, 0)
if (ncall.eq.1 .and. kstack.eq.1) then
call plcol(kcoltext)
call pllab(ix,iy," ")
chsizsav = chsiz
chsiz = 0.8 * chsiz
call plschr (zero4, chsiz)
call plmtex ('t',dispm, posm, justm, khead1)
if (kdatewrite.gt.0)
& call plmtex ('t',dispm2, posm2, justm2, fdate())
endif
if (ncall.ne.1) call pllab(ix,iy," ")
endif
c---------------------------------------------------
c
if (kover.eq.1 .and.
& ncall.eq. nxplots*nyplots*kstack) ncall = 0
c
if (kover.gt.1 .and. ncall.eq.kover) ncall = 0
c
c Draw the line(s) through the data or draw points
c jcolodfdata solid line or points
c 2 Aug 2002
kcol4dat = jcolodfdata(numrun)
kcol4b = jcolodfbayes(numrun)
kcol4th = jcolodftheory(numrun)
kcol4box = kboxcol(numrun,1)
iytest = iy(1:13)
if (iytest.eq.'Bayes / Exp' .or.
& iytest.eq.'<Bayes / Exp>') kcol4dat = kcol4b
c
if (kdebug > 1) print *, ' numrun, kcolorodf(numrun) =',
& numrun, kcolorodf(numrun)
call plcol(kcol4dat)
c
kdash = jdashodf(numrun)
call setlinestyle ! Set line style for data
c
call plwid(klinewid)
if (kdebug > 0) print *,
&' pdover: n x4 y4 yeb4lo yeb4hi y4temp y4temp1'
do n=1,kpts
kpdex = kplo+n-1
x4(n) = x(kpdex)
y4(n) = y(kpdex)
y4err(n) = yerr(kpdex)
y4temp(n) = temp(kpdex)
y4temp1(n) = temp1(kpdex)
yeb4lo(n) = y4(n) - y4err(n)
yeb4hi(n) = y4(n) + y4err(n)
c
if (kdebug > 1) print '(i4, f10.3, 5f12.6)', n, x4(n), y4(n),
& yeb4lo(n), yeb4hi(n),y4temp(n), y4temp1(n)
enddo
c----------------------------------------------------------------------
if(kdebug.gt.0) print *,' pdover: kpts, kaverage =',kpts,kaverage
if (kaverage > 1)
& call paverage (kpts, x4, y4, y4temp, y4temp1,y4err,yeb4lo,yeb4hi)
if(kdebug > 0) print *,' pdover: after paverage: kpts =', kpts
c----------------------------------------------------------------------
if (kdebug > 1) print *,
&' pdover: n x4 y4 yeb4lo yeb4hi'
if (kyaxlog > 0) then
do n=1,kpts
if (kxaxlog.gt.0) then
if (x4(n).le.0.0) x4(n) = 1.0e-8
x4(n) = alog10(x4(n))
endif
yyylo = y4(n) - y4err(n)
yyyhi = y4(n) + y4err(n)
if (y4(n).le.0.0) y4(n) = 1.0e-22
if (y4temp(n).le.0.0) y4temp(n) = 1.0e-22
if (y4temp1(n).le.0.0) y4temp1(n)= 1.0e-22
y4(n) = alog10(y4(n))
y4temp(n) = alog10(y4temp(n))
y4temp1(n)= alog10(y4temp1(n))
yeb4lo(n) = -100.
yeb4hi(n) = -100.
if (yyylo.gt.0.0) yeb4lo(n) = dlog10(yyylo)
if (yyyhi.gt.0.0) yeb4hi(n) = dlog10(yyyhi)
if(kdebug > 1) print '(i4, f10.3, 5f12.6)',
& n, x4(n), y4(n),yeb4lo(n), yeb4hi(n)
enddo
endif
if (kxaxlog > 0. and. kyaxlog.le.0) then
do n=1,kpts
if (x4(n).le.0.0) x4(n) = 1.0e-8
x4(n) = alog10(x4(n))
enddo
endif
c----------------------------------------------------------------------
ncode = kpoints
if (kerrbar > 0) then
if (kpoints.gt.0) call plpoin (kpts,x4,y4, ncode)
ytick = 0.0
call plsmin (ytick, ytick)
call plerry (kpts,x4,yeb4lo, yeb4hi)
call plsmin (ytick, yone)
endif
if (kpoints.gt.0) then
call plpoin (kpts,x4,y4, ncode)
else if (kpoints.eq.0) then
c
x4his(2) = 0.5*(x4(1) + x4(2)) ! histogram
x4his(1) = 2.0*x4(1) - x4his(2)
y4his(1) = y4(1)
y4his(2) = y4(1)
do n=2,kpts
ip1 = 2*n -1
x4his(ip1) = x4his(ip1-1)
if (n.eq.kpts) then
x4his(ip1+1) = 1.5*x4(n) - 0.5*x4(n-1)
else
x4his(ip1+1) = 0.5*(x4(n) + x4(n+1))
endif
y4his(ip1) = y4(n)
y4his(ip1+1) = y4(n)
enddo
kpts2 = 2*kpts
call plline (kpts2,x4his,y4his)
endif
34 continue
iytest = iy(1:13)
if (kdebug > 1) print *, ' iytest = ', iytest
if (iytest.eq.'Bayes / Exp' .or.
& iytest.eq.'<Bayes / Exp>') then
if (kdebug > 0) print * ,' drawing cyan line for Bayes / Exp'
call drawcyanline (xmin, xmax)
endif
if (kcurves.ge.2) then
c
call plcol(kcol4b)
c 14 Aug 2002
c
if (kdedet > 0 .and. kcurves.eq.2) call plcol(kcol4th)
c
if (numrun.eq.12) call plcol(1) ! Red solid line
c
kdash = jdashbayes(numrun)
if (kdedet > 0) kdash = jdashtheory(numrun)
call setlinestyle ! Set line style for theory or Bayes
cccccccc call pllsty(1)
call plwid(klinewid)
if (ktheoryhist > 0) then
call plot_hist (kpts, x4, y4temp)
else
call plline (kpts,x4,y4temp)
endif
endif
if (kcurves.ge.3) then
c yellow dashed line
call plcol(kcol4th)
call pllsty(klinewid)
c 13 Oct 2004
kdash = jdashtheory(numrun)
call setlinestyle ! Set line style for
theory
c
call plwid(klinewid)
call plline (kpts,x4,y4temp1)
endif
if (kdebug > 0) print *, 'iy =',iy
if (numrun.eq.12 .and. iy.ne.'Theory2 / Theory1'
& .and. iy.ne.'Bayes2 / Bayes1'
& .and. iy.ne.'Theory2 - Theory1'
& .and. iy.ne.'Bayes2 - Bayes1'
& .and. iy.ne.'Data2 - Data1'
& .and. iy.ne.'Data2 / Data1'
& .and. iy.ne.'<Data2> / <Data1>'
& .and. iy.ne.'<Data2> - <Data1>')
& then
c Overlay Run 1 & Run 2
c Cyan dotted line
call plcol(11)
call pllsty(2)
call plline (kpts,x4,y4temp1)
c
c * * * * * * * * * * * * * * * * * * * * Legend
c
x4(1) = xmin + 0.02 *(xmax-xmin)
x4(2) = xmin + 0.06 *(xmax-xmin)
y4(1) = ymax - 0.04 *(ymax-ymin)
y4(2) = y4(1)
npts2 = 2
if (kdebug > 0) print *,'legend: ', x4(1),x4(2),y4(1),y4(2)
call plline (npts2,x4,y4)
c Red solid line
y4(1) = y4(1) + 0.02 *(ymax-ymin)
y4(2) = y4(1)
call plcol(1)
call pllsty(1)
call plline (npts2,x4,y4)
c
endif
c
c * * * * * * * * * * * * * * * * * * * * Annotate with Figure caption
c
if (kcaption > 0) then
if (nxplots*nyplots ==1 .and. nstack < 2 .or.
& nxplots*nyplots > 1 .and. ncall == 0) then
chsiz = capsize
call plschr (zero4, chsiz)
call plcol(13) ! magenta
dispcaparg = dispcap
if (nyplots > 1) dispcaparg = dispcaparg -1.0
if (nyplots > 2) dispcaparg = dispcaparg -1.0
if (nxplots < 2)
& call plmtex ('b',dispcaparg, poscap, justcap, caption_string)
if (kdebug > 0)
& print "('dispcaparg, poscap, justcap',/3f8.3,i10)",
& dispcaparg, poscap, justcap
endif
endif
c....................................................................
c....................................................................
if (kstack.eq.1 .or. nstack.eq.kstack) then
c
c * * * * * * * * * * * * * * * * * * * * * * * Annotate with CHISQ/N
c
if (ichisqn.ne.0 .and. kchiwrite.gt.0) then
chsiz = 0.4 ! character size
if (nxplots*nyplots.eq.1) chsiz = 0.5
call plschr (zero4, chsiz)
call plcol(3) ! Green
call plmtex ('t',dispchi, poschi, justchi, kchisqn)
endif
c
c
c * * * * * * * * * * * * * * * * * * * * * * * Annotate with kaverage
c
if (kaverage > 1 .and. kavgwrite.ne.0) then
c character size
chsiz = 0.4
if (nxplots*nyplots.eq.1) chsiz = 0.5
call plschr (zero4, chsiz)
call plcol(13) ! magenta
write (unit=kavgtext, fmt="('<', i2, '>')") kaverage
call plmtex ('t',dispavg, posavg, justavg, kavgtext)
endif
c
if (knormwrite > 0 .and. kodfwrite > 0) then
c * * * * * * * * * * * * * * * * * * * * * * * Annotate with ANORMODF
c
chsiz = 0.4 ! character size
if (nxplots*nyplots.eq.1) chsiz = 0.5
call plschr (zero4, chsiz)
call plcol(13) ! magenta
if (ncall.eq.1) then
dispnorm = -44.
write (unit=kekevtext,fmt="(' ANORM')")
call plmtex ('t',dispnorm, posnorm, justnorm, kekevtext)
endif
write (unit=kekevtext,fmt='(f7.3)') anormodf(numrun)
dispnorm = dispnorm - 2.0
call plmtex ('t',dispnorm, posnorm, justnorm, kekevtext)
c
endif
c
if (kodfwrite.gt.0) then
c * * * * * * * * * * * * * * * * * * * * * * * Write odf name(s)
if (kdebug > 0) then
print "(' ------ odf names after CHISQ:',
&/' ncall nstack kstack numrun klength klength2',6i8)",
& ncall, nstack, kstack, numrun, klength, klength2
print *, ' chsiz = ', chsiz
endif
if (klength > 40 .or. klength2 > 40) then
chsizsav = chsiz
chsiz = 0.4
endif
call plschr (zero4, chsiz)
call plcol(1) ! red
if (numrun > 0) then
chsizsav = chsiz
chsiz = 0.4
call plschr (zero4, chsiz)
if (numrun.eq.12) then
call plmtex ('t',dispodf,posodf, justodf, ' ')
posodf = posodf + 0.07
call plmtex ('t',dispodf,posodf, justodf, odfnames(1))
posodf = posodf - 0.07
else
call plmtex ('t',dispodf, posodf, justodf, odfname)
endif
endif
if (numrun.eq.12) then
call plcol(11) ! cyan
chsizsav = chsiz
chsiz = 0.4
call plschr (zero4, chsiz)
call plmtex ('t',dispodf2,posodf2,justodf2, '- - -')
posodf2 = posodf2 + 0.07
call plmtex ('t',dispodf2,posodf2,justodf2,odfnames(2))
posodf2 = posodf2 - 0.07
endif
chsiz = chsizsav
call plschr (zero4, chsiz)
endif
c * * * * * * * * * * * * * * * * * * * * * * * end of odf name write
endif
c......................................................................
if (kdebug > 2) print *, ' keres nglo nghi ngbig'
if (kdebug > 2) print *, keres, nglo,nghi,ngbig
c write resonance energies on plots
if (keres.ne.0) then
ysamlo = - 88.
if (kstack.gt.1) ysamlo = ysamlo/float(kstack)
if (ngbig > 0) then
npts2 = 2
chsiz = 0.4
cccccccc call pllsty(2)
call plstyl (1, 600, 600)
call plschr (zero4, chsiz)
k = 0
cccccccc nglo = 1
nghi = ngbig
if (keres.eq.2) then
nglo = nslbw
nghi = nglo
print *, 'nslbw =', nslbw
endif
if (kdebug.gt.2) print * ,' ng jr samer(jr,ng)'
do ng=nglo,nghi
do jr=1, ires(ng)
k = k + 1
ersave(k) = samer(jr,ng)
if (kdebug.gt.2) print 44, ng, jr, samer(jr,ng)
44 format (2i5, 5f11.4)
enddo
enddo
kup = k
ccccc ysam = 0.9 * ymax4
ysam = -2.0
indx(1) = 1
if (kup.gt.1) call indexx (kup, ersave, indx)
call plcol(13) ! magenta
if (kdebug > 2) then
print * ,' k k samerout'
print *, 'rplnum: samerout disp pos just ksamerout'
endif
do k=1,kup
samerout = ersave(indx(k))
if (samerout.lt.xmin4) cycle
if (samerout.gt.xmax4) ysam = ysam - 2.0
if (ysam.ge.ysamlo)
& call rplnum (chsiz, xmin4,xmax4, ymin4,ymax4, samerout, ysam)
c
c Vertical markers
y4(1) = ymin4
y4(2) = ymin4 + 0.06 * (ymax4 - ymin4)
x4(1) = samerout
x4(2) = x4(1)
call plline (npts2, x4, y4)
enddo
c
if (keres.ne.2 .and. nyplots < 3) then
c Group numbers above vertical markers
if (nyplots == 1) then
deltax = 0.0
ysam = -81.7
if (kstack > 1)
& ysam = -85.7
else
deltax = 0.007 * (xmax4-xmin4)
ysam = -57.5
endif
if (kcaption > 0) ysam = ysam + 5.2/float(kstack)
do ng=nglo,nghi
do jr=1, ires(ng)
kgroup = ng
samerout = samer(jr,ng) + deltax
if (samerout.lt.xmin4) cycle
if (samerout.gt.xmax4) cycle
call iplnum (chsiz, xmin4,xmax4, ymin4,ymax4,
& samerout, kgroup, ysam)
enddo
enddo
endif
c
if (kgamres.gt.0) then
c......................................................................
c Resonance widths below x axis
if (kdebug.gt.2) then
print *,'rwidplnum: samerout ysam reswidth ksamerout disp pos'
print*, ' ng jr gamr2(jr,ng) samgamr2(jr,ng) reswidth'
endif
c---
if (nyplots.lt.2) call widthlab (kgamres, numpeaks)
c---
do ng=nglo,nghi
do jr=1, ires(ng)
kgroup = ng
samerout = samer(jr,ng)
if (samerout.lt.xmin4) cycle
if (samerout.gt.xmax4) cycle
c neutron widths (magenta)
call plcol(13)
ysam = -95.5
if (kcaption > 0) ysam = ysam + 5.2
reswidth = samgamr(jr,ng)
call rwidplnum (chsiz, xmin4,xmax4, ymin4,ymax4,
& samerout, reswidth, ysam)
ysam = -97.2
if (kcaption > 0) ysam = ysam + 5.2
reswidth = samgamr2(jr,ng)
if (reswidth .ne. 0.0)
& call rwidplnum (chsiz, xmin4,xmax4, ymin4,ymax4,
& samerout, reswidth, ysam)
if (numpeaks.gt.0) then
c total widths at top of plot for PKFIND
ysam = 7.8
reswidth = gamr(jr,ng)
if (reswidth .ne. 0.0)
& call rwidplnum (chsiz, xmin4,xmax4, ymin4,ymax4,
& samerout, reswidth, ysam)
cycle
endif
if (kgamres.lt.2) cycle
c......................................................................
c % diff in widths at top of plot
ysam = 7.8
reswidth = samgamr(jr,ng)
c (CYAN)
call plcol(11)
if (reswidth .ne. zero4) then
if (gamr(jr,ng) .ne. zero4)
& reswidth = 100.* (samgamr(jr,ng)/gamr(jr,ng) - yone)
if (abs(reswidth) .gt. 0.1)
& call rwidplnum (chsiz, xmin4,xmax4, ymin4,ymax4,
& samerout, reswidth, ysam)
endif
ysam = 6.1
reswidth = samgamr2(jr,ng)
call plcol(11)
if (reswidth .ne. zero4) then
if (gamr2(jr,ng) .ne. zero4)
& reswidth = 100.* (samgamr2(jr,ng)/gamr2(jr,ng) - yone)
if (abs(reswidth) > 0.1)
& call rwidplnum (chsiz, xmin4,xmax4, ymin4,ymax4,
& samerout, reswidth, ysam)
if (kdebug > 0) print '(2i5, 3f12.4)',
& ng, jr, gamr2(jr,ng), samgamr2(jr,ng), reswidth
endif
c ()
call plcol(13) ! magenta
c END % diff in widths
enddo
enddo
c......................................................................
endif
endif
endif
c......................................................................
c annotate with text string
c
if (nxplots*nyplots.lt.3) then
chsiz = 0.6
else
chsiz = 0.5
endif
call plschr (zero4, chsiz)
c
call textwrite (ncallann, kcoltext, chsiz)
c------------------------------------------------------------------------
c......................................................................
c DRAW BOX
c
kstr = numbox(ncallann)
if (kstr < 1) return
c
do k=1,kstr
karrow = kboxarrow(ncallann,k)
kcol4box = kboxcol(ncallann,k)
call plcol(kcol4box)
kdash = kboxdash(ncallann,k)
call setlinestyle
npts2 = 2
if (kdebug > 0) print *,'BOX: ncallann:', ncallann
if (karrow > 0) then
kkcall = k
call drawarrow (ncallann, karrow, kkcall)
else
x4(1) = xbox1(ncallann,k)
x4(2) = xbox2(ncallann,k)
y4(1) = ybox1(ncallann,k)
y4(2) = y4(1)
if (kdebug > 0) print *,'BOX: ', x4(1),x4(2),y4(1),y4(2)
call plline (npts2,x4,y4)
x4(1) = xbox2(ncallann,k)
x4(2) = x4(1)
y4(1) = ybox1(ncallann,k)
y4(2) = ybox2(ncallann,k)
call plline (npts2,x4,y4)
x4(1) = xbox2(ncallann,k)
x4(2) = xbox1(ncallann,k)
y4(1) = ybox2(ncallann,k)
y4(2) = y4(1)
call plline (npts2,x4,y4)
x4(1) = xbox1(ncallann,k)
x4(2) = x4(1)
y4(1) = ybox2(ncallann,k)
y4(2) = ybox1(ncallann,k)
call plline (npts2,x4,y4)
endif
enddo
c------------------------------------------------------------------------
return
end
c***************************************************************************
subroutine drawarrow (ncallann, karrow, k) ! draw arrow
c
use rsap_constants
use rsap_BLK2
use rsap_BLK6
use rsap_BLK12
integer*4 ncallann, karrow, ksymarrow, k
integer*4 :: npts2 = 2
real*4 x4(4), y4(4)
real*4 x4j1, y4j1, x4j2, y4j2
real*8 dxarrow, dyarrow, bigrarrow, tarrow, theta, beta
real*8 :: facpsi = 1.d0, aspectrat = 1.38d0, facplus, facminus
real*8 :: rarrow = 0.014d0, sizearrow = 0.01d0
real*8 :: cospsi = 0.8660254d0, sinpsi = 0.5d0
c
x4(1) = xbox1(ncallann,k)
x4(2) = xbox2(ncallann,k)
y4(1) = ybox1(ncallann,k)
y4(2) = ybox2(ncallann,k)
c
x4j1 = xbox1(ncallann,k)
x4j2 = xbox2(ncallann,k)
y4j1 = ybox1(ncallann,k)
y4j2 = ybox2(ncallann,k)
if (kdebug > 0) print *,'drawarrow: ', x4(1),x4(2),y4(1),y4(2)
cccccccc call plline (npts2,x4,y4)
call pljoin (x4j1,y4j1, x4j2,y4j2)
if (karrow < 1) return
dxarrow = sizearrow * (xmax - xmin)
dyarrow = sizearrow * (ymax - ymin)
c
if (y4j1 .eq. y4j2) then ! horizontal line
theta = zero
if(x4j2 < x4j1) then
x4j1 = x4j2 + dxarrow
else
x4j1 = x4j2 - dxarrow
endif
y4j1 = y4j2 + 0.5*dyarrow
call pljoin (x4j1,y4j1, x4j2,y4j2)
y4j1 = y4j2 - 0.5*dyarrow
call pljoin (x4j1,y4j1, x4j2,y4j2)
else if (x4j1 .eq. x4j2) then
theta = 90.d0 ! vertical line
if(y4j2 < y4j1) then
y4j1 = y4j2 + dyarrow*aspectrat
else
y4j1 = y4j2 - dyarrow*aspectrat
endif
x4j1 = x4j2 + 0.5*dxarrow/aspectrat
call pljoin (x4j1,y4j1, x4j2,y4j2)
x4j1 = x4j2 - 0.5*dxarrow/aspectrat
call pljoin (x4j1,y4j1, x4j2,y4j2)
else
theta = datan((y4j2-y4j1) * dxarrow /(dyarrow * (x4j2-x4j1)))
theta = datan2((y4j2-y4j1)/dyarrow , (x4j2-x4j1)/dxarrow)
bigrarrow = ((y4j2-y4j1)/(ymax-ymin))**2
& + ((x4j2-x4j1)/(xmax-xmin))**2
bigrarrow = dsqrt(bigrarrow)
tarrow = bigrarrow**2 + rarrow**2 - two*rarrow*bigrarrow*cospsi
tarrow = dsqrt(tarrow)
beta = dasin(rarrow * sinpsi / tarrow)
beta = theta - beta
if (kdebug > 0) then
print '("drawarrow: x1,y1,x2,y2 =",4f10.5)', (x4(j),y4(j),j=1,2)
print '(" theta bigrarrow tarrow beta")'
print '(4f9.5)', theta, bigrarrow, tarrow, beta
endif
x4(3) = tarrow * cos(beta)
y4(3) = tarrow * sin(beta)
x4(4) = tarrow * cos(two*theta - beta)
y4(4) = tarrow * sin(two*theta - beta)
cccccccc print '(" x3 y3 x4 y4")'
cccccccc print '(4f9.5)', x4(3), y4(3), x4(4), y4(4)
do j=3,4
x4(j) = x4(1) + x4(j) * (xmax-xmin)
y4(j) = y4(1) + y4(j) * (ymax-ymin)
enddo
cccccccc print '(" x3 y3 x4 y4")'
cccccccc print '(4f9.3)', x4(3), y4(3), x4(4), y4(4)
c
y4j1 = y4(3) ; x4j1 = x4(3)
call pljoin (x4j1,y4j1, x4j2,y4j2)
y4j1 = y4(4) ; x4j1 = x4(4)
call pljoin (x4j1,y4j1, x4j2,y4j2)
endif
return
end
c***************************************************************************
subroutine setlinestyle ! Set PLPLOT line style
c
use rsap_BLK2
use rsap_BLK12
if (kdash.eq.0) then
call pllsty(1)
else if (kdash < 0) then
jdotlength = 100*iabs(kdash)
call plstyl (1, jdotlength,jdotlength)
else
call pllsty(kdash)
endif
return
end
c***************************************************************************
subroutine caption_write (device) ! Annotate with Figure caption
c
use rsap_pdo_constants
use rsap_viewport
use rsap_BLK3
use rsap_BLK12
character*8 device
c
print *, 'caption_write: ', device, kcaption, nxplots
if (kcaption < 1 .or. nxplots < 2) return
call pladv(0)
print *, 'caption_write: pladv called'
call set_viewport
call plvpor (xvmin, xvmax, yvmin,yvmax)
print *, 'caption_write: plvpor called'
c
chsiz = capsize
call plschr (zero4, chsiz)
call plcol(13) ! magenta
dispcaparg = dispcap
if (nyplots > 1) dispcaparg = dispcaparg -1.0
if (nyplots > 2) dispcaparg = dispcaparg -1.0
call plmtex ('b',dispcaparg, poscap, justcap, caption_string)
if (kdebug.gt.0)
& print "('dispcaparg, poscap, justcap',/3f8.3,i10)",
& dispcaparg, poscap, justcap
return
end
c***************************************************************************
subroutine int_sig_over_e (numodf, nrecdo, yplot, timep)
c
use rsap_constants
use rsap_BLK2
use rsap_BLK3
use rsap_BLK4
use rsap_BLK5
use rsap_BLK6
use rsap_BLK7
use rsap_BLK8
use rsap_BLK9
use rsap_BLK12
use rsap_yaxislabels
c
implicit real*8 (a-h,o-z)
implicit integer*8 (i-n)
integer*8 numodf, nrecdo
real*8 yplot(1), timep(1)
c
elow = xmin ; ehigh = xmax
if (kdebug > 0) then
print *, ' elow ehigh big(1,1) big(nrecdo,1)'
print *, elow, ehigh, big(1,1,numodf), big(nrecdo,1,numodf)
endif
if (elow.lt.big(1,1,numodf)) then
elow = big(1,1,numodf)
print *, ' elow set to min ODF energy: ', elow
endif
if (ehigh.gt.big(nrecdo,1,numodf) .or. ehigh.lt.elow) then
ehigh = big(nrecdo,1,numodf)
print *, ' ehigh set to max ODF energy: ', ehigh
endif
print 1500, elow, ehigh
1500 format(' integrate dE * sig / E from ',f12.6,' to',f12.6,' keV')
c
energyscale = one
if (kplev.eq.1) energyscale = 1000.d0
c
if (kdebug.gt.0) print *,
&' j energy(keV) sigtheory(b) dE(keV) dE*sig/E integral'
sum = zero
sumtrap = zero
c trapezoidal rule
do j=1,nrecdo
evalue = big(j,1,numodf)
if (evalue .lt. elow) cycle
if (evalue .lt. zero) cycle
if (evalue .gt. ehigh) exit
sigtheory = big(j,4,numodf)
if (j.eq.1) then
deltae = big(2,1,numodf) - big(1,1,numodf)
else
if (j.eq.nrecdo) then
deltae = big(nrecdo,1,numodf) - big(nrecdo-1,1,numodf)
else
deltae = 0.5d0 * (big(j+1,1,numodf) - big(j-1,1,numodf))
endif
endif
if (deltae.lt.zero) then
print *, ' deltae < 0, energy ignored ',evalue
cycle
endif
timep(j) = energyscale* evalue
sigtheory = big(j,4,numodf)
yplot(j) = deltae * sigtheory / evalue
sum = sum + yplot(j)
temp (j) = sum
if (j.lt.nrecdo) then
deltaetrap = big(j+1,1,numodf) - big(j,1,numodf)
fac = big(j+1,4,numodf)/big(j+1,1,numodf)
& + big( j,4,numodf)/big( j,1,numodf)
yplot(j) = deltaetrap * fac * 0.5d0
sumtrap = sumtrap + deltaetrap * fac * 0.5d0
temp (j) = sumtrap
endif
if (kdebug > 0) print '(i6, 1p6e14.4)',
& j, evalue, sigtheory, deltae, yplot(j), sum, sumtrap
enddo
print *, ' integral = ', sum, ' barns'
return
end
--
Summary: gfortran installed on Mac OS X using fink commander.
Product: gcc
Version: unknown
Status: UNCONFIRMED
Severity: critical
Priority: P2
Component: fortran
AssignedTo: unassigned at gcc dot gnu dot org
ReportedBy: sayerro at ornl dot gov
CC: gcc-bugs at gcc dot gnu dot org
http://gcc.gnu.org/bugzilla/show_bug.cgi?id=20430
