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