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

580 lines
12 KiB
Fortran

subroutine fit_file(calstep, filename)
! --------------------------------------
implicit none
include 'fit.inc'
real calstep
character filename*(*)
real step, s, xx
character fname*132
integer l,i,jset,j,iostat
real xm1, xm2, a1, a2
real fifun, voigt ! function
if (iscx .eq. 0) then
xm1=xval(nxmin)
xm2=xval(nxmax)
do i=nxmin,nxmax
if (xval(i) .lt. xm1) xm1=xval(i)
if (xval(i) .gt. xm2) xm2=xval(i)
enddo
a2=0.5/(nxmax+2-nxmin)
a1=1+a2
xbeg=xm1*a1-xm2*a2
xend=xm2*a1-xm1*a2
endif
step=(xend-xbeg)/200.0
if (filename .eq. ' ') then
write (isyswr,'(/X,A,$)') 'Output file name:'
read (isysrd, '(A)', end=999,err=999) fname
if (fname .eq. ' ') fname='out'
if (nu .gt. 0) then
write (isyswr,'(/X,A,F10.3,A,$)')
1 'X-Step for calculated data (',abs(step),'):'
read (isysrd, '(F20.0)', err=999,end=999) s
if (s .ne. 0) step=sign(s,step)
endif
else
fname=filename
if (calstep .ne. 0) then
step=sign(calstep,step)
endif
endif
call save_delimiter(char(9))
do jset=1,nset
if (nu .eq. 0 .and. nset .le. 1) then ! no calculation, single dataset
call sys_parse(fname, l, fname, '.obs', 0)
else
call sys_parse(fname, l, '.obs', fname, 0)
if (nset .gt. 1 .and. l .lt. len(fname)) then
call cvt_real_str(fname(l+1:), j, float(jset), 0,0,6,3)
l=l+j
endif
endif
call sys_open(2, fname, 'w', iostat)
if (iostat .ne. 0) then
print *,'Cannot open ',fname(1:l)
return
endif
write(isyswr, '(/x,2a)') 'Observed data: ', fname(1:l)
do i=nxmin,nxmax
if (iset(i) .eq. jset) then
call save_fill(xval(i))
call save_fill(YVAL(i))
call save_fill(sig(i))
call save_wrt(2)
endif
enddo
close(2)
enddo
if (nu .le. 0) return
call sys_parse(fname, l, '.cal', fname, 0)
call sys_open(2, fname, 'w', iostat)
if (iostat .ne. 0) then
print *,'Can not open ',fname(1:l)
return
endif
write(isyswr, '(x,2a)') 'Calculated data: ', fname(1:l)
xx=xbeg
do while (xx .le. xend+step*1e-3)
call save_fill(xx)
actset=1
call save_fill(fifun(xx))
if (ififu .eq. 1) then
call save_fill(u(1)+u(2)*(xx-u(3))) ! Background
do i=3,nu,5
call save_fill(u(1)+u(2)*(xx-u(3)) ! Peaks
1 +voigt(xx-u(i), u(i+3), u(i+4))*u(i+2))
enddo
elseif (ififu .eq. 7 .and. nset .gt. 1) then
do actset=2,nset
call save_fill(fifun(xx))
enddo
endif
call save_wrt(2)
xx=xx+step
enddo
close(2)
return
999 print *,'input error'
end
subroutine fit_save(filename)
! -----------------------------
implicit none
include 'fit.inc'
character filename*(*)
integer lunit
character fname*132/' '/, line*132, filtyp*8
character utit*16, parname*16
integer i,l,j,nx1,nx2,ispec,lun,iostat,nu_user
real sbuf(32)
integer fit_userini
external fit_userini
external fit_wrapper
call save_delimiter(',')
if (fname .le. ' ') then
fname='last.fit3'
endif
if (filename .ne. ' ') then
line=filename
else
call sys_parse(fname, l, fname, '.FIT3', 0)
write(isyswr, '(/X,3A,$)') 'Filename [',fname(1:l),']: '
read(isysrd, '(a)', end=999,err=999) line
endif
call sys_parse(line, l, line, fname, 0)
lun=2
call sys_parse(line, l, line, ' ', 0)
call sys_open(lun, line, 'w', iostat)
if (iostat .ne. 0) then
print *,'can not open ',line(1:l)
return
endif
write(isyswr, '(/1x,2a)') 'Save parameters and data on '
1 ,line(1:l)
write(2,'(A)') 'FitSave 3.8'
call save_fill(float(nu))
call save_fill(float(ififu))
call dat_gettyp(filtyp)
call save_str(filtyp)
call save_fill(float(nxmin))
call save_fill(float(nxmax))
if (iscx .eq. 0) then
xbeg=1
xend=1
endif
if (iscy .eq. 0) then
ybeg=1
yend=1
endif
call save_fill(xbeg)
call save_fill(xend)
call save_fill(ybeg)
call save_fill(yend)
if (nu .gt. 0) then
call save_fill(amin)
else
call save_fill(0.0)
endif
call save_wrt(lun)
do i=1,nu
if (psho(i) .ne. ' ') then
call str_trim(psho(i),psho(i),l)
call save_str(psho(i)(1:l)//':'//pnam(i))
else
call save_str(pnam(i))
endif
call save_fill(u(i))
call save_fill(werr(i))
call save_fill(alim(i))
call save_fill(blim(i))
call save_fill(float(lcode(i)))
call save_fill(cfac(i))
call save_fill(coff(i))
call save_fill(float(icsw(i)))
call save_fill(float(icto(i)))
call save_fill(werrs(i))
call save_wrt(lun)
enddo
if (ififu .eq. 7) then
write(lun,'(a)') usertit
endif
call sym_put_real('Monitor', ymon)
call sym_list(lun, 0, 3, ' ')
do i=1,min(nstyl,32)
sbuf(i)=styl(i)
enddo
call save_array('Style', sbuf, nstyl)
call save_wrt(lun)
if (keepregion) then
call save_array('RegX1', regx1, nregion)
call save_wrt(lun)
call save_array('RegX2', regx2, nregion)
call save_wrt(lun)
call save_array('RegY1', regy1, nregion)
call save_wrt(lun)
call save_array('RegY2', regy2, nregion)
call save_wrt(lun)
endif
call str_trim(fillis, fillis, l)
write(lun,'(3a)') 'Filelist=''',fillis(1:l),''''
write(lun,*)
do i=1,npkt
call save_fill(xval(i))
call save_fill(YVAL(i))
call save_fill(sig(i))
call save_fill(rmon(i))
call save_fill(float(iset(i)))
call save_wrt(lun)
enddo
close(lun)
return
999 print *,'input error'
return
entry fit_load_nopn(lunit)
! --------------------------
lun=lunit
rewind lun
inquire(lun, name=line)
call sys_parse(line, l, line, '.', 0)
goto 20
entry fit_load(filename)
! ------------------------
call sys_parse(line, l, filename, '.FIT3', 0)
lun=2
call sys_open(lun, line(1:l), 'r', iostat)
if (iostat .ne. 0) goto 29
20 write(isyswr, '(/x,2a)') 'Load parameters and data from '
1 ,line(1:l)
line=' '
read(lun, '(A)',err=23,end=23) line
23 if (line(1:8) .ne. 'FitSave ') then
write(isyswr,*) 'Unknown file format'
close(lun)
return
endif
call sym_purge(1)
line(1:1)=' '
if (line(9:11) .gt. '1.0' .and. line(9:11) .lt. '3.4') then
read(lun, '(A)',err=23,end=23) filnam
endif
nu=0
ififu=0
xbeg=1
xend=1
ybeg=1
yend=1
filtyp=' '
if (line(9:11) .lt. '3.2') then
read(lun, *, err=23,end=23) nu,ififu,ispec,nx1,nx2
1,xbeg,xend,ybeg,yend
if (ispec .eq. 1) then
filtyp='IN3'
elseif (ispec .ne. 0) then
filtyp='LNS'
endif
else
read(lun, *, err=23,end=23) nu,ififu,filtyp,nx1,nx2
1 ,xbeg,xend,ybeg,yend
endif
if (filtyp .ne. ' ') call dat_settyp(filtyp)
if (line(9:11) .lt. '3.0') then
ififu=ififu-1
if (ififu .le. 4 .or. ififu .eq. 9) then
write(isyswr,*) 'Incompatible version (older than 3.0)'
nu=0
ififu=8
close(lun)
return
endif
endif
do i=1,nu
u(i)=0
werr(i)=0
alim(i)=0
blim(i)=0
lcode(i)=0
cfac(i)=0
coff(i)=0.0
icsw(i)=0
icto(i)=0
werrs(i)=0
if (line(9:11) .ge. '3.8') then
read(lun,*,err=27,end=27) parname, u(i), werr(i)
1 , alim(i), blim(i), lcode(i), cfac(i), coff(i)
1 , icsw(i), icto(i), werrs(i)
else
read(lun,*,err=27,end=27) parname, u(i), werr(i)
1 , alim(i), blim(i), lcode(i), cfac(i)
1 , icsw(i), icto(i), werrs(i)
endif
j=index(parname,':')
pnam(i)=parname(j+1:)
if (j .gt. 1) then
psho(i)=parname(1:j-1)
else
psho(i)=' '
endif
enddo
if (ififu .eq. 7 .and. line(9:11) .ge. '3.5') then
read(lun,'(a)',err=23,end=23) utit
if (line(9:11) .lt. '3.7') read(lun,*,err=23,end=23)
if (usernp .lt. nu .or. utit .ne. usertit) goto 25
do i=1,nu
if (pnam(i) .ne. userpar(i)(1:8)) goto 25
enddo
goto 26
25 do i=nu+1,usernp
pnam(i)=' '
enddo
do i=usernp+1,nu
userpar(i)=' '
enddo
write(isyswr,*)
write(isyswr,*) 'Present function data file'
write(isyswr,*) '---------------------------'
write(isyswr,*) usertit,' ',utit
write(isyswr,*)
do i=1,max(nu,usernp)
write(isyswr,'(x,a,10x,a)') userpar(i)(1:8),pnam(i)
pnam(i)=userpar(i)
psho(i)=usersho(i)
enddo
nu=usernp
write(isyswr,'(/x,a/)')
1 'User function mismatch, check Parameters for validity'
26 continue
endif
if (line(9:11) .lt. '3.3') then
close(lun)
if (filnam .eq. ' ') then
write(isyswr,*) 'No datafile'
else
load_state=3
endif
else
if (line(9:11) .lt. '3.4') then
j=0
ymon=0
read(lun, *, err=27,end=27) j,ymon
else
call sym_read(lun, fit_wrapper)
call fit_get_real('Monitor', ymon)
keepregion=(nregion .ne. 0)
ymon0=ymon
endif
i=0
30 i=i+1
if (i .gt. maxdat) then
print *,'too many data points --> truncated'
npkt=maxdat
goto 39
endif
31 read(lun,'(a)',err=27,end=39) line
j=index(line,'/') ! for compatibility with versions 3.3 and older
if (j .gt. 0) line(j:)=' ' ! "
! read(line,'(bn,4f20.0,i20)',err=37,end=37)
! 1 xval(i),YVAL(i),sig(i),rmon(i),j
read(line,*,err=37,end=37) xval(i),YVAL(i),sig(i),rmon(i),j
iset(i)=max(1,j)
npkt=i
goto 30
37 print *,'error at point ',i
goto 31
39 if (ymon .eq. 0) ymon=rmon(1)
close(lun)
endif
wavlen=0
call sym_get_real('lambda', wavlen)
temp=0
call sym_get_real('Temp', temp)
dtemp=0
call sym_get_real('dTemp', dtemp)
call sym_get_str('Title', l, itit)
if (nx1 .lt. npkt .and. nx2 .le. npkt .and. nx2 .gt. nx1) then
nxmin=nx1
nxmax=nx2
else
nxmin=1
nxmax=npkt
endif
nset=0
do i=nxmin,nxmax
nset=max(iset(i),nset)
enddo
if (ififu .eq. 7) then
nu_user=fit_userini(7)
endif
nfcn=0
call sym_list(isyswr, 1, 1, ' ')
if (load_state .eq. 3) return
call fit_set(0,0.0,0.0,0.0,0.0)
call fit_print(1)
call fit_scale(xbeg,xend,ybeg,yend)
return
27 write(isyswr,*) 'Error in ',filename
close(lun)
return
29 write(isyswr, *) 'Can not open ',filename
end
subroutine save_fill(val)
implicit none
real val, values(*)
integer lun, nsize
character str*(*), delim*1
character line*132
integer i,j,l/0/,ll
save line,delim,l,ll
if (delim .eq. char(9)) then
ll=7
else
ll=12
endif
if (l+ll .ge. len(line)) return
call cvt_real_str(line(l+1:l+ll), j, val, 0,0,6,3)
l=l+j
line(l+1:l+1)=delim
l=l+1
return
entry save_str(str)
if (l+len(str)+3 .ge. len(line)) return
line(l+1:l+1)=''''
l=l+1
call str_trim(line(l+1:l+len(str)), str, j)
if (j .eq. 0) j=1
l=l+j
line(l+1:l+1)=''''
l=l+1
line(l+1:l+1)=delim
l=l+1
return
entry save_array(str, values, nsize)
if (l+len(str)+1 .ge. len(line)) goto 39
call str_trim(line(l+1:l+len(str)), str, j)
if (j .eq. 0) j=1
l=l+j
line(l+1:l+1)='='
l=l+1
do i=1,min(32,nsize)
if (l+9 .ge. len(line)) goto 39
call cvt_real_str(line(l+1:l+7), j, values(i), 0,0,6,3)
l=l+j
line(l+1:l+1)=' '
l=l+1
enddo
line(l:l)=delim
return
39 print *,str,' truncated'
return
entry save_wrt(lun)
l=l-1
write(lun, '(a)') line(1:l)
l=0
return
entry save_delimiter(str)
delim=str
end
subroutine fit_wrapper(str, value, putval)
include 'fit.inc'
character str*(*)
real value
external putval
integer ns/0/, nx1/0/, nx2/0/, ny1/0/, ny2/0/
if (str .eq. ' ') then ! reset
ns=0
nx1=0
nx2=0
ny1=0
ny2=0
return
endif
if (str(1:min(len(str),9)) .eq. 'Filelist=') then
fillis=str(10:)
return
endif
if (str .eq. 'Style') then
if (ns .lt. maxset) then
ns=ns+1
styl(ns)=nint(value)
endif
nstyl=ns
return
endif
if (len(str) .eq. 5 .and. str(1:min(len(str),3)) .eq. 'Reg') then
if (str .eq. 'RegX1') then
if (nx1 .lt. maxregion) then
nx1=nx1+1
regx1(nx1)=value
endif
else if (str .eq. 'RegX2') then
if (nx2 .lt. maxregion) then
nx2=nx2+1
regx2(nx2)=value
endif
else if (str .eq. 'RegY1') then
if (ny1 .lt. maxregion) then
ny1=ny1+1
regy1(ny1)=value
endif
else if (str .eq. 'RegY2') then
if (ny2 .lt. maxregion) then
ny2=ny2+1
regy2(ny2)=value
endif
else
goto 9
endif
nregion=min(nx1,nx2,ny1,ny2)
return
endif
9 call putval(str, value)
end