Skip to content
Snippets Groups Projects
Commit a368d07b authored by Andreas Korpi-Lagg's avatar Andreas Korpi-Lagg
Browse files

sync

parent 78df26cb
No related branches found
No related tags found
No related merge requests found
...@@ -4,115 +4,115 @@ ...@@ -4,115 +4,115 @@
!call_pikaia. Only during a simulated measurement this value should be !call_pikaia. Only during a simulated measurement this value should be
!considered! !considered!
subroutine read_prefilter(file,prefilter_wlerr,wlin,nwl,prefilter,verbose) subroutine read_prefilter(file,prefilter_wlerr,wlin,nwl,prefilter,verbose)
use all_type use all_type
use interpol use interpol
use random use random
character(LEN=maxstr) :: file,line,word(2) character(LEN=maxstr) :: file,line,word(2)
real(8),dimension(maxwl) :: wlin real(8),dimension(maxwl) :: wlin
real(8),dimension(maxcwl) :: wl real(8),dimension(maxcwl) :: wl
real(4),dimension(maxcwl) :: val real(4),dimension(maxcwl) :: val
integer(4) :: err,cindex,sidx,j,idx(1),nw integer(4) :: err,cindex,sidx,j,idx(1),nw
integer(2) iw,verbose,nwl,ncwl,i integer(2) iw,verbose,nwl,ncwl,i
integer(4) :: iu integer(4) :: iu
real(4) :: prefilter_wlerr,wlerr real(4) :: prefilter_wlerr,wlerr
type(prefiltertyp) :: prefilter type(prefiltertyp) :: prefilter
prefilter%nwl=0 prefilter%nwl=0
prefilter%doprefilter=0 prefilter%doprefilter=0
if (trim(file(1:1)).eq.' ') return if (trim(file(1:1)).eq.' ') return
if (verbose.ge.1) write(*,*) 'Read prefilter file '//trim(file) if (verbose.ge.1) write(*,*) 'Read prefilter file '//trim(file)
open(iu,file=trim(file),iostat=err,action='READ',form='FORMATTED') open(iu,file=trim(file),iostat=err,action='READ',form='FORMATTED')
if (err.ne.0) then if (err.ne.0) then
write(*,*) 'Could not open prefilter file: ',trim(file) write(*,*) 'Could not open prefilter file: ',trim(file)
write(*,*) '****************************' write(*,*) '****************************'
write(*,*) 'ERROR: no prefilter applied!' write(*,*) 'ERROR: no prefilter applied!'
write(*,*) '****************************' write(*,*) '****************************'
prefilter%doprefilter=0 prefilter%doprefilter=0
return return
end if end if
iw=1 iw=1
do while (err.eq.0) do while (err.eq.0)
read(iu, fmt='(a)',iostat=err) line read(iu, fmt='(a)',iostat=err) line
cindex=index(trim(line),';')-1 cindex=index(trim(line),';')-1
if ((cindex.lt.0).and.(len(trim(line)).gt.0)) then if ((cindex.lt.0).and.(len(trim(line)).gt.0)) then
j=1 j=1
sidx=1 sidx=1
word='' word=''
!split input line into words separated by spaces !split input line into words separated by spaces
do while ((sidx.gt.0).and.(j.le.2).and.len(trim(line)).gt.0) do while ((sidx.gt.0).and.(j.le.2).and.len(trim(line)).gt.0)
sidx=index(line,' ') sidx=index(line,' ')
if (sidx.gt.1) then if (sidx.gt.1) then
word(j)=line(1:sidx-1) word(j)=line(1:sidx-1)
line=trim(line(sidx+1:)) line=trim(line(sidx+1:))
j=j+1 j=j+1
else else
if (sidx.eq.1) line=trim(line(sidx+1:)) if (sidx.eq.1) line=trim(line(sidx+1:))
end if end if
end do end do
nw=j-1 nw=j-1
sidx=index(trim(line),' ') sidx=index(trim(line),' ')
if (iw.gt.maxcwl) then if (iw.gt.maxcwl) then
write(*,*) 'Max. # of WL points in prefilter-file& write(*,*) 'Max. # of WL points in prefilter-file&
& exceeded: ', iw,maxcwl & exceeded: ', iw,maxcwl
write(*,*) 'Change MAXCWL in all_type.f90 and maxpar.pro' write(*,*) 'Change MAXCWL in all_type.f90 and maxpar.pro'
write(*,*) '****************************' write(*,*) '****************************'
write(*,*) 'ERROR: no prefilter applied!' write(*,*) 'ERROR: no prefilter applied!'
write(*,*) '****************************' write(*,*) '****************************'
prefilter%doprefilter=0 prefilter%doprefilter=0
return return
else else
if (nw.ge.2) then if (nw.ge.2) then
read(word(1),*) wl(iw) read(word(1),*) wl(iw)
read(word(2),*) val(iw) read(word(2),*) val(iw)
iw=iw+1 iw=iw+1
end if end if
end if
end if end if
end if end do
end do ncwl=iw-1
ncwl=iw-1 close(iu)
close(iu)
if (prefilter_wlerr.gt.0) then
if (prefilter_wlerr.gt.0) then CALL random_seed()
CALL random_seed() wlerr=gennor(0.,prefilter_wlerr)
wlerr=gennor(0.,prefilter_wlerr) if (verbose.ge.1) write(*,'(a,2f7.3)') &
if (verbose.ge.1) write(*,'(a,2f7.3)') & 'Adding WL error (\AA) to prefilter (normal noise): ',&
'Adding WL error (\AA) to prefilter (normal noise): ',& prefilter_wlerr,wlerr
prefilter_wlerr,wlerr wl=wl+wlerr
wl=wl+wlerr endif
endif
prefilter%doprefilter=1
prefilter%doprefilter=1
!the WL bins for the prefilter must
!the WL bins for the prefilter must ! be the
! be the !same as for the observations.
!same as for the observations.
! manual interpolation routine (to be
! manual interpolation routine (to be ! consistent with the fortran version)
! consistent with the fortran version) if ((minval(wlin(1:nwl)).lt.minval(wl(1:ncwl))).or.&
if ((minval(wlin(1:nwl)).lt.minval(wl(1:ncwl))).or.& (maxval(wlin(1:nwl)).ge.maxval(wl(1:ncwl)))) then
(maxval(wlin(1:nwl)).ge.maxval(wl(1:ncwl)))) then write(*,*) 'Prefilter curve does not cover the full WL range of the'//&
write(*,*) 'Prefilter curve does not cover the full WL range of the'//& ' observations. Please extend the WL coverage for the prefilter curve.'
' observations. Please extend the WL coverage for the prefilter curve.' stop
stop end if
end if prefilter%wl(1:nwl)=wlin(1:nwl)
prefilter%wl(1:nwl)=wlin(1:nwl) do i=1,nwl
do i=1,nwl idx=minloc(abs(wlin(i)-wl(1:ncwl)))
idx=minloc(abs(wlin(i)-wl(1:ncwl))) prefilter%val(i)=(val(idx(1)+1)-val(idx(1)))/(wl(idx(1)+1)-wl(idx(1)))* &
prefilter%val(i)=(val(idx(1)+1)-val(idx(1)))/(wl(idx(1)+1)-wl(idx(1)))* & (wlin(i)-wl(idx(1)))+val(idx(1))
(wlin(i)-wl(idx(1)))+val(idx(1)) end do
end do
!normalize prefilter function to it's max
!normalize prefilter function to it's max ! prefilter%val(1:nwl)= prefilter%val(1:nwl)/maxval(prefilter%val(1:nwl))
! prefilter%val(1:nwl)= prefilter%val(1:nwl)/maxval(prefilter%val(1:nwl)) prefilter%nwl=nwl
prefilter%nwl=nwl
end subroutine read_prefilter
end subroutine read_prefilter
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment