Files
fit/gen/dat_utils.f
2022-08-19 15:22:33 +02:00

627 lines
12 KiB
Fortran

subroutine dat_filelist(file_list)
include 'dat_utils.inc'
character*(*) file_list
lpos=0
filelist=file_list
end
subroutine dat_init_list
include 'dat_utils.inc'
call dat_init
lout=0
end
subroutine dat_read_next(putval, nmax, nread, xx, yy, ss, ww)
external putval
integer nread
integer nmax ! max. number of points
real xx(nmax) ! x-values
real yy(nmax) ! y-values
real ss(nmax) ! sigma
real ww(nmax) ! weights
include 'dat_utils.inc'
logical silent
call dat_init
call dat_get_silent(silent)
10 if (lpos .ge. len(filelist)) then
nread=-1
else
if (nmax .eq. 0) then
print *,'no more datapoints accepted'
nread=0
else
if (silent) call dat_silent
call dat_open_next(filelist, lpos, outlist, lout
1 , putval, nmax, nread, xx, yy, ss, ww)
endif
if (nread .le. 0) then
call sym_purge(1) ! a little hack
goto 10
endif
endif
end
subroutine dat_get_filelist(file_list)
character*(*) file_list
include 'dat_utils.inc'
if (lout .lt. 2) then
file_list=' '
else
file_list=outlist(1:lout-1)
endif
end
subroutine dat_read_all(filename, putval, nread, nmax,xx,yy,ss,ww)
character filename*(*)
external putval
integer nread
integer nmax ! max. number of points
real xx(nmax) ! x-values
real yy(nmax) ! y-values
real ss(nmax) ! sigma
real ww(nmax) ! weights
include 'dat_utils.inc'
integer n,k
call dat_init
lpos=0
lout=0
k=0
1 if (k .ge. nmax) then
print *,'no more datapoints accepted'
n=0
else
call dat_open_next(filename, lpos, outlist, lout
1 , putval, nmax-k, n, xx(k+1),yy(k+1),ss(k+1),ww(k+1))
endif
if (n .le. 0) goto 1
k=k+n
if (lpos .lt. len(filename)) goto 1
if (lout .gt. 1) then
print *,'files read: ',outlist(1:lout-1)
nread=k
else
print *,'no files read'
nread=-1
endif
end
subroutine dat_put_filelist(filename)
character filename*(*)
include 'dat_utils.inc'
integer n
external dat_put_none
call dat_init
lpos=0
lout=0
1 call dat_open_next(filename, lpos, outlist, lout
1 , dat_put_none, 0, n, 0.0, 0.0, 0.0, 0.0) ! check syntax only
if (lpos .lt. len(filename)) goto 1
end
subroutine dat_next_filename(filename, len_name)
character filename*(*)
integer len_name
integer n
include 'dat_utils.inc'
external dat_put_none
if (lpos .ge. len(filelist)) then
len_name=0
return
endif
call dat_open_next(filelist, lpos, outlist, lout
1 , dat_put_none, 0, n, 0.0, 0.0, 0.0, 0.0) ! syntax only
call dat_get_filename(filename, len_name)
end
subroutine dat_group(show, putval)
integer show
external putval
call putval('ShowLevel', float(show))
end
subroutine dat_put_none(name, value)
c
c dummy routine
c
character name*(*)
real value
end
subroutine dat_put_one(name, value)
c put a real value: NAME must not contain '='
c
c put a string value: NAME must contain '='
c syntax name=text, VALUE must be 0.0
character name*(*), text*(*)
real value
integer i
character name0*32/' '/, text0*80
save name0, text0
i=index(name, '=')
if (i .gt. 1 .and. i .lt. len(name)) then
name0=name(1:i-1)
text0=name(i+1:)
endif
return
entry dat_get_one(name, text)
name=name0
text=text0
name0=' '
end
subroutine dat_insert_year(str, year, ok)
character str*(*)
integer year
logical ok
integer i,y,m,d
ok=.false.
1 i=index(str,'%%')
if (i .ne. 0) then
if (year .eq. 0) then
call sys_date(y,m,d)
else
y=year
endif
if (str(i:min(i+3,len(str))) .eq. '%%%%') then
write(str(i:i+3), '(i4.4)') mod(y,10000)
else
write(str(i:i+1), '(i2.2)') mod(y,100)
endif
ok=.true.
goto 1
endif
end
subroutine dat_set_options(options)
character options*(*)
character defopt*256/' '/
integer leng/0/, deflen/1/
leng=0
call dat_insert_options(options, leng)
if (defopt(1:deflen) .ne. ' ') then
call dat_insert_options(defopt(1:deflen), leng)
endif
return
entry dat_def_options(options)
! set default options
call str_trim(defopt, options, deflen)
return
entry dat_get_def_options(options)
options=defopt
end
subroutine dat_insert_options(options, leng)
character name*(*), value*(*), options*(*)
integer leng
character opt*256/' '/
integer idx/0/
integer nam1(32), nam2(32)
integer val1(32), val2(32)
integer cnt(32)
character nam*16
integer i,j,nopt,mopt,found,named
integer start,ende
save mopt,nopt,nam1,nam2,val1,val2,cnt
if (leng .eq. 0) then
opt=' '
leng=1
named=0
nopt=0
mopt=0
else
mopt=nopt
endif
start=0
ende=0
do while (ende .ge. 0)
if (nopt .ge. 32) then
print *,'too many options, truncated'
return
endif
call str_split(options, ',', start, ende)
if (start .le. ende) then
call str_first_nonblank(options(start:ende), j)
if (j .gt. 0) start=start+j-1
endif
if (start .le. ende) then
i=index(options(start:ende), '=')
if (i .eq. 0) then
if (options(start:ende) .ne. ' ' .or. start .gt. 1) then
if (named .ne. 0) then
print *,'unnamed options must appear first: ',options
endif
nopt=nopt+1
nam1(nopt)=1
nam2(nopt)=1
val1(nopt)=leng+1
call str_append(opt, leng, options(start:ende))
val2(nopt)=leng
! print *,'unnamed "',options(start:ende),'"'
endif
else
named=1
if (i .eq. 1) then
print *,'syntax error in options: ',options(start:ende)
else
nopt=nopt+1
nam1(nopt)=leng+1
call str_append(opt, leng, options(start:start+i-2))
nam2(nopt)=leng
call str_upcase(opt(nam1(nopt):leng),opt(nam1(nopt):leng))
if (start+i .le. ende) then
val1(nopt)=leng+1
call str_append(opt, leng, options(start+i:ende))
val2(nopt)=leng
! print *,'named "',options(start:ende),'"'
else
val1(nopt)=1
val2(nopt)=1
! print *,'named empty "',options(start:ende),'"'
endif
endif
endif
endif
enddo
if (leng .eq. len(opt)) then
print *,'too many options, truncated'
endif
return
entry dat_start_options
if (mopt .eq. 0) mopt=nopt
do i=1,nopt
cnt(i)=0
enddo
idx=0
return
entry dat_str_option(name, value) ! for dat_xxx
found=0
idx=idx+1
if (idx .le. mopt .and. nam2(idx) .eq. 1) then
value=opt(val1(idx):val2(idx))
cnt(idx)=cnt(idx)+1
found=1
endif
call str_upcase(nam, name)
do i=1,nopt
if (opt(nam1(i):nam2(i)) .eq. nam) then
if (found .ne. 0) then
if (value .ne. opt(val1(i):val2(i)) .and. i .le. mopt) then
cnt(i)=cnt(i)+1
endif
RETURN
endif
value=opt(val1(i):val2(i))
cnt(i)=cnt(i)+1
found=1
endif
enddo
return
entry dat_end_options
do i=1,mopt
if (cnt(i) .eq. 0) then
if (nam2(i) .eq. 1) then
print *,'superflous option: ', opt(val1(i):val2(i))
else
print *,'unknown option: '
1 ,opt(nam1(i):nam2(i)), '=', opt(val1(i):val2(i))
endif
elseif (cnt(i) .eq. 2) then
print *,'ambigous option: '
1 ,opt(nam1(i):nam2(i)), '=', opt(val1(i):val2(i))
endif
enddo
return
end
subroutine dat_set_index(idx)
integer idx, flg
integer indx/0/, flag/0/
indx=idx
flag=1
return
entry dat_get_index(idx) ! for dat_xxx, call only when index is supported
idx=indx
flag=2
return
entry dat_next_index(flg)
if (indx .gt. 0) then
indx=indx+1
flg=flag
else
flg=flag
endif
end
subroutine dat_real_option(name, value)
character name*(*)
integer ivalue
real value
character str*64
integer l
str=' '
call dat_str_option(name, str)
if (str .eq. ' ') RETURN
value=0.0
read(str, *, err=9, end=9) value
RETURN
9 call str_trim(str, str, l)
print *,'option ',name,'=',str(1:l),' must be a real number'
RETURN
entry dat_int_option(name, ivalue)
str=' '
call dat_str_option(name, str)
if (str .eq. ' ') RETURN
ivalue=0
read(str, *, err=19, end=19) ivalue
RETURN
19 call str_trim(str, str, l)
print *,'option ',name,'=',str(1:l),' must be an integer'
RETURN
end
subroutine dat_use_mon(imon)
!
! convention: 0: auto, 1: standard, 2: primary beam, 3: auxiliary, 4: time, 5: proton current/reactor power
!
integer imon
integer used_mon/1/
used_mon=imon
return
entry dat_used_mon(imon)
imon=used_mon
end
subroutine dat_use_axis(xaxis, yaxis)
character xaxis*(*), yaxis*(*)
character xax*64/' '/, yax*64/' '/
xax=xaxis
yax=yaxis
entry dat_used_axis(xaxis, yaxis)
xaxis=xax
yaxis=yax
end
subroutine dat_ask_filelist(filelist, text)
character filelist*(*), text*(*)
character spec*32
integer high, ls
call dat_getdef(spec)
call str_trim(spec, spec, ls)
high=0
call dat_get_high(high)
20 print *
if (ls .ne. 0) then
if (high .eq. 0) then
print '(X,3A)',
1 'Default for numors: ', spec(1:ls)
else
print '(X,3A,i6,a)',
1 'Default for numors: ', spec(1:ls),' (highest: ',high,')'
endif
endif
c call dat_gettyp(spec)
c if (spec .ne. ' ') then
c call str_trim(spec, spec, ls)
c print '(x,3a)', 'Default type for files (', spec(1:ls), ')'
c endif
if (text .ne.' ') print *,text
print '(/X,A,$)', 'Filename(s) or numor(s) (? to get help): '
read(*, '(a)',err=91,end=91) filelist
if (filelist .eq. '?') then
print '(X,A)'
1, 'Valid examples:'
1, ' '
1, '1001 a numor of the default instrument'
1, 'TASP/1997/1234 numor 1234 of year 1997 of TASP'
1, 'IN3/923-925,927 numors 923,924,925 and 927 of IN3'
1, 'TEST.FIT3 or any other file specification'
goto 20
endif
91 end
subroutine dat_powder_trf(xin,xout,yfact)
real xin
real xout
real yfact
real lambda_i
character xaxis*(*), oaxis*(*)
real pi
parameter (pi=3.14159265)
real lambda/0.0/
integer mode/0/ ! 1: 2theta to q, 2: 2theta to d
real sint
character up*16, oup*16
real sind, asind, x
sind(x)=sin(x*3.14159265/180.0)
asind(x)=asin(x)*180./3.14159265
yfact=1.0 ! the feature to change y was not appreciated
if (mode .eq. 0) then
xout=xin
RETURN
endif
if (mode .le. 2) then ! 2t -> q
xout=max(1.0e-6,min(xin,179.999))
if (mode .eq. 1) then
xout=(4*pi)*sind(xin/2)/lambda
! yfact=1.0/cosd(xin/2) ! removed
elseif (mode .eq. 2) then ! 2t -> d
sint=sind(xin/2)
if (sint .eq. 0) sint=1e-3
xout=lambda/2/sint
! yfact=tand(xin/2)*sint ! removed
endif
elseif (mode .eq. 3) then ! q -> 2t
xout=2*asind(max(-1.0,min(xin*lambda/(4*pi), 1.0)))
elseif (mode .eq. 6) then ! d -> 2t
xout=2*asind(min(1.0,lambda/2/max(1e-6,xin)))
elseif (mode .eq. 5 .or. mode .eq. 7) then ! q <-> d
xout=(2*pi)/max(1e-6,xin)
endif
RETURN
entry dat_powder_init(lambda_i, oaxis, xaxis)
mode=0
call str_upcase(oup, oaxis)
call str_upcase(up, xaxis)
if (oup .eq. 'Q') then
mode=3
else if (oup .eq. 'D') then
mode=6
else if (oup .ne. ' ' .and. oup(1:1) .ne. '2'
1 .and. oup(1:1) .ne. 'T') then
print *,'cannot transform from ',oaxis
xaxis=' '
goto 9
endif
if (up .eq. 'Q') then
mode=mode+1
elseif (up .eq. 'D') then
mode=mode+2
elseif (up .eq. ' ' .or. up(1:1) .eq. '2') then
xaxis='2-Theta'
else
print *,'cannot transform to ',xaxis
xaxis=' '
endif
if (lambda_i .eq. 0 .and. mode .gt. 0 .and.
1 mode .le. 3 .or. mode .eq. 6) then
print *,'no wavelength lambda defined'
xaxis=' '
mode=0
endif
9 lambda=lambda_i
end
subroutine dat_get_datanumber(file, numor)
character file*(*)
integer numor
integer luntmp, iostat
call sys_get_lun(luntmp)
call sys_open(luntmp, file, 'r', iostat)
if (iostat .eq. 0) then
read(luntmp, *, err=18,end=18) numor
18 close(luntmp)
endif
call sys_free_lun(luntmp)
end