191 lines
5.4 KiB
Fortran
191 lines
5.4 KiB
Fortran
program deteff
|
|
|
|
* 1. linearer Fit (A+xB) einer Detektoreichmessung *
|
|
* 2. Normierung der Daten mit der verfeinerten Funktion *
|
|
* 3. Speicherung im Format fuer deteff.dat *
|
|
* 27.4.99 keller, Aenderung 15.7.99,28.9.99 *
|
|
* Aenderung 3.11.99: schneidet Punkte die 0 sind weg (HRPT!) *
|
|
* 7.3.00, 24.5.02, 9.8.02 *
|
|
* *
|
|
* > myfit deteff.f *
|
|
* > cp a.out deteff *
|
|
|
|
implicit none
|
|
|
|
external FIT_LIN_FUN
|
|
|
|
real par(16),err(16),int(2000),tth(2000),par1,ti,step,tf,lim
|
|
integer n,i,nret,l,count
|
|
character*36 spec
|
|
character*5 flag,flag2
|
|
character*78 title
|
|
|
|
call fit_init
|
|
|
|
! Function title and parameter names
|
|
call fit_userfun('STRAIGHT LINE', fit_lin_fun) ! function title, function
|
|
call fit_userpar('Bg(0)') ! first parameter: background at zero
|
|
call fit_userpar('dBg/dX') ! second parameter: slope
|
|
|
|
call fit_dat_options('cal=0') ! turns off calibration of raw data
|
|
|
|
do i=1,2000
|
|
int(i)=0.
|
|
enddo
|
|
|
|
|
|
call sys_getenv('dat_defspec',spec)
|
|
if (spec.eq.' ') spec='DMC'
|
|
call str_trim(spec,spec,l) ! Laenge bestimmen (l)
|
|
|
|
100 write(*,*)
|
|
write(*,'(x,a)')'DMC = 1 / HRPT = 2'
|
|
write(*,'(x,3a,$)')'Instrument (default: ',spec(1:l),'): '
|
|
read(*,'(a)')flag ! ^ schreibt spec von
|
|
! Zeichen 1 bis l
|
|
|
|
call str_upcase(flag,flag) ! schreibt Inhalt von flag gross
|
|
if (flag.eq.'1') flag='DMC'
|
|
if (flag.eq.'2') flag='HRPT'
|
|
if (flag.ne.' ') then
|
|
if (flag.ne.'DMC' .and. flag.ne.'HRPT') goto 100
|
|
spec=flag
|
|
endif
|
|
|
|
call sys_setenv('dat_defspec',spec)
|
|
|
|
|
|
C WRITE (*,'(X,A,$)') 'Calibration data files (e.g.: DMC/1999/
|
|
C A4-34) : '
|
|
C READ (*,'(A)') FILES
|
|
|
|
C CALL FIT_DAT(FILES)
|
|
|
|
C call fit_dat(' ')
|
|
C call fit_merge(0.02)
|
|
call fit_dat_merge(' ',0.02) ! ersetzt fit_dat und fit_merge
|
|
call fit_get_array('X',tth,2000,nret)
|
|
call fit_get_array('Y',int,2000,nret)
|
|
|
|
write(*,*)
|
|
write(*,'(x,a,$)')'Title of detector efficiency file: '
|
|
read(*,'(a)')title
|
|
|
|
if (title.eq.' ') then
|
|
if (nret.eq.400) title='DMC detector efficiency file'
|
|
if (nret.eq.1600) title='HRPT detector efficiency file'
|
|
endif
|
|
|
|
if (nret.ne.400) then
|
|
if (nret.ne.1600) then
|
|
write(*,*)
|
|
write(*,*)
|
|
write(*,*)
|
|
write(*,*)
|
|
write(*,*)' **************** WARNING ***************'
|
|
write(*,*)
|
|
write(*,*)' The number of data points does not match '
|
|
write(*,*)' the number of detectors of DMC or HRPT'
|
|
write(*,*)
|
|
write(*,*)
|
|
write(*,'(a,$)')' continue? (y/N): '
|
|
read(*,'(a)')flag2
|
|
call str_upcase(flag2,flag2)
|
|
if (flag2.ne.'Y') then
|
|
write(*,'(a)')' program aborted'
|
|
goto 102
|
|
endif
|
|
endif
|
|
endif
|
|
|
|
write(*,*)
|
|
write(*,'(x,a,$)')'slope=0 ? (Y/n): '
|
|
read(*,'(a)') flag2
|
|
call str_upcase(flag2,flag2)
|
|
|
|
TI=TTH(1)
|
|
TF=TTH(NRET)
|
|
STEP=(TTH(NRET)-TTH(1))/(NRET-1.)
|
|
|
|
PAR1=0
|
|
DO I=1,NRET
|
|
PAR1=PAR1+INT(I)
|
|
ENDDO
|
|
PAR(1)=PAR1/NRET ! starting value for Bg(0)
|
|
PAR(2)=0. ! starting value for slope
|
|
ERR(1)=1.
|
|
ERR(2)=1.
|
|
|
|
lim=par(1)/4. ! points below lim are excluded!
|
|
|
|
call fit_exclude(0.,0.,0.,lim)
|
|
|
|
CALL FIT_FUN(7,2,PAR,ERR)
|
|
if (flag2 .ne. 'N') call fit_fix(2)
|
|
CALL FIT_FIT(0)
|
|
CALL FIT_FIT(0)
|
|
CALL FIT_FIT(0)
|
|
CALL FIT_GET_ARRAY('P',PAR,2,N) ! reads fit parameters
|
|
CALL FIT_GET_ARRAY('E',ERR,2,N) ! reads parameter errors
|
|
|
|
OPEN (1,FILE='deteff.dat',STATUS='UNKNOWN')
|
|
|
|
WRITE (1,'(2A)')'#',TITLE
|
|
if (nret.gt.999) write(1,'(I4)')nret
|
|
if (nret.lt.1000) write(1,'(I3)')nret
|
|
|
|
count=0 ! count: # of excluded points
|
|
|
|
DO 101 I=1,NRET
|
|
|
|
if (int(i).lt.lim) then
|
|
write(1,'(F8.6)')0.0
|
|
count=count+1
|
|
else
|
|
WRITE(1,'(F8.6)')INT(I)/(PAR(1)+((I-1)*STEP+TI)*PAR(2))
|
|
endif
|
|
|
|
101 CONTINUE
|
|
|
|
CLOSE(1)
|
|
|
|
call str_trim(spec,spec,l)
|
|
WRITE(*,*)
|
|
WRITE(*,*)
|
|
WRITE(*,*)
|
|
WRITE(*,*)
|
|
WRITE(*,*)' New detector efficiency file: deteff.dat'
|
|
WRITE(*,'(3a)')' Copy this file into the directory /data/
|
|
Alnslib/data/',spec(1:l),'/<year>/'
|
|
WRITE(*,*)
|
|
WRITE(*,*)' Normalization function: P(1)+2th*P(2)'
|
|
WRITE(*,'(A,F12.5,A,F10.5,A,F10.5,A,F9.5,A)')' P(1)='
|
|
A,PAR(1),'(',ERR(1),'), P(2)=',PAR(2),'(',ERR(2),')'
|
|
WRITE(*,*)
|
|
if (count.gt.0) write (*,*)' Number of excluded points: ',count
|
|
WRITE(*,*)
|
|
|
|
102 END
|
|
|
|
|
|
real function fit_lin_fun(x,p,n,mode,cinfo)
|
|
! -------------------------------------------
|
|
|
|
implicit none
|
|
|
|
real x ! x-value
|
|
integer n ! number of parameters
|
|
real p(n) ! parameters
|
|
integer mode ! mode
|
|
integer cinfo ! calculation information (see below)
|
|
|
|
if (mode .eq. 0) then
|
|
|
|
! Define here your own function
|
|
|
|
fit_lin_fun=p(1)+x*p(2)
|
|
|
|
endif
|
|
end
|
|
|