259 lines
5.9 KiB
Fortran
259 lines
5.9 KiB
Fortran
subroutine dat_rita
|
|
c -------------------
|
|
|
|
external dat_rita_desc
|
|
external dat_rita_opts
|
|
external dat_rita_read
|
|
|
|
integer dtype/0/
|
|
|
|
call dat_init_desc(dtype, dat_rita_desc)
|
|
call dat_init_opts(dtype, dat_rita_opts)
|
|
call dat_init_read(dtype, dat_rita_read)
|
|
end
|
|
|
|
|
|
subroutine dat_rita_desc(text)
|
|
! ------------------------------
|
|
character*(*) text ! (out) description
|
|
|
|
! type description
|
|
! ----------------------------------
|
|
text='RITA TASCOM (RITA) format'
|
|
end
|
|
|
|
|
|
subroutine dat_rita_opts
|
|
! ------------------------
|
|
print '(x,a)'
|
|
1,'x,y,mon: xaxis,yaxis,monitor to be choosen'
|
|
end
|
|
|
|
|
|
subroutine dat_rita_read
|
|
1 (lun, forced, nread, putval, nmax, xx, yy, ss, ww)
|
|
! ----------------------------------------------------
|
|
implicit none
|
|
|
|
integer lun ! (in) logical unit number (file will be closed if successful)
|
|
integer forced ! 0: read only if type is sure; 1: forced read
|
|
integer nread ! (out) >=0: = number of points read, file closed
|
|
! -1: not correct type, file rewinded
|
|
! -2: correct type, but unreadable, file rewinded
|
|
external putval ! (in) subroutine to put name/value pairs.
|
|
! for numeric data: call putval('name', value) ! value must be real
|
|
! for character data: call putval('name=text', 0.0)
|
|
integer nmax ! max. number of points
|
|
real xx(*) ! x-values
|
|
real yy(*) ! y-values
|
|
real ss(*) ! sigma
|
|
real ww(*) ! weights (original monitor)
|
|
|
|
! local
|
|
integer mcol
|
|
parameter (mcol=64)
|
|
|
|
character line*1024, xaxis*8, yaxis*8, monam*8, name*8, col1*8
|
|
character fih*132
|
|
real values(64), r, f, s, y, ymon
|
|
real tt, tm, ts, td, xmon
|
|
integer i,j,l,ncol,ccol,pcol,xcol,tcol
|
|
integer iostat
|
|
|
|
nread=0
|
|
10 read(lun,'(a)',err=98,end=98) line
|
|
if (line .eq. ' ') goto 10
|
|
|
|
if (forced .le. 0) then
|
|
if (line(1:4).ne.'#fdt') goto 100
|
|
endif
|
|
|
|
xaxis=' '
|
|
yaxis=' '
|
|
monam='MON'
|
|
|
|
call dat_start_options
|
|
call dat_str_option('x', xaxis)
|
|
call dat_str_option('y', yaxis)
|
|
call dat_str_option('mon', monam)
|
|
|
|
call str_upcase(xaxis, xaxis)
|
|
call str_upcase(yaxis, yaxis)
|
|
call str_upcase(monam,monam)
|
|
|
|
call dat_group(1, putval)
|
|
call putval('Instrument=RITA',0.0)
|
|
|
|
11 continue
|
|
call str_trim(line, line, l)
|
|
l=max(10,min(l,len(line)))
|
|
|
|
if (line(1:4) .eq. '#fdt') then
|
|
call dat_group(2, putval)
|
|
i=index(line(10:l), ' ')+9
|
|
call putval('OrigFile='//line(10:i), 0.0)
|
|
call dat_group(1, putval)
|
|
call putval('Date='//line(i+1:l), 0.0)
|
|
else if (line(1:4) .eq. '#txt') then
|
|
call dat_group(1, putval)
|
|
call putval('Title='//line(10:l), 0.0)
|
|
else if (line(1:4) .eq. '#cmd') then
|
|
call dat_group(2, putval)
|
|
call putval('cmd='//line(10:l), 0.0)
|
|
else if (line(1:4) .eq. '#fih') then
|
|
fih=line(5:)
|
|
else if (line(1:4) .eq. '#fhp') then
|
|
call dat_group(2, putval)
|
|
do i=1,mcol
|
|
values(i)=0.0
|
|
enddo
|
|
ncol=0
|
|
i=1
|
|
read(line(5:), *, err=15,end=15) values
|
|
15 call str_get_elem(fih, i, name)
|
|
if (name .ne. ' ') then
|
|
ncol=ncol+1
|
|
call putval(name, values(ncol))
|
|
goto 15
|
|
endif
|
|
else if (line(1:4) .eq. '#fip') then
|
|
line(1:4)=' '
|
|
goto 12
|
|
else if (line(1:4) .eq. '#plp') then
|
|
continue ! ignore
|
|
else if (line(1:4) .eq. '#pls') then
|
|
continue ! ignore
|
|
else if (line(1:4) .eq. '#fmp') then
|
|
continue ! ignore
|
|
else if (line(1:4) .eq. '#plv') then
|
|
i=5
|
|
call str_get_elem(line, i, name)
|
|
if (xaxis .eq. ' ') xaxis=name
|
|
call str_get_elem(line, i, name)
|
|
if (yaxis .eq. ' ') yaxis=name
|
|
else if (line(1:4) .eq. '#com') then
|
|
continue
|
|
elseif (line(1:l) .ne. ' ') then
|
|
print *,'unknown header line: ',line(1:l)
|
|
endif
|
|
read(lun, '(a)') line
|
|
goto 11
|
|
|
|
12 i=1
|
|
line(len(line):len(line))=' '
|
|
ncol=0
|
|
ccol=0
|
|
pcol=0
|
|
tcol=0
|
|
xcol=0
|
|
col1=' '
|
|
ymon=0
|
|
31 do while (line(i:i) .eq. ' ')
|
|
i=i+1
|
|
if (i .gt. len(line)) goto 39
|
|
enddo
|
|
l=i
|
|
do while (line(i:i) .ne. ' ')
|
|
i=i+1
|
|
enddo
|
|
ncol=ncol+1
|
|
if (ncol .eq. 1) col1=line(l:i)
|
|
if (line(l:i) .eq. yaxis .or. yaxis .eq. ' '
|
|
1 .and. line(l:i) .eq. 'I') ccol=ncol
|
|
if (line(l:i) .eq. monam) pcol=ncol
|
|
if (line(l:i) .eq. xaxis) xcol=ncol
|
|
if (line(l:i) .eq. 'TT') tcol=ncol
|
|
if (tcol .eq. 0 .and. line(min(i,l+1):i) .eq. 'TEM') tcol=ncol
|
|
goto 31
|
|
|
|
39 if (yaxis .eq. ' ') yaxis='I'
|
|
if (ccol .eq. 0) then
|
|
print *,'no values found for ',yaxis
|
|
goto 99
|
|
endif
|
|
if (pcol .eq. 0) then
|
|
xmon=0
|
|
read(monam, *, iostat=iostat) xmon
|
|
if (xmon .le. 0) then
|
|
if (monam .eq. ' ') monam='Monitor'
|
|
print *,'no values found for ',monam,' take 1'
|
|
xmon=1
|
|
endif
|
|
endif
|
|
if (xcol .eq. 0) then
|
|
if (xaxis .ne. ' ') then
|
|
print *,'no values found for ',xaxis,', take ',col1
|
|
endif
|
|
xcol=1
|
|
xaxis=col1
|
|
endif
|
|
|
|
call dat_group(2, putval)
|
|
call putval('XAxis='//xaxis, 0.0)
|
|
if (yaxis .eq. 'I') then
|
|
call putval('YAxis=Intensity', 0.0)
|
|
else
|
|
call putval('YAxis='//yaxis, 0.0)
|
|
endif
|
|
call dat_group(1, putval)
|
|
|
|
l=min(mcol,max(xcol,pcol,ccol,tcol))
|
|
f=1.0
|
|
tm=0
|
|
ts=0
|
|
|
|
4 read(lun,*,err=19,end=9) (values(j),j=1,l)
|
|
if (nread .ge. nmax) goto 29
|
|
y=values(ccol)
|
|
if (pcol .eq. 0) then
|
|
r=xmon
|
|
else
|
|
r=values(pcol)
|
|
endif
|
|
if (ymon .eq. 0) then
|
|
ymon=r
|
|
if (r .eq. 0) r=1.
|
|
endif
|
|
if (r .le. 0.0) r=ymon
|
|
f=ymon/r
|
|
if (f .le. 0.0) f=1.0
|
|
if (y .gt. 0) then
|
|
s=sqrt(y) ! statistical error of detector
|
|
else
|
|
s=1
|
|
endif
|
|
nread=nread+1
|
|
xx(nread)=values(xcol)
|
|
yy(nread)=y*f
|
|
ss(nread)=s*f
|
|
ww(nread)=r
|
|
if (tcol .ne. 0) then
|
|
tt=values(tcol)
|
|
td=(tt-tm)/nread
|
|
tm=tm+td ! mean temp.
|
|
ts=ts+(tt-tm)**2+td*td*(nread-1) ! sum of (temp(i)-mean)**2
|
|
c print *,'temp',tt,tm,ts
|
|
endif
|
|
goto 4
|
|
|
|
9 close(lun)
|
|
call putval('Monitor', ymon)
|
|
if (tcol .ne. 0) then
|
|
call putval('Temp', tm)
|
|
if (nread .gt. 1) call putval('dTemp', sqrt(ts/(nread-1)))
|
|
endif
|
|
return
|
|
|
|
19 print *,'DAT_RITA: error at point ',nread
|
|
goto 4
|
|
29 print *,'DAT_RITA: too many points'
|
|
goto 100
|
|
98 if (forced .le. 0) goto 100
|
|
99 print *,'DAT_RITA: error during read'
|
|
rewind lun
|
|
nread=-2
|
|
return
|
|
100 nread=-1
|
|
rewind lun
|
|
end
|