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

201 lines
5.9 KiB
Fortran

subroutine FIT_PRINT(IKODE)
c ---------------------------
include 'fit.inc'
integer ikode
integer i,l, n1, n2, j, k, ll
real dx, al, ba, du1, du2, pos, yi, dpos, rms, xm, drms, dyi
real fact
character txt*32,pname*8,line*80
real uc
if (nu .eq. 0) return
if (nfree .le. 0) nfree=1
if (ikode .eq. 2) goto 400
if (ikode .ge. 1) call fit_confun_print
if (isyswr .ne. 0) then
write (isyswr,1000) ififu,amin,nfree,nfcn,sigma/max(1e-10,amin)
if (isw(2) .ge. 1) then
write (isyswr,1001) 'Error'
else
write (isyswr,1001) 'Step '
endif
endif
! calculate errors of corr. parameters.
do j=1,ncor
i=icord(j)
fact=cfac(i)
k=icto(i)
if (icsw(i) .lt. 0) then
if (ififu .eq. 1) then
if (u(k) .ne. 0) fact=u(i)/u(k) ! Error of Max.Int from Int.Int and vice versa
endif
endif
werr(i) = werr(k)*fact
enddo
uc=up*max(1.0,amin)/nfree
do i= 1, nu
20 l = lcorsp(i)
if (l .ne. 0) then ! variable parameter. calculate external error if v exists
if (isw(2) .ge. 1) then
dx = sqrt(abs(v(l,l)*uc))
if (lcode(i) .gt. 1) then
al = alim(i)
ba = blim(i) - al
du1 = al + 0.5 *(sin(x(l)+dx) +1.0) * ba - u(i)
du2 = al + 0.5 *(sin(x(l)-dx) +1.0) * ba - u(i)
if (dx .gt. 1.0) du1 = ba
dx = 0.5 * (abs(du1) + abs(du2))
endif
werr(i) = dx
endif
if (isyswr .ne. 0) then
call cvt_real_str(line( 1:14), ll, u(i), 14, 3, 6, 0)
call cvt_real_str(line(15:28), ll, werr(i), 14, 3, 6, 0)
call cvt_real_str(line(29:42), ll, x(l), 14, 3, 6, 0)
call cvt_real_str(line(43:56), ll, dirin(l),14, 3, 6, 0)
write (isyswr,1002) i,psho(i),pnam(i),line(1:56)
1002 FORMAT (1X,I3,2X,A4,1X,A8,2X,a)
if (lcode(i) .gt. 1) then
if (abs(cos(x(l))) .lt. 0.001) then
write(isyswr,1004) alim(i), blim(i)
1 , ' --- Parameter is at limit ---'
else
write(isyswr,1004) alim(i), blim(i)
endif
endif
endif
elseif (ikode .gt. 0) then ! fixed / corr. parameter. print only if ikode .gt.0
if (icsw(i) .ne. 0) then
pname=psho(icto(i))
if (pname .eq. ' ') pname=pnam(icto(i))
if (icsw(i) .gt. 0) then
if (coff(i) .eq. 0.0) then
write(txt, '(A,f10.3,2a)',err=390)
1 '=', cfac(i), '*', pname
else
write(txt, '(A,f10.3,3a,g12.5)',err=390)
1 '=', cfac(i), '*', pname
1 ,'+', coff(i)
endif
else
werr(i)=0
write(txt, '(6A)') 'correlated to ',pname(1:3)
1 ,', G',pname(2:3),' and L',pname(2:3)
endif
else
txt='fixed'
endif
390 if (isyswr .ne.0) then
call cvt_real_str(line( 1:14), ll, u(i), 14, 3, 6, 0)
write (isyswr,1003) i,psho(i),pnam(i),line(1:14),txt
1003 FORMAT (1X,I3,2X,A4,1X,A8,2X,a,3x,a)
endif
endif
if (ififu .eq. 1 .and. mod(i,5) .eq. 2 .and. isyswr .ne.0)
1 write(isyswr, *)
enddo
if (ikode.lt.1 .and.isw(2).lt.1 .or. isyswr .eq. 0) return
400 YINTEG = 0.
DYINTEG = 0.
ni=0
IF (IFIFU .EQ. 6) THEN
if (nxmin .ge. nxmax) return
call fit_sort(nxmin,nxmax)
N1=NXMIN+1
do while (abs(xval(n1)-u(3)) .gt. u(4)/2
1 .and. u(3) .gt. xval(n1) .and. n1 .lt. nxmax)
n1=n1+1
enddo
N2=NXMAX-1
do while (abs(xval(n2)-u(3)) .gt. u(4)/2
1 .and. u(3) .lt. xval(n2) .and. n2 .gt. nxmin)
n2=n2-1
enddo
if (n2 .le. n1) return
! calculate center of gravity
POS=0
DO 401,I=N1,N2
DX = ABS( XVAL(I+1)-XVAL(I-1) )/2. ! weight of point I
YI = (YVAL(I)-U(1)-U(2)*(XVAL(I)-U(3)))*DX ! intensity area (minus background)
YINTEG=YINTEG+YI
DYI=DX*DX*(SIG(I)**2+WERR(1)*WERR(1)+(WERR(2)*XVAL(I))**2) ! error^2 of intensity
DYINTEG=DYINTEG+DYI
POS=POS+(XVAL(I)-U(3))*YI ! sum of positions
401 CONTINUE
DYINTEG=SQRT(DYINTEG)
if (YINTEG .eq. 0) YINTEG=1
POS=POS/YINTEG+U(3) ! mean position
! calculate mean square deviation and error of center of gravity
DPOS=0
RMS=0
DO 402,I=N1,N2
XM = XVAL(I)-POS
DX = ABS( XVAL(I+1)-XVAL(I-1) )/2.
YI = (YVAL(I)-U(1)-U(2)*(XVAL(I)-U(3)))*DX
DYI=DX*DX*(SIG(I)**2+WERR(1)*WERR(1)+(WERR(2)*XVAL(I))**2)
DPOS=DPOS+XM*XM*DYI ! sum for error of position
RMS=RMS+XM*XM*YI ! sum of deviation^2 form mean pos.
402 CONTINUE
RMS=ABS(RMS/YINTEG)
! calculate error of mean square deviation
DRMS=0
DO 403,I=N1,N2
DX = ABS( XVAL(I+1)-XVAL(I-1) )/2.
DYI=DX*DX*(SIG(I)**2+WERR(1)*WERR(1)+(WERR(2)*XVAL(I))**2)
DRMS=DRMS+((XVAL(I)-POS)**2-RMS)*DYI ! sum for error of root mean square deviation
403 CONTINUE
U(5)=POS
U(6)=SQRT(RMS*8*LOG(2.0))
WERR(5)=SQRT(ABS(DPOS))/YINTEG
WERR(6)=SQRT(ABS(DRMS*8*LOG(2.0)))/YINTEG
ni=3
u(7)=YINTEG
werr(7)=dYINTEG
if (ikode .eq. 2) return
WRITE (ISYSWR,'(5X,A,5X,2F14.5)') 'Mean Posit.',u(5),werr(5)
WRITE (ISYSWR,'(5X,A,5X,2F14.5)') 'Mean Width ',U(6),WERR(6)
elseif (ififu .eq. 1) then
DO 409 I=NXMIN+1,NXMAX-1
XM = XVAL(I)-U(3)
DX = ABS( XVAL(I+1)-XVAL(I-1) )/2.
YI = (YVAL(I)-U(1)-U(2)*XM)*DX
YINTEG=YINTEG+YI
DYI=DX*DX*(SIG(I)**2+WERR(1)*WERR(1)+(WERR(2)*XM)**2)
DYINTEG=DYINTEG+DYI
409 CONTINUE
DYINTEG=SQRT(DYINTEG)
if (nu .lt. maxext) then
u(nu+1)=YINTEG
werr(nu+1)=dYINTEG
ni=1
endif
ENDIF
if (ikode .eq. 2) return
if (YINTEG .ne. 0 .and. dYINTEG .ne. 0)
1 WRITE (ISYSWR,1007) YINTEG,DYINTEG
990 WRITE (ISYSWR,1005) UP
RETURN
1000 FORMAT(/' ... fcn',I2,' Chi**2 =',F12.5,2X,'nf =',i6,i8
1 ,' Calls , edm =',E9.2,' ...')
1001 FORMAT (6X,' Parameter Value ',A,' ',
1 ' Intern.Value Int.Step Size')
1004 FORMAT (1X,6X,'Limits:', 2F14.5,A)
1005 FORMAT (' ...... Errors correspond to function change of ',
1 F10.4,' ......')
1007 FORMAT (5X,'Integr.Int. Exp:',2F14.3)
1008 FORMAT (5X,'Max.Intens. ',I2,' :',2F14.3)
END