1!fast disturbance smoother for simulation 2subroutine smoothsimfast(yt, ymiss, timevar, zt, ht,tt, rtv,qt,a1, ft,kt,& 3finf, kinf, dt, jt, p, m, n,r,epshat,etahat,rt0,rt1,needeps) 4 5 implicit none 6 7 logical, intent(in) :: needeps 8 integer, intent(in) :: p, m, r,n,dt,jt 9 integer :: t, i 10 integer, intent(in), dimension(n,p) :: ymiss 11 integer, intent(in), dimension(5) :: timevar 12 double precision, intent(in), dimension(n,p) :: yt 13 double precision, intent(in), dimension(p,m,(n-1)*timevar(1)+1) :: zt 14 double precision, intent(in), dimension(p,p,(n-1)*timevar(2)+1) :: ht 15 double precision, intent(in), dimension(m,m,(n-1)*timevar(3)+1) :: tt 16 double precision, intent(in), dimension(m,r,(n-1)*timevar(4)+1) :: rtv 17 double precision, intent(in), dimension(r,r,(n-1)*timevar(5)+1) :: qt 18 double precision, intent(in), dimension(m) :: a1 19 double precision, intent(in), dimension(p,n) :: ft,finf 20 double precision, intent(in), dimension(m,p,n) :: kt,kinf 21 double precision, intent(inout), dimension(p,n) :: epshat 22 double precision, intent(inout), dimension(r,n) :: etahat 23 double precision, dimension(p,n) :: vt 24 double precision, dimension(m) :: at 25 double precision, dimension(m,m) :: im 26 double precision, intent(inout), dimension(m) :: rt0,rt1 27 double precision :: lik 28 lik = 0.0d0 29 at = a1 30 if(dt .GT. 0) then 31 !diffuse filtering begins 32 do t = 1, dt - 1 33 call dfilter1stepnv(ymiss(t,:),yt(t,:),transpose(zt(:,:,(t-1)*timevar(1)+1)),& 34 tt(:,:,(t-1)*timevar(3)+1),at,vt(:,t),ft(:,t),kt(:,:,t),& 35 finf(:,t),kinf(:,:,t),p,m,p,lik) 36 end do 37 38 t = dt 39 call dfilter1stepnv(ymiss(t,:),yt(t,:),transpose(zt(:,:,(t-1)*timevar(1)+1)),& 40 tt(:,:,(t-1)*timevar(3)+1),at,vt(:,t),ft(:,t),kt(:,:,t),& 41 finf(:,t),kinf(:,:,t),p,m,jt,lik) 42 !non-diffuse filtering begins 43 if(jt .LT. p) then 44 call filter1stepnv(ymiss(t,:),yt(t,:),transpose(zt(:,:,(t-1)*timevar(1)+1)),& 45 tt(:,:,(t-1)*timevar(3)+1),at,vt(:,t),ft(:,t),kt(:,:,t),p,m,jt,lik) 46 end if 47 end if 48 !Non-diffuse filtering continues from t=d+1, i=1 49 do t = dt + 1, n 50 call filter1stepnv(ymiss(t,:),yt(t,:),transpose(zt(:,:,(t-1)*timevar(1)+1)),& 51 tt(:,:,(t-1)*timevar(3)+1),at,vt(:,t),ft(:,t),kt(:,:,t),p,m,0,lik) 52 end do 53 54 55 !smoothing begins 56 57 im = 0.0d0 58 do i = 1, m 59 im(i,i) = 1.0d0 60 end do 61 62 rt0 = 0.0d0 63 64 do t = n, dt+1, -1 65 call smooth1step(ymiss(t,:), transpose(zt(:,:,(t-1)*timevar(1)+1)), ht(:,:,(t-1)*timevar(2)+1), & 66 tt(:,:,(t-1)*timevar(3)+1), rtv(:,:,(t-1)*timevar(4)+1), qt(:,:,(t-1)*timevar(5)+1), vt(:,t), & 67 ft(:,t),kt(:,:,t), im,p,m,r,1,rt0,etahat(:,t),epshat(:,t),needeps) 68 end do 69 70 if(dt .GT. 0) then 71 t = dt 72 if(jt .LT. p) then 73 call smooth1step(ymiss(t,:), transpose(zt(:,:,(t-1)*timevar(1)+1)), ht(:,:,(t-1)*timevar(2)+1), & 74 tt(:,:,(t-1)*timevar(3)+1), rtv(:,:,(t-1)*timevar(4)+1), qt(:,:,(t-1)*timevar(5)+1), vt(:,t), & 75 ft(:,t),kt(:,:,t), im,p,m,r,jt+1,rt0,etahat(:,t),epshat(:,t),needeps) 76 end if 77 rt1 = 0.0d0 78 call dsmooth1step(ymiss(t,:), transpose(zt(:,:,(t-1)*timevar(1)+1)), ht(:,:,(t-1)*timevar(2)+1), & 79 tt(:,:,(t-1)*timevar(3)+1), rtv(:,:,(t-1)*timevar(4)+1), qt(:,:,(t-1)*timevar(5)+1), vt(:,t), & 80 ft(:,t),kt(:,:,t), im,p,m,r,jt,rt0,rt1,finf(:,t),kinf(:,:,t),etahat(:,t),epshat(:,t),needeps) 81 do t = (dt - 1), 1, -1 82 call dsmooth1step(ymiss(t,:), transpose(zt(:,:,(t-1)*timevar(1)+1)), ht(:,:,(t-1)*timevar(2)+1), & 83 tt(:,:,(t-1)*timevar(3)+1), rtv(:,:,(t-1)*timevar(4)+1), qt(:,:,(t-1)*timevar(5)+1), vt(:,t), & 84 ft(:,t),kt(:,:,t), im,p,m,r,p,rt0,rt1,finf(:,t),kinf(:,:,t),etahat(:,t),epshat(:,t),needeps) 85 end do 86 end if 87 88end subroutine smoothsimfast 89