da_get_innov_vector_ssmi_rv.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_ssmi_rv (it, grid, ob, iv)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 integer, intent(in) :: it ! External iteration.
10 type(domain), intent(in) :: grid ! first guess state.
11 type(y_type), intent(in) :: ob ! Observation structure.
12 type(iv_type), intent(inout) :: iv ! O-B structure.
13
14 integer :: n ! Loop counter.
15 integer :: i, j ! Index dimension.
16 real :: dx, dxm ! Interpolation weights.
17 real :: dy, dym ! Interpolation weights.
18 real :: model_tpw ! Model value tpw at oblocation
19 real :: model_speed ! Model value speed at oblocation
20 integer :: itpw, itpwf, ispeed, ispeedf
21
22 if (trace_use) call da_trace_entry("da_get_innov_vector_ssmi_rv")
23
24 itpw = 0 ; itpwf = 0 ; ispeed = 0 ; ispeedf = 0
25 do n=iv%info(ssmi_rv)%n1,iv%info(ssmi_rv)%n2
26
27 ! compute innovation vector
28 ! =========================
29
30 ! Obs coordinates on model grid
31
32 ! TPW
33
34 i = iv%info(ssmi_rv)%i(1,n)
35 j = iv%info(ssmi_rv)%j(1,n)
36 dx = iv%info(ssmi_rv)%dx(1,n)
37 dy = iv%info(ssmi_rv)%dy(1,n)
38 dxm = iv%info(ssmi_rv)%dxm(1,n)
39 dym = iv%info(ssmi_rv)%dym(1,n)
40
41 iv % ssmi_rv(n) % tpw % inv = 0.0
42 if (abs(ob % ssmi_rv(n) % tpw - missing_r) > 1.0 .and. iv % ssmi_rv(n) % tpw % qc >= obs_qc_pointer) then
43 model_tpw = dym*(dxm*grid%xb%tpw(i,j) + dx*grid%xb%tpw(i+1,j)) + dy *(dxm*grid%xb%tpw(i,j+1) + dx*grid%xb%tpw(i+1,j+1))
44 iv % ssmi_rv(n) % tpw % inv = ob % ssmi_rv(n) % tpw - model_tpw
45 end if
46
47 ! surface wind speed
48
49 iv % ssmi_rv(n) % speed % inv = 0.0
50 if (abs(ob % ssmi_rv(n) % speed - missing_r) > 1.0 .and. &
51 iv % ssmi_rv(n) % speed % qc >= obs_qc_pointer) then
52
53 model_speed = dym*(dxm*grid%xb%speed(i,j ) + dx*grid%xb%speed(i+1,j )) &
54 + dy *(dxm*grid%xb%speed(i,j+1) + dx*grid%xb%speed(i+1,j+1))
55 iv % ssmi_rv(n) % speed % inv = ob % ssmi_rv(n) % speed - model_speed
56 end if
57 end do
58
59 !------------------------------------------------------------------
60 ! Perform optional maximum error check:
61 !------------------------------------------------------------------
62
63 if (check_max_iv) call da_check_max_iv_ssmi_rv(iv, it, itpw,itpwf, ispeed,ispeedf)
64
65 if (rootproc .and. check_max_iv_print) then
66 write(unit = check_max_iv_unit, fmt ='(A,i5,A)') &
67 'For outer iteration ',it,&
68 ' Total Rejections for SSMI(TPW and SPEED) follows:'
69 write(unit = check_max_iv_unit, fmt = '(/,4(2(A,I6),/))') &
70 'Number of failed ssmi tpw observations: ',itpwf, ' on ',itpw, &
71 'Number of failed ssmi speed observations: ',ispeedf,' on ',ispeed
72 end if
73
74 if (trace_use) call da_trace_exit("da_get_innov_vector_ssmi_rv")
75
76 end subroutine da_get_innov_vector_ssmi_rv
77
78