201                      do j=1,ntpts
203                         read(567,*) (BCt(j,n),n=1,4)
251                do j=1,ntpts
310       integer  iBC(numnp),nlast,i,j,nper   local
321          do j=2,nlast   !loop to find the interval that we are in
323             if(BCtptr(i)%p(j,4).gt.t) then  
326       wr=(t-BCtptr(i)%p(j-1,4))/(BCtptr(i)%p(j,4)-BCtptr(i)%p(j-1,4))
327                BC(nbct(i),3:5)= BCtptr(i)%p(j-1,1:3)*(one-wr) 
328      &                        + BCtptr(i)%p(j,1:3)*wr
364             do j=1,numDataImp
365                read(817,*) (ValueImpt(j,n,k),n=1,2) ! n=1 time, 2 value
385       integer nlast, nper, k, j , jstep  local
393          do j=2,nlast   !loop to find the interval that we are in
395             if(ValueImpt(j,1,k).gt.ptime) then  ! this is upper bound, j-1 is lower
396                wr=(ptime-ValueImpt(j-1,1,k))
397      &             / ( ValueImpt(j,1,k)-ValueImpt(j-1,1,k) )
398                ValueListImp(jstep,k)= ValueImpt(j-1,2,k)*(one-wr) 
399      &                        + ValueImpt(j,2,k)*wr
416       integer nptf, cutfreq, j, k, m, s, Filtime(nptf)  local
422       do j=1,nptf
423          Filtime(j) = j-1
424          Window(j) = 0.42+0.5*cos(2*pi*Filtime(j)/nptf)
425      &              +0.08*cos(4*pi*Filtime(j)/nptf)
426          Sinc(j) = sin(pi*Filtime(j)*windK/nptf)/sin(pi*Filtime(j)/nptf) 
430       do j=1,nptf     
431          FilterSW(j) = Window(nptf+1-j)*Sinc(nptf+1-j) !filter for convolution
436          do j=1,nptf
437             s=modulo(m-nptf+j,nptf)
442      &              +FilterSW(j)*DataHist(s,:)/nptf !filter convolution
468             do j=1,3
469                read(818,*) ValueListRCR(j,k) ! reads Rp,C,Rd
471             do j=1,numDataRCR
472                read(818,*) (ValuePdist(j,n,k),n=1,2) ! n=1 time, 2 value
489       integer nlast, nper, k, j  local
497          do j=2,nlast   !loop to find the interval that we are in
499             if(ValuePdist(j,1,k).gt.ptime) then  ! this is upper bound, j-1 is lower
500                wr=(ptime-ValuePdist(j-1,1,k))
501      &             / ( ValuePdist(j,1,k)-ValuePdist(j-1,1,k) )
502                Pdist(k)= ValuePdist(j-1,2,k)*(one-wr) 
503      &                        + ValuePdist(j,2,k)*wr
527         do j=1,nTimePoint+1
528             pressHist(k) = pressHist(k) + QHist(j,k)*betas(j,k)