719 lines
15 KiB
Fortran
719 lines
15 KiB
Fortran
subroutine fit_dat(filenames)
|
|
|
|
character*(*) filenames
|
|
|
|
integer iret, fit_dat_opt
|
|
|
|
iret=fit_dat_opt(filenames, 0, 0.0) ! option 0: standard
|
|
end
|
|
|
|
|
|
subroutine fit_link(filenames)
|
|
|
|
character*(*) filenames
|
|
|
|
integer iret, fit_dat_opt
|
|
|
|
iret=fit_dat_opt(filenames, 1, 0.0) ! option 1: link to existing data
|
|
end
|
|
|
|
|
|
subroutine fit_dat_merge(filenames, step)
|
|
|
|
character*(*) filenames
|
|
real step
|
|
|
|
integer iret, fit_dat_opt
|
|
|
|
iret=fit_dat_opt(filenames, 2, step) ! option 2: merge read data
|
|
end
|
|
|
|
|
|
subroutine fit_open(filenames)
|
|
|
|
character*(*) filenames
|
|
|
|
integer iret, fit_dat_opt
|
|
|
|
iret=fit_dat_opt(filenames, 4, 0.0) ! option 4: <return> reads last.fit3
|
|
end
|
|
|
|
|
|
subroutine fit_dat_next(filenames, flag)
|
|
|
|
character*(*) filenames
|
|
integer flag ! input: 0: start read, 1: read next file, 2: read next range
|
|
! output: 0: end of list reached, 1: o.k., no range, 2: o.k., next range
|
|
|
|
|
|
call fit_dat_next_opt(filenames, flag, 0, 0.0)
|
|
end
|
|
|
|
|
|
subroutine fit_dat_next_opt(filenames, flag, opt, step)
|
|
|
|
character*(*) filenames
|
|
integer flag ! input: 0: start read, 1: read next file, 2: read next range
|
|
! output: 0: end of list reached, 1: o.k., no range, 2: o.k., next range
|
|
integer opt ! 0: DAT, 1: LINK, 2: DAT&MERGE
|
|
real step
|
|
|
|
integer lopt, fit_dat_opt
|
|
integer iret
|
|
integer idx/0/
|
|
|
|
logical silent, silent0
|
|
common/fit_dat_com/silent
|
|
|
|
if (flag .eq. 0) then
|
|
call dat_set_index(1)
|
|
lopt=8 ! read 1st
|
|
else if (flag .eq. 1) then
|
|
call dat_set_index(1)
|
|
lopt=16 ! read next file
|
|
else
|
|
idx=idx+1
|
|
lopt=24 ! read again
|
|
endif
|
|
silent0=silent
|
|
iret=fit_dat_opt(filenames, lopt+opt, step)
|
|
if (iret .eq. 0 .and. flag .eq. 2) then
|
|
call dat_set_index(1)
|
|
silent=silent0
|
|
iret=fit_dat_opt(filenames, 16+opt, step) ! read next file
|
|
endif
|
|
if (iret .eq. 0) then
|
|
flag=0
|
|
else
|
|
call dat_next_index(flag)
|
|
endif
|
|
end
|
|
|
|
|
|
integer function fit_dat_opt(filenames, opt, step)
|
|
! --------------------------------------------------
|
|
|
|
implicit none
|
|
|
|
include 'fit.inc'
|
|
|
|
character filenames*(*)
|
|
integer opt
|
|
real step
|
|
|
|
integer n_mul
|
|
parameter (n_mul=256)
|
|
character filelist*(maxflen), sample*80
|
|
integer i,j,l,ll,nread,n1
|
|
integer merge, list_mode
|
|
integer fidx, sets/0/
|
|
logical init, init_dat/.true./, purged
|
|
integer nmax, ncol, nset0
|
|
external fit_putval, fit_fit3_read
|
|
real stp, t1, dt1, w1, tmean, wmean, ds, upar
|
|
|
|
logical silent
|
|
common/fit_dat_com/silent
|
|
data silent/.false./
|
|
|
|
! setup different flags
|
|
init=mod(opt,2) .eq. 0 ! bit0=0: purge old data
|
|
merge=mod(opt,4)/2 ! bit1=1: merge data ?
|
|
if (merge .ne. 0) stp=step ! step for merge
|
|
load_state=mod(opt,8)/4 ! bit2: do a load for FitSave files
|
|
list_mode=mod(opt,32)/8 ! bit3+4: 1: read first file, 2: read next file (ignore filenames argument), 3: read again
|
|
|
|
fit_dat_opt=1
|
|
if (list_mode .gt. 0) then
|
|
if (list_mode .ne. 3) then
|
|
call dat_init_list
|
|
if (list_mode .eq. 1) then
|
|
call dat_filelist(filenames)
|
|
endif
|
|
endif
|
|
goto 15
|
|
endif
|
|
fidx=1
|
|
11 if (filenames(fidx:fidx) .eq. '+') then
|
|
fidx=fidx+1
|
|
if (load_state .eq. 0) init=.false.
|
|
goto 11
|
|
elseif (filenames(fidx:fidx) .eq. '&') then
|
|
fidx=fidx+1
|
|
merge=1
|
|
goto 11
|
|
endif
|
|
|
|
if (init_dat) then
|
|
call dat_init
|
|
call dat_fit3_replace(fit_fit3_read)
|
|
init_dat=.false.
|
|
endif
|
|
if (filenames(fidx:) .eq. ' ') then
|
|
if (load_state .eq. 1) then
|
|
call dat_ask_filelist(filelist
|
|
1 , 'Hit <RETURN> to continue with last data')
|
|
else
|
|
call dat_ask_filelist(filelist, ' ')
|
|
endif
|
|
if (filelist .eq. ' ') then
|
|
if (load_state .eq. 1) then
|
|
filelist='last.fit3'
|
|
else
|
|
fit_dat_opt=0
|
|
goto 91
|
|
endif
|
|
endif
|
|
fidx=1
|
|
12 if (filelist(fidx:fidx) .eq. '+') then
|
|
fidx=fidx+1
|
|
if (load_state .eq. 0) init=.false.
|
|
goto 12
|
|
elseif (filelist(fidx:fidx) .eq. '&') then
|
|
fidx=fidx+1
|
|
merge=1
|
|
goto 12
|
|
endif
|
|
call str_trim(filelist,filelist(fidx:),ll)
|
|
else
|
|
call str_trim(filelist,filenames(fidx:),ll)
|
|
endif
|
|
|
|
|
|
if (init) then
|
|
sets=0
|
|
call dat_init_list
|
|
else
|
|
call dat_put_filelist(fillis)
|
|
endif
|
|
|
|
call dat_filelist(filelist(1:ll))
|
|
filesave=filelist(1:ll)
|
|
|
|
15 if (init) then
|
|
npkt=0
|
|
call sym_purge(1)
|
|
temp=0
|
|
dtemp=0
|
|
wtemp=0
|
|
endif
|
|
|
|
purged=.false.
|
|
n1=npkt+1
|
|
30 if (silent) call dat_silent
|
|
call fit_dat_init_table
|
|
call fit_ncurves(0)
|
|
if (list_mode .ne. 3) then
|
|
call dat_read_next(fit_putval, maxdat-npkt, nread
|
|
1 , xval(npkt+1), YVAL(npkt+1), sig(npkt+1), rmon(npkt+1))
|
|
else
|
|
call dat_read_again(fit_putval, maxdat-npkt, nread
|
|
1 , xval(npkt+1), YVAL(npkt+1), sig(npkt+1), rmon(npkt+1))
|
|
endif
|
|
if (nread .gt. 0) then
|
|
if (npkt .ne. 0) then
|
|
if (ymon0 .ne. ymon) then
|
|
if (ymon .eq. 0) then
|
|
ymon=ymon0
|
|
elseif (ymon0 .eq. 0) then
|
|
if (ymon .ne. 0) then
|
|
do i=npkt+1,npkt+nread
|
|
rmon(i)=ymon
|
|
enddo
|
|
endif
|
|
else
|
|
do i=npkt+1,npkt+nread
|
|
YVAL(i)=YVAL(i)*ymon/ymon0
|
|
sig(i)=sig(i)*ymon/ymon0
|
|
enddo
|
|
endif
|
|
endif
|
|
else
|
|
nset=0
|
|
ymon=ymon0
|
|
endif
|
|
|
|
w1=0
|
|
j=1
|
|
nset0=nset
|
|
nmax=npkt+nread
|
|
do while (npkt .lt. nmax)
|
|
ncol=cols(j)
|
|
do i=0,min(nmax-npkt,rows(j)*ncol-1)
|
|
if (npkt .lt. nmax) then
|
|
npkt=npkt+1
|
|
iset(npkt)=nset+mod(i,ncol)+1
|
|
w1=w1+rmon(npkt) ! weight for temperature stddev
|
|
endif
|
|
enddo
|
|
nset=nset+ncol
|
|
j=min(nmult,j+1)
|
|
enddo
|
|
t1=0
|
|
dt1=0
|
|
call sym_get_real('Temp', t1)
|
|
if (t1 .eq. 0) then
|
|
call sym_get_real('sTemp', t1) ! this is for old DMC/HRPT files
|
|
endif
|
|
call sym_get_real('dTemp', dt1)
|
|
|
|
if (t1 .ne. 0 .and. w1 .ne. 0) then ! calculate merged temperature stddev
|
|
wmean=wtemp+w1
|
|
tmean=(temp*wtemp+t1*w1)/wmean
|
|
dtemp=sqrt((wtemp*(dtemp**2+(temp-tmean)**2)
|
|
1 + w1*( dt1**2+ (t1-tmean)**2))/wmean )
|
|
temp=tmean
|
|
wtemp=wmean
|
|
endif
|
|
|
|
if (nset .gt. maxset) goto 31
|
|
|
|
do i=1,npd
|
|
upar=none
|
|
call sym_get_real(userpar(usernp+i), upar)
|
|
if (upar .eq. none) then
|
|
if (pdpar(i, nset) .eq. none .and.
|
|
1 userpar(usernp+i) .ne. 'Numor' .and.
|
|
1 userpar(usernp+i) .ne. 'Temp') then
|
|
print *,'parameter ',userpar(usernp+i),' not found'
|
|
pdpar(i, nset)=0.0
|
|
endif
|
|
else
|
|
if (userpar(usernp+i) .eq. 'Numor') then
|
|
if (nset0+1 .eq. nset) then
|
|
ds=0.0
|
|
call sym_get_real('Pal', ds)
|
|
else
|
|
ds=0.1
|
|
endif
|
|
do j=nset0+1,nset
|
|
pdpar(i, j)=upar+ds
|
|
ds=ds+0.1
|
|
enddo
|
|
else
|
|
do j=nset0+1,nset
|
|
pdpar(i, j)=upar
|
|
enddo
|
|
endif
|
|
endif
|
|
enddo
|
|
|
|
31 continue
|
|
ISW(1) = 0
|
|
ISW(2) = 0
|
|
ISW(3) = 0
|
|
NFCN = 1
|
|
SIGMA = 0.
|
|
|
|
if (ififu .eq. 0) ififu=8 ! set to plot only mode
|
|
c print '(x,a,i7)','Npkt=',nread
|
|
call sym_show(1)
|
|
call meta_put('Npkt', float(nread))
|
|
call meta_set_format('Npkt', 5003)
|
|
call meta_set_format('Numor', 5003)
|
|
call meta_set_format('Counts', 5003)
|
|
if (.not. silent) call sym_list(isyswr,2,1,' ')
|
|
sets=sets+1
|
|
if (merge .ne. 0 .and. sets .gt. 1) then
|
|
if (npkt .gt. n1 .and.
|
|
1 (npkt .gt. maxdat/2 .or. nset .ge. maxset)) then ! probably not enough space -> merge
|
|
nxmin=n1
|
|
nxmax=npkt
|
|
call fit_merge(stp)
|
|
merge=1
|
|
else
|
|
merge=2
|
|
endif
|
|
endif
|
|
|
|
if (keepregion) then
|
|
if (npkt .gt. n1 .and. npkt .gt. maxdat/2) then ! probably not enough space -> purge
|
|
call fit_restore_region(.true.)
|
|
purged=.true.
|
|
endif
|
|
endif
|
|
|
|
if (list_mode .eq. 0) then
|
|
if (npkt .lt. maxdat) goto 30
|
|
print *,'no more datapoints accepted'
|
|
endif
|
|
else if (list_mode .ne. 0) then
|
|
fit_dat_opt=0
|
|
endif
|
|
|
|
if (merge .eq. 2 .and. npkt .ge. n1) then ! merge (if not yet done)
|
|
nxmin=n1
|
|
nxmax=npkt
|
|
call fit_merge(stp)
|
|
endif
|
|
if (keepregion) then
|
|
call fit_restore_region(purged)
|
|
if (purged) then
|
|
print *,'points outside window/region were purged'
|
|
endif
|
|
else
|
|
call fit_win(1.0,1.0) ! disable window
|
|
endif
|
|
|
|
if (load_state .gt. 1) then
|
|
if (load_state .eq. 3) then ! for FitSave 3.2 and older
|
|
call dat_init_list
|
|
call dat_filelist(filnam)
|
|
load_state=0
|
|
goto 30
|
|
endif
|
|
else
|
|
call dat_get_filelist(fillis)
|
|
call str_trim(fillis, fillis, l)
|
|
if (.not. silent) then
|
|
print *
|
|
if (list_mode .eq. 0) then
|
|
if (fillis(1:l) .eq. ' ') then
|
|
print *,'no file read'
|
|
else
|
|
print *,'Files read: ',fillis(1:l)
|
|
endif
|
|
endif
|
|
endif
|
|
|
|
wavlen=0
|
|
call sym_get_real('lambda', wavlen)
|
|
call sym_get_str('Title', ll, itit)
|
|
sample=' '
|
|
call sym_get_str('sample', l, sample)
|
|
if (sample(1:l) .ne. ' '
|
|
1 .and. index(itit(1:ll),sample(1:l)) .eq. 0) then
|
|
if (ll+l .lt. len(itit)) then
|
|
if (ll+l+7 .lt. len(itit)) then
|
|
itit(ll+1:)=' Sample='
|
|
ll=ll+8
|
|
else
|
|
ll=ll+1
|
|
endif
|
|
itit(ll+1:)=sample
|
|
endif
|
|
endif
|
|
|
|
call sym_show(1)
|
|
call sym_put_real('Temp', temp)
|
|
call sym_put_real('dTemp', dtemp)
|
|
call sym_put_real('Monitor', ymon)
|
|
|
|
if (npkt .gt. 0) call fit_set(0,0.0,0.0,0.0,0.0)
|
|
endif
|
|
91 silent=.false.
|
|
end
|
|
|
|
|
|
subroutine fit_dat_silent
|
|
|
|
logical silent
|
|
common/fit_dat_com/silent
|
|
|
|
silent=.true.
|
|
end
|
|
|
|
|
|
subroutine fit_datraw
|
|
! ---------------------
|
|
print *,'Subroutine FIT_DATRAW no longer available'
|
|
end
|
|
|
|
|
|
subroutine fit_spec
|
|
! -------------------
|
|
print *,'Subroutine FIT_SPEC is no longer available'
|
|
END
|
|
|
|
|
|
|
|
subroutine fit_dat_put(mode, xx, nx, yy, ny, ss, ns, ww, nw) !!
|
|
!! ------------------------------------------------------------
|
|
|
|
! preconditions: (ny>0) and (nx=ny or nx=1 or nx=2) and (ns=ny or ns=1) and (nw=ny or nw=1)
|
|
! nx=1: x-values are xx(1),xx(1)+1,xx(1)+2,...,xx(1)+(nx-1)
|
|
! nx=2: x-values are xx(1),xx(2),...,xx(1)+(xx(2)-xx(1))*(nx-1)
|
|
! ns=1: ss(1) <> 0.0: sigma values are ss(1)
|
|
! ss(1) = 0.0: sigma values are sqrt(max(1.0,yy(i)))
|
|
! nw=1: weights are ww(1)
|
|
|
|
include 'fit.inc'
|
|
|
|
integer mode !! mode=0: purge before, mode=1: link new dataset, mode=2: link to existing dataset
|
|
integer nx, ny, ns, nw !!
|
|
real xx(nx), yy(ny), ss(ns), ww(nw) !!
|
|
|
|
integer i,j,nread
|
|
real x0,xs
|
|
|
|
if (mode .eq. 0) then
|
|
npkt=0
|
|
call sym_purge(1)
|
|
endif
|
|
|
|
if (ny .lt. 0) stop 'FIT_DAT_INT: illegal ny'
|
|
|
|
nread=ny
|
|
if (npkt+ny .gt. maxdat) then
|
|
print *,'Too many datapoints: truncated'
|
|
nread=maxdat-npkt
|
|
endif
|
|
|
|
j=npkt
|
|
do i=1,nread
|
|
j=j+1
|
|
YVAL(j)=yy(i)
|
|
enddo
|
|
|
|
j=npkt
|
|
if (nx .eq. ny) then
|
|
do i=1,nread
|
|
j=j+1
|
|
xval(j)=xx(i)
|
|
enddo
|
|
else
|
|
x0=xx(1)
|
|
if (nx .eq. 1) then
|
|
xs=1
|
|
elseif (nx .eq. 2) then
|
|
xs=(xval(nx)-x0)/(nx-1)
|
|
else
|
|
stop 'FIT_DAT_INT: illegal nx'
|
|
endif
|
|
x0=x0-xs
|
|
do i=nx+1,nread
|
|
j=j+1
|
|
xval(j)=x0+xs*i
|
|
enddo
|
|
endif
|
|
|
|
j=npkt
|
|
if (ns .eq. ny) then
|
|
do i=1,nread
|
|
j=j+1
|
|
sig(j)=ss(i)
|
|
enddo
|
|
elseif (ns .eq. 1) then
|
|
xs=ss(1)
|
|
if (xs .eq. 0.0) then
|
|
do i=1,nread
|
|
j=j+1
|
|
sig(j)=sqrt(max(1.0,YVAL(j)))
|
|
enddo
|
|
else
|
|
do i=1,nread
|
|
j=j+1
|
|
sig(j)=xs
|
|
enddo
|
|
endif
|
|
else
|
|
stop 'FIT_DAT_INT: illegal ns'
|
|
endif
|
|
|
|
j=npkt
|
|
if (nw .eq. ny) then
|
|
ymon0=ww(1)
|
|
do i=1,nread
|
|
j=j+1
|
|
rmon(j)=ww(i)
|
|
enddo
|
|
elseif (nw .eq. 1) then
|
|
ymon0=ww(1)
|
|
do i=1,nread
|
|
j=j+1
|
|
rmon(j)=ymon0
|
|
enddo
|
|
else
|
|
stop 'FIT_DAT_INT: illegal nw'
|
|
endif
|
|
|
|
if (npkt .ne. 0) then
|
|
if (mode .eq. 1) nset=nset+1
|
|
if (ymon0 .ne. ymon) then
|
|
if (ymon .eq. 0) then
|
|
ymon=ymon0
|
|
elseif (ymon0 .eq. 0) then
|
|
if (ymon .ne. 0) then
|
|
do i=npkt+1,npkt+nread
|
|
rmon(i)=ymon
|
|
enddo
|
|
endif
|
|
else
|
|
do i=npkt+1,npkt+nread
|
|
YVAL(i)=YVAL(i)*ymon/ymon0
|
|
sig(i)=sig(i)*ymon/ymon0
|
|
enddo
|
|
endif
|
|
endif
|
|
else
|
|
nset=1
|
|
ymon=ymon0
|
|
endif
|
|
do i=npkt+1,npkt+nread
|
|
iset(i)=nset
|
|
enddo
|
|
npkt=npkt+nread
|
|
|
|
ISW(1) = 0
|
|
ISW(2) = 0
|
|
ISW(3) = 0
|
|
NFCN = 1
|
|
SIGMA = 0.
|
|
|
|
nxmin=1
|
|
nxmax=npkt
|
|
end
|
|
|
|
|
|
subroutine fit_putval(str, val)
|
|
|
|
character*(*) str
|
|
real val
|
|
|
|
include 'fit.inc'
|
|
|
|
if (str .eq. 'ShowLevel') then
|
|
call sym_show(nint(val))
|
|
return
|
|
elseif (str .eq. 'Monitor' .or. str .eq. 'monitor') then
|
|
ymon0=val
|
|
endif
|
|
call meta_put(str, val)
|
|
end
|
|
|
|
|
|
subroutine fit_fit3_read(lun,forced,nread,putval,nmax,xx,yy,ss,ww)
|
|
! ------------------------------------------------------------------
|
|
implicit none
|
|
|
|
include 'fit.inc'
|
|
|
|
integer lun, forced, nread, nmax
|
|
external putval
|
|
real xx(nmax), yy(nmax), ss(nmax), ww(nmax)
|
|
|
|
character str*8
|
|
|
|
if (load_state .eq. 1) then
|
|
read(lun, '(a)') str
|
|
rewind lun
|
|
if (str .eq. 'FitSave') then
|
|
load_state=2
|
|
call fit_load_nopn(lun)
|
|
nread=0
|
|
return
|
|
endif
|
|
endif
|
|
call dat_fit3_read(lun, forced, nread, putval, nmax,xx,yy,ss,ww)
|
|
end
|
|
|
|
|
|
subroutine fit_dat_options(options)
|
|
|
|
include 'fit.inc'
|
|
|
|
character*(*) options
|
|
|
|
character typ*16, desc*256
|
|
integer done
|
|
integer l
|
|
|
|
if (options .eq. '?' .or. options .eq. ' ') then
|
|
call dat_gettyp(typ)
|
|
call str_trim(typ, typ, l)
|
|
print *
|
|
done = 1
|
|
call dat_desc_opt(done, typ(1:l))
|
|
if (done .eq. 0) then
|
|
print '(x,2a)','No options available for filetype ',typ(1:l)
|
|
else
|
|
call dat_get_def_options(desc)
|
|
if (desc .ne. ' ') then
|
|
call str_trim(desc, desc, l)
|
|
print '(x,2a)','Actual options: ',desc(1:l)
|
|
endif
|
|
if (options .eq. ' ') then
|
|
print *
|
|
print '(x,a)','Syntax: option1=value1,option2=value2 ...'
|
|
print '(x/x,a,$)','Enter default options: '
|
|
read(*,'(a)',end=9,err=9) desc
|
|
call dat_def_options(desc)
|
|
endif
|
|
endif
|
|
else
|
|
call dat_def_options(options)
|
|
endif
|
|
9 end
|
|
|
|
|
|
subroutine fit_dat_table(itbl, icols, irows)
|
|
! to be called within dat_..._read routines
|
|
! table number itbl has icols columns and irows rows
|
|
include 'fit.inc'
|
|
|
|
integer itbl, irows, icols
|
|
integer i
|
|
|
|
if (itbl .gt. nmult) then
|
|
do i=nmult+1,itbl
|
|
rows(i)=0
|
|
cols(i)=1
|
|
enddo
|
|
nmult=itbl
|
|
endif
|
|
rows(itbl)=irows
|
|
cols(itbl)=icols
|
|
end
|
|
|
|
|
|
subroutine fit_dat_init_table
|
|
|
|
include 'fit.inc'
|
|
|
|
nmult=1
|
|
rows(1)=maxdat
|
|
cols(1)=1
|
|
end
|
|
|
|
|
|
subroutine fit_dat_pdp_idx(name, ipdp)
|
|
! to be called within dat_..._read routines
|
|
include 'fit.inc'
|
|
|
|
character name*(*)
|
|
integer ipdp
|
|
|
|
integer i
|
|
character unam*32, upar*32
|
|
|
|
call str_upcase(unam, name)
|
|
do i=1,npd
|
|
call str_upcase(upar, userpar(usernp+i))
|
|
if (upar .eq. unam) then
|
|
ipdp=i
|
|
return
|
|
endif
|
|
enddo
|
|
if (usernp .lt. maxext) then
|
|
call fit_userpdp(name)
|
|
ipdp=npd
|
|
else
|
|
ipdp=0
|
|
endif
|
|
end
|
|
|
|
|
|
subroutine fit_dat_pdp_set(ipdp, idx, value)
|
|
! to be called within dat_..._read routines
|
|
|
|
include 'fit.inc'
|
|
|
|
integer ipdp, idx
|
|
real value
|
|
|
|
if (ipdp .gt. 0 .and. ipdp .le. npd .and.
|
|
1 idx .gt. 0 .and. nset+idx .le. maxset) then
|
|
pdpar(ipdp, nset+idx)=value
|
|
endif
|
|
end
|