Page MenuHomeHEPForge
Files F11222285

interpolate_TS.f90
No OneTemporary
Actions

interpolate_TS.f90
View Options

!------------------------------------------------------------
program interpolate
! Compile with
! gfortran -C interpolate_TS.f90 -o interpolate_TS.exe
!
! This program turns an input table with different mass separations
! into a table with equal mass separation. It interpolates the columns
! linearly. NOTE: delta_mh must be a divisor of all entries in the first column.
!------------------------------------------------------------
implicit none
!------------------------------------------------------------
integer :: i,j,k,m,x
integer :: skip,ncol
integer :: l_in,l_out,f_in,f_out, step
character(LEN=100) :: filename
double precision :: delta_mh, sep
double precision, allocatable :: mh_XS_in(:,:),mh_XS_in_even(:,:),mh_XS_out(:,:)
character(LEN=150), allocatable :: description(:)
character(LEN=4) :: ending
character(LEN=2) :: a
!------------------------------------------------------------
filename='../HiggsBounds_KW/Expt_tables/CMStables/11014_CMS_H-WW_1.55fb-1' !main part of filename (without file extension)
! filename='sample/interpolate_TS_demo' !main part of filename (without file extension)
! filename='../../docs/Expt/2011134_ATLAS_WW-lnln_files/observed_fig14a'
ending='.txt' !file extension i.e. will read interpol_demo.txt
skip=5 !number of lines to be skipped before data starts
ncol=3 !number of columns of data
delta_mh=5. !Higgs mass separation for output
!------------------------------------------------------------
if((ncol.lt.2).or.(ncol.gt.100))stop 'wrong number of columns'
if(skip.lt.0)stop 'need to set skip'
write(*,*)'starting '//trim(filename)//'...'
! open input and output files
f_in=20
f_out=f_in*10
open(f_in,file=trim(filename)//ending)
open(f_out,file=trim(filename)//'_interpol'//ending)
x=getfilelength(f_in)
! find length of input file
l_in=x-skip
write(*,*)'skip=',skip
write(*,*)'l_in=',l_in
allocate(mh_XS_in(l_in,ncol))
allocate(mh_XS_in_even(l_in,ncol))
! read in input file
if(skip.ne.0)then
allocate(description(skip))
do i=1,skip
read(f_in,'(a)')description(i)
enddo
endif
do i=1,l_in
read(f_in,*)(mh_XS_in(i,j),j=1,ncol)
enddo
do i=1,l_in
do m=1,ncol
if(mh_XS_in(i,m).lt.0.0D0)stop 'negative input data'
if(m.ne.1) then
mh_XS_in_even(i,m)=mh_XS_in(i,m)
endif
enddo
if(mod(mh_XS_in(i,1),delta_mh).lt.delta_mh/2.) then
mh_XS_in_even(i,1)=int(mh_XS_in(i,1)/delta_mh)*delta_mh
else
mh_XS_in_even(i,1)=(int(mh_XS_in(i,1)/delta_mh)+1)*delta_mh
endif
enddo
l_out=int((mh_XS_in_even(l_in,1)-mh_XS_in_even(1,1))/delta_mh)+1
allocate(mh_XS_out(l_out,ncol))
k=0
do i=1,l_in-1
sep=mh_XS_in_even(i+1,1)-mh_XS_in_even(i,1)
step=int(sep/delta_mh)
do j=1,step
k=k+1
mh_XS_out(k,1) = mh_XS_in_even(i,1)+(j-1)*delta_mh
do m=2,ncol
mh_XS_out(k,m) = mh_XS_in_even(i,m)+(j-1)*(mh_XS_in_even(i+1,m)-mh_XS_in_even(i,m))/dble(step)
enddo
enddo
enddo
do m=1,ncol
mh_XS_out(l_out,m) = mh_XS_in_even(l_in,m)
enddo
write(a,'(I2)')ncol-1
do i=1,l_out
write(f_out,'(1f15.2,'//trim(a)//'g15.7)')(mh_XS_out(i,j),j=1,ncol)
enddo
!********************************************
deallocate(mh_XS_in)
deallocate(mh_XS_in_even)
deallocate(mh_XS_out)
close(f_in)
close(f_out)
contains
!-------------------------------------------------------------
function getfilelength(fileid)
!-------------------------------------------------------------
implicit none
!-------------------------------------output
integer fileid
!-----------------------------------internal
integer :: getfilelength
integer :: n,ios,m
character(len=100) :: temp
!-------------------------------------------
n= 0
do !this will count the number of lines which end in a newline character
read(fileid,*,iostat=ios)
if(ios.lt.0) exit
if(ios.gt.0)then
write(*,*)'error in input file',fileid ; call flush(6)
exit
endif
n = n + 1
enddo
rewind(fileid)
m = 0 ; temp="" !this will count the number of lines in the file, including the last one,
do ! even if it doesn't end in a newline character
read(fileid,'(a)',iostat=ios)temp
if(ios.lt.0) exit
m = m + 1
enddo
rewind(fileid)
if(m.ne.n)then
write(*,*)'error in input file: file id=',fileid ; call flush(6)
stop'Error: file needs to end with a newline character (see standard output for file id)'
endif
getfilelength=n
end function getfilelength
end program interpolate

File Metadata

Mime Type
text/plain
Expires
Wed, May 14, 11:42 AM (10 h, 42 m)
Storage Engine
blob
Storage Format
Raw Data
Storage Handle
5111491
Default Alt Text
interpolate_TS.f90 (4 KB)

Event Timeline