! $Id: muri_abi_ocean_lut_mod.f90 2424 2017-12-18 14:30:36Z awalther $ ! ! 2019 Feb 6 : fixed scalong of aot_550nm with factor 10 module muri_abi_ocean_lut_mod use cx_sds_type_definitions_mod, only: use cx_sds_io_mod !, only: & ! cx_sds_finfo & ! , cx_sds_read !, cx_sds_read_6d_real use lib_array, only: interp1d use aw_lib_array, only: interp4d implicit none integer, parameter :: N_BANDS = 6 integer :: N_OPT integer, parameter :: N_FINE = 4 integer, parameter :: N_COARSE = 5 integer :: n_ws integer :: n_mode integer :: n_sol integer :: n_sat integer :: n_azi type muri_ocean_lut_type logical :: is_read = .false. integer :: n_opt real, allocatable :: sol(:) real, allocatable :: sat(:) real, allocatable :: azi(:) real, allocatable :: ws(:) real, allocatable :: aot_550nm (:) real, allocatable :: app_refl(:,:,:,:,:,:,:) !real, allocatable :: aot_aer(:,:,:,:,:,:,:) real, allocatable :: aot_aer_fine(:,:,:) real, allocatable :: aot_aer_coarse(:,:,:) real :: refl_fine(6,8,4) real :: refl_coarse(6,8,5) contains procedure ::read_lut => muri_lut_type__read_lut procedure ::sub_table => muri_lut_type__sub_table end type muri_ocean_lut_type type ( muri_ocean_lut_type),save :: ocean_lut contains ! ! ! subroutine muri_lut_type__read_lut(this, sol, sat, azi, ws, path) use lib_array, only:interp1d class(muri_ocean_lut_type ) :: this real, intent(in) :: sol real, intent(in) :: sat real, intent(in) :: azi real, intent(in) :: ws character(len = *) , intent(in), optional :: path character(len = 400) :: path_local character (len = 400)::lut_file character (len = 400)::lut_file_1 integer :: istatus integer :: ftype integer :: nsds integer :: natt character ( len = MAXNCNAM), allocatable :: sds_name(:) character ( len = MAXNCNAM), allocatable :: att_name(:) real,allocatable :: sol_zen_ang(:,:) real,allocatable :: sat_zen_ang(:,:) logical :: file_exists real, allocatable :: temp_2d_real(:,:) real, allocatable :: temp_6d_real(:,:,:,:,:,:) integer :: band character :: band_string integer :: shp_6d(6) integer :: i_ws,i_ws_a(1) integer :: i_opt,i_mode if ( present(path)) then path_local = trim(path) else path_local = trim('/apollo/cloud/scratch/mino/MURI_aerosol_LUT/') end if if ( this % is_read) return lut_file = trim(path_local)//trim('/ABI_Ocean_Aerosol_LUT_RB2_v212.hdf') !print*,'lut_file ', lut_file ! the following geo parameters are the same for all bands. ! (ABI band 1 is not using in over ocean retrieval) istatus = cx_sds_read ( trim(lut_file),'Solar_Zenith_Angles', temp_2d_real) allocate ( this %sol(size(temp_2d_real(:,1)) ), source = temp_2d_real(:,1)) istatus = cx_sds_read ( trim(lut_file),'View_Zenith_Angles',temp_2d_real ) allocate ( this %sat(size(temp_2d_real(:,1))), source = temp_2d_real(:,1)) istatus = cx_sds_read ( trim(lut_file),'Relative_Azimuth_Angles', temp_2d_real) allocate ( this %azi(size(temp_2d_real(:,1))), source = temp_2d_real(:,1)) !print*,'azi',this%azi istatus = cx_sds_read ( trim(lut_file),'Wind_Speed', temp_2d_real) allocate ( this %ws(size(temp_2d_real(:,1))), source = temp_2d_real(:,1)) !print*,'ws',this%ws istatus = cx_sds_read ( trim(lut_file),'AOT_at_550nm',temp_2d_real) !print*,temp_2d_real ! - add scale factor Jan 2019 AW allocate ( this %aot_550nm(size(temp_2d_real(1,:))), source = temp_2d_real(1,:)) ! new orientation this % aot_550nm(:) = temp_2d_real (1,:) /100. ! -- new scale factor with v212 print*,this % aot_550nm do band = 1,N_BANDS write ( band_string, "(i1)") band lut_file = trim(path_local)//trim('/ABI_Ocean_Aerosol_LUT_RB'//band_string//'_v212.hdf') !print*,lut_file INQUIRE(file = lut_file,EXIST=file_exists) if ( .not. file_exists) then print*,'MURI LUT file not there stopping' print*,'CLAVR-x was searching at ',lut_file stop end if istatus = cx_sds_read ( trim(lut_file), 'Apparent_Reflectance_ocean' , temp_6d_real) if ( band .eq. 1) then if ( .not. allocated( this%app_refl)) then shp_6d = shape(temp_6d_real) !print*,'shp_6d',shp_6d n_sol= shp_6d(1) n_sat = shp_6d(2) n_azi = shp_6d(3) n_opt = shp_6d(4) this % n_opt = n_opt n_ws = shp_6d(5) n_mode = shp_6d(6) ! !allocate ( this%app_refl(n_sol,n_sat,n_azi,N_bands,N_opt,N_mode)) allocate ( this%app_refl(n_sol,n_sat,n_azi,n_ws,N_bands,N_opt,N_mode)) end if this%app_refl = -999. end if !- closest wind speed ! - 16 Juy 2019 AW !i_ws_a = minloc ( abs( ws - this % ws)) !this%app_refl(:,:,:,band,:,:) = 0.0001 * temp_6d_real(:,:,:,:,i_ws_a(1),:) ! in this new version, we will interpolate do i_opt=1,n_opt do i_mode=1,n_mode do i_ws=1,n_ws this%app_refl(:,:,:,i_ws,band,i_opt,i_mode) = 0.0001 * temp_6d_real(:,:,:,i_opt,i_ws,i_mode) end do end do end do end do !print*,'shape this%app_refl',shape(this%app_refl) this % is_read = .true. end subroutine muri_lut_type__read_lut ! ! ! subroutine muri_lut_type__sub_table ( this, sol,sat,azi,ws ) use aw_lib_array use lib_array, only:interp1d class ( muri_ocean_lut_type ) :: this real, intent(in) :: sol real, intent(in) :: sat real, intent(in) :: azi real, intent(in) :: ws integer,parameter :: idp = selected_int_kind(13) integer,parameter :: sp = selected_real_kind(p=6,r=37) integer,parameter :: dp = selected_real_kind(p=15,r=307) real ,allocatable :: temp_4d(:,:,:,:) integer :: i_band, i_mode, i_opt integer :: i_azi,i_sat integer :: pos_sol integer :: pos_sat integer :: pos_azi real :: tmp_3d(10,15,21),tmp_2d(15,21),tmp_1d(21) real :: tmp_0d !print*,'azi',this%azi ! if(size(this%sol).gt.1) then do i_band=1,6 do i_opt=1,8 do i_mode=1,4 this % refl_fine(i_band,i_opt,i_mode)=interp4d(this%sol,this%sat & ,this%azi, this%ws & ,this%app_refl(:,:,:,:,i_band,i_opt,i_mode) & ,sol,sat,azi,ws & , bounds_error = .false., FILL_VALUE = -999.) end do end do end do do i_band=1,6 do i_opt=1,8 do i_mode=1,5 this % refl_coarse(i_band,i_opt,i_mode)=interp4d(this%sol,this%sat & ,this%azi, this%ws & ,this%app_refl(:,:,:,:,i_band,i_opt,i_mode+4) & ,sol,sat,azi,ws & , bounds_error = .false., FILL_VALUE = -999.) end do end do end do ! end if ! call dcomp_interpolation_weight(size( this % sol) , sol , this % sol & ! &, near_index = pos_sol ) ! call dcomp_interpolation_weight(size( this % sat) , sat , this % sat & ! &, near_index = pos_sat ) ! call dcomp_interpolation_weight(size( this % azi) , azi , this % azi & ! &, near_index = pos_azi ) ! this % refl_fine = this % app_refl(pos_sol,pos_sat,pos_azi,:,:,1:4) ! this % refl_coarse = this % app_refl(pos_sol,pos_sat,pos_azi,:,:,5:9) end subroutine muri_lut_type__sub_table end module muri_abi_ocean_lut_mod