var/da/da_radiance/da_write_oa_rad_ascii.inc

   1 subroutine da_write_oa_rad_ascii (it, ob, iv, re )
   2
   3    !---------------------------------------------------------------------------
   4    ! Purpose: write out OMB and OMA vector structure for radiance data.
   5    !---------------------------------------------------------------------------
   6
   7    implicit none
   8
   9    integer      ,     intent(in)  :: it       ! outer loop count
  10    type (y_type),     intent(in)  :: ob       ! Observation structure.
  11    type (iv_type),    intent(in)  :: iv       ! O-B structure.
  12    type (y_type),     intent(in)  :: re       ! O-A structure.
  13
  14    integer                        :: n        ! Loop counter.
  15    integer                        :: i, k     ! Index dimension.
  16    integer                        :: nlevelss ! Number of obs levels.
  17
  18    integer            :: ios, oma_rad_unit
  19    character(len=filename_len)  :: filename
  20    character(len=7)   :: surftype
  21    integer            :: ndomain
  22    logical            :: amsr2
  23
  24    if (trace_use) call da_trace_entry("da_write_oa_rad_ascii")
  25
  26    write(unit=message(1),fmt='(A)') 'Writing radiance OMA ascii file'
  27    call da_message(message(1:1))
  28
  29    do i = 1, iv%num_inst
  30       if (iv%instid(i)%num_rad < 1) cycle
  31
  32       ! count number of obs within the proc_domain
  33       !---------------------------------------------
  34       ndomain = 0
  35       do n =1,iv%instid(i)%num_rad
  36          if (iv%instid(i)%info%proc_domain(1,n)) then
  37             ndomain = ndomain + 1
  38          end if
  39       end do
  40       if (ndomain < 1) cycle
  41
  42       amsr2 = index(iv%instid(i)%rttovid_string,'amsr2') > 0
  43
  44       write(unit=filename, fmt='(i2.2,a,i4.4)') it,'_oma_'//trim(iv%instid(i)%rttovid_string)//'.', myproc
  45
  46       call da_get_unit(oma_rad_unit)
  47       open(unit=oma_rad_unit,file=trim(filename),form='formatted',iostat=ios)
  48       if (ios /= 0) then
  49          call da_error(__FILE__,__LINE__, &
  50             (/"Cannot open oma radiance file"//filename/))
  51       end if
  52       write(unit=oma_rad_unit,fmt='(a,a,i7,a,i5,a)') trim(iv%instid(i)%rttovid_string), &
  53                            ' number-of-pixels : ', ndomain, &
  54                            ' channel-number-of-each-pixel : ', iv%instid(i)%nchan, &
  55                            ' index-of-channels : '
  56       write(unit=oma_rad_unit,fmt='(10i5)') iv%instid(i)%ichan
  57       if ( amsr2 ) then
  58          write(unit=oma_rad_unit,fmt='(a)') ' pixel-info : i date scanpos landsea_mask  elv lat lon  satzen satazi solzen solazi clw'
  59       else
  60          write(unit=oma_rad_unit,fmt='(a)') ' pixel-info : i date scanpos landsea_mask  elv lat lon  satzen satazi solzen solazi'
  61       end if
  62       write(unit=oma_rad_unit,fmt='(a)') ' xb-surf-info : i t2m mr2m(ppmv) u10 v10 ps ts smois tslb snowh isflg &
  63                        & soiltyp vegtyp vegfra elev clwp'
  64       ndomain = 0
  65       do n=1,iv%instid(i)%num_rad
  66          if (iv%instid(i)%info%proc_domain(1,n)) then
  67             ndomain=ndomain+1
  68             if ( amsr2 ) then !write out clw
  69                write(unit=oma_rad_unit,fmt='(a,i7,2x,a,i6,i3,f6.0,6f8.2,f8.3)') 'INFO : ', ndomain, &
  70                                    iv%instid(i)%info%date_char(n), &
  71                                    iv%instid(i)%scanpos(n),        &
  72                                    iv%instid(i)%landsea_mask(n),   &
  73                                    iv%instid(i)%info%elv(n),       &
  74                                    iv%instid(i)%info%lat(1,n),     &
  75                                    iv%instid(i)%info%lon(1,n),     &
  76                                    iv%instid(i)%satzen(n),         &
  77                                    iv%instid(i)%satazi(n),         &
  78                                    iv%instid(i)%solzen(n),         &
  79                                    iv%instid(i)%solazi(n),         &
  80                                    iv%instid(i)%clw(n)
  81             else !no clw info
  82                write(unit=oma_rad_unit,fmt='(a,i7,2x,a,i6,i3,f6.0,6f8.2)') 'INFO : ', ndomain, &
  83                                    iv%instid(i)%info%date_char(n), &
  84                                    iv%instid(i)%scanpos(n),        &
  85                                    iv%instid(i)%landsea_mask(n),   &
  86                                    iv%instid(i)%info%elv(n),       &
  87                                    iv%instid(i)%info%lat(1,n),     &
  88                                    iv%instid(i)%info%lon(1,n),     &
  89                                    iv%instid(i)%satzen(n),         &
  90                                    iv%instid(i)%satazi(n),         &
  91                                    iv%instid(i)%solzen(n),         &
  92                                    iv%instid(i)%solazi(n)
  93             end if
  94             select case (iv%instid(i)%isflg(n))
  95             case (0) ;
  96                surftype = ' SEA : '
  97             case (1) ;
  98                surftype = ' ICE : '
  99             case (2) ;
 100                surftype = 'LAND : '
 101             case (3) ;
 102                surftype = 'SNOW : '
 103             case (4) ;
 104                surftype = 'MSEA : '
 105             case (5) ;
 106                surftype = 'MICE : '
 107             case (6) ;
 108                surftype = 'MLND : '
 109             case (7) ;
 110                surftype = 'MSNO : '
 111             end select
 112             write(unit=oma_rad_unit,fmt='(a,i7,9f10.2,3i3,f8.3,f10.2,f8.3,f15.5)') surftype, n, &
 113                              iv%instid(i)%t2m(n), &
 114                              iv%instid(i)%mr2m(n),   &
 115                              iv%instid(i)%u10(n), &
 116                              iv%instid(i)%v10(n),  &
 117                              iv%instid(i)%ps(n),  &
 118                              iv%instid(i)%ts(n),  &
 119                              iv%instid(i)%smois(n),  &
 120                              iv%instid(i)%tslb(n),  &
 121                              iv%instid(i)%snowh(n), &
 122                              iv%instid(i)%isflg(n), &
 123                              nint(iv%instid(i)%soiltyp(n)), &
 124                              nint(iv%instid(i)%vegtyp(n)), &
 125                              iv%instid(i)%vegfra(n), &
 126                              iv%instid(i)%elevation(n), &
 127                              iv%instid(i)%clwp(n), &
 128                              iv%instid(i)%cloud_frac(n)
 129
 130             write(unit=oma_rad_unit,fmt='(a)') 'OBS  : '
 131             write(unit=oma_rad_unit,fmt='(10f11.2)') ob%instid(i)%tb(:,n)
 132             write(unit=oma_rad_unit,fmt='(a)') 'BAK  : '
 133             write(unit=oma_rad_unit,fmt='(10f11.2)') iv%instid(i)%tb_xb(:,n)
 134             write(unit=oma_rad_unit,fmt='(a)') 'IVBC : '
 135             write(unit=oma_rad_unit,fmt='(10f11.2)')  iv%instid(i)%tb_inv(:,n)
 136             write(unit=oma_rad_unit,fmt='(a)') 'OMA  : '
 137             write(unit=oma_rad_unit,fmt='(10f11.2)')  re%instid(i)%tb(:,n)
 138             write(unit=oma_rad_unit,fmt='(a)') 'EMS  : '
 139             write(unit=oma_rad_unit,fmt='(10f11.2)')  iv%instid(i)%emiss(1:iv%instid(i)%nchan,n)
 140             write(unit=oma_rad_unit,fmt='(a)') 'ERR  : '
 141             write(unit=oma_rad_unit,fmt='(10f11.2)') iv%instid(i)%tb_error(:,n)
 142             write(unit=oma_rad_unit,fmt='(a)') 'QC   : '
 143             write(unit=oma_rad_unit,fmt='(10i11)') iv%instid(i)%tb_qc(:,n)
 144
 145             if (write_profile) then
 146                nlevelss  = iv%instid(i)%nlevels
 147                if ( rtm_option == rtm_option_rttov ) then
 148 #ifdef RTTOV
 149                   ! first, write RTTOV levels
 150                   write(unit=oma_rad_unit,fmt='(a)') 'RTM_level pres(mb) T(k) Q(ppmv)'
 151                   do k = 1, nlevelss
 152                      write(unit=oma_rad_unit,fmt='(i3,f10.2,f8.2,e11.4)') &
 153                         k, &                             ! RTTOV levels
 154                         coefs(i) % coef % ref_prfl_p(k) , &
 155                         iv%instid(i)%t(k,n) , &
 156                         iv%instid(i)%mr(k,n)
 157                   end do  ! end loop RTTOV levels
 158                   ! second, write WRF model levels
 159                   write(unit=oma_rad_unit,fmt='(a)') &
 160                      'WRF_level pres(mb) T(k) q(g/kg) clw(g/kg) rain(g/kg)'
 161                   do k=kts,kte
 162                      write(unit=oma_rad_unit,fmt='(i3,f10.2,f8.2,3e11.4)') &
 163                         k,  &                     ! WRF model levels
 164                         iv%instid(i)%pm(k,n) , &
 165                         iv%instid(i)%tm(k,n) , &
 166                         iv%instid(i)%qm(k,n)*1000 , &
 167                         iv%instid(i)%qcw(k,n)*1000.0, &
 168                         iv%instid(i)%qrn(k,n)*1000.0
 169                   end do   ! end loop WRF model levels
 170 #endif
 171                end if  ! end if rtm_option_rttov
 172
 173                if ( rtm_option == rtm_option_crtm ) then
 174 #ifdef CRTM
 175                   write(unit=oma_rad_unit,fmt='(a)') &
 176                      'level fullp(mb) halfp(mb) t(k) q(g/kg) water(mm) ice(mm) rain(mm) snow(mm) graupel(mm) hail(mm)'
 177                   if (crtm_cloud) then
 178                      do k=1,iv%instid(i)%nlevels-1
 179                         write(unit=oma_rad_unit,fmt='(i3,2f10.2,f8.2,13f8.3)') &
 180                            k,  &
 181                            iv%instid(i)%pf(k,n), &
 182                            iv%instid(i)%pm(k,n), &
 183                            iv%instid(i)%tm(k,n), &
 184                            iv%instid(i)%qm(k,n), &
 185                            iv%instid(i)%qcw(k,n), &
 186                            iv%instid(i)%qci(k,n), &
 187                            iv%instid(i)%qrn(k,n), &
 188                            iv%instid(i)%qsn(k,n), &
 189                            iv%instid(i)%qgr(k,n), &
 190                            iv%instid(i)%qhl(k,n), &
 191                            iv%instid(i)%rcw(k,n), &
 192                            iv%instid(i)%rci(k,n), &
 193                            iv%instid(i)%rrn(k,n), &
 194                            iv%instid(i)%rsn(k,n), &
 195                            iv%instid(i)%rgr(k,n), &
 196                            iv%instid(i)%rhl(k,n)
 197                      end do ! end loop profile
 198                   else  ! no cloud
 199                      do k=1,iv%instid(i)%nlevels-1
 200                         write(unit=oma_rad_unit,fmt='(i3,2f10.2,f8.2,7f8.3)') &
 201                            k,  &
 202                            iv%instid(i)%pf(k,n), &
 203                            iv%instid(i)%pm(k,n), &
 204                            iv%instid(i)%tm(k,n), &
 205                            iv%instid(i)%qm(k,n), &
 206                            0.0, &
 207                            0.0, &
 208                            0.0, &
 209                            0.0, &
 210                            0.0, &
 211                            0.0
 212                      end do ! end loop profile
 213                   end if ! end if crtm_cloud
 214 #endif
 215                end if  ! end if crtm_option
 216
 217             end if   ! end if write_profile
 218          end if    ! end if proc_domain
 219       end do     ! end do pixels
 220       close(unit=oma_rad_unit)
 221       call da_free_unit(oma_rad_unit)
 222    end do    !! end do instruments
 223
 224    if (trace_use) call da_trace_exit("da_write_oa_rad_ascii")
 225
 226 end subroutine da_write_oa_rad_ascii
 227
 228