da_get_innov_vector_ssmi_rv.inc
References to this file elsewhere.
1 subroutine da_get_innov_vector_ssmi_rv( it, xb, ob, iv)
2
3 !-----------------------------------------------------------------------
4 ! Purpose: TBD
5 !-----------------------------------------------------------------------
6
7 implicit none
8
9 integer, intent(in) :: it ! External iteration.
10 type(xb_type), intent(in) :: xb ! first guess state.
11 type(y_type), intent(in) :: ob ! Observation structure.
12 type(ob_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 (iv % num_ssmi_retrieval < 1) return
23
24 if (trace_use) call da_trace_entry("da_get_innov_vector_ssmi_rv")
25
26 itpw = 0 ; itpwf = 0 ; ispeed = 0 ; ispeedf = 0
27 do n=1, iv % num_ssmi_retrieval
28
29 ! COMPUTE inNOVATION VECTOR
30 ! =========================
31
32 ! Obs coordinates on model grid
33
34 ! TPW
35
36 i = iv%ssmi_retrieval(n)%loc%i
37 j = iv%ssmi_retrieval(n)%loc%j
38 dx = iv%ssmi_retrieval(n)%loc%dx
39 dy = iv%ssmi_retrieval(n)%loc%dy
40 dxm = iv%ssmi_retrieval(n)%loc%dxm
41 dym = iv%ssmi_retrieval(n)%loc%dym
42
43
44 iv % ssmi_retrieval(n) % tpw % inv = 0.0
45 if (abs(ob % ssmi_retrieval(n) % tpw - missing_r) > 1. .and. &
46 iv % ssmi_retrieval(n) % tpw % qc >= obs_qc_pointer) then
47 model_tpw = dym*(dxm*xb%tpw(i,j) + dx*xb%tpw(i+1,j)) + &
48 dy *(dxm*xb%tpw(i,j+1) + dx*xb%tpw(i+1,j+1))
49
50 iv % ssmi_retrieval(n) % tpw % inv = &
51 ob % ssmi_retrieval(n) % tpw - model_tpw
52 end if
53
54 ! SURFACE WinD SPEED
55
56 iv % ssmi_retrieval(n) % speed % inv = 0.0
57 if (abs(ob % ssmi_retrieval(n) % speed - missing_r) > 1. .and. &
58 iv % ssmi_retrieval(n) % speed % qc >= obs_qc_pointer) then
59
60 model_speed = dym*(dxm*xb%speed(i,j ) + dx*xb%speed(i+1,j )) &
61 + dy *(dxm*xb%speed(i,j+1) + dx*xb%speed(i+1,j+1))
62 iv % ssmi_retrieval(n) % speed % inv = &
63 ob % ssmi_retrieval(n) % speed - model_speed
64 end if
65
66 !------------------------------------------------------------------
67 ! Perform optional maximum error check:
68 !------------------------------------------------------------------
69
70 if (check_max_iv) then
71 call da_check_max_iv_ssmi_rv(it, iv % ssmi_retrieval(n),&
72 itpw,itpwf, ispeed,ispeedf)
73 end if
74 end do
75
76 if (rootproc .and. check_max_iv_print) then
77 write(unit = check_max_iv_unit, fmt ='(A,i5,A)') &
78 'For outer iteration ',it,&
79 ' Total Rejections for SSMI(TPW and SPEED) follows:'
80 write(unit = check_max_iv_unit, fmt = '(/,4(2(A,I6),/))') &
81 'Number of failed ssmi tpw observations: ',itpwf, ' on ',itpw, &
82 'Number of failed ssmi speed observations: ',ispeedf,' on ',ispeed
83 end if
84
85 if (trace_use) call da_trace_exit("da_get_innov_vector_ssmi_rv")
86
87 end subroutine da_get_innov_vector_ssmi_rv
88
89