! $Id$ ! module sat_mod use sensor_mod use date_tools_mod implicit none integer , dimension(16) , private , parameter :: & modis_chn_list = [ 8 , 9 , 10 , 12 , 1 , 15 , 2 , 5 , 26 , 6 , 7 , 20 , 22 , 29 , 31 , 32 ] type sat_data_type real, allocatable, dimension(:,:) :: rad real, allocatable, dimension(:,:) :: ref real, allocatable, dimension(:,:) :: bt end type sat_data_type type geo_data_type real, allocatable, dimension(:,:) :: lat real, allocatable, dimension(:,:) :: lon real, allocatable, dimension(:,:) :: sol_zen real, allocatable, dimension(:,:) :: sat_zen real, allocatable, dimension(:,:) :: rel_az real, allocatable, dimension(:,:) :: glint_zen end type geo_data_type type sat_type character (len = 255) :: file character (len = 255) :: l1b_path integer :: ETsensor integer :: num_pix integer :: num_elem integer :: num_seg integer :: num_elem_per_seg_general integer :: num_elem_this_seg integer :: idx_seg logical , dimension(40) :: chan_on type ( sat_data_type) , dimension(40) :: data type ( geo_data_type) :: geo type ( date_type ) :: start_time type ( date_type ) :: end_time type ( sensor_type) :: sensor_const contains procedure :: initialize_sat procedure :: set_dimensions procedure :: read_l1b procedure :: deallocate_all end type sat_type contains subroutine initialize_sat ( this , conf , i_file) use config_mod class ( sat_type) :: this type ( conf_user_opt_type ) , intent ( in ) :: conf integer , intent (in) :: i_file this % file = conf % file % infile (i_file) this % ETsensor = conf % file % ETsensor (i_file) this % l1b_path = conf % file % l1b_path associate ( sensor => this % sensor_const) call sensor % populate(this%file, this%l1b_path) end associate call this % set_dimensions print*,' == Sensor is ===> ',this%sensor_const%sensor_name end subroutine initialize_sat ! ! ! subroutine set_dimensions ( this ) use config_mod class ( sat_type) :: this integer :: year, month , day , start_hour, start_minute integer :: start_sec, end_hour , end_minute, end_sec, orbit character ( len =100)::infile associate ( time_0 => this % start_time , & time_1 => this % end_time ) select case ( this % ETsensor ) case (ETsensor_aqua) case (ETsensor_viirs) ! --- Read data from the file name infile = this % file read(Infile(12:15), fmt="(I4)") year read(Infile(16:17), fmt="(I2)") month read(Infile(18:19), fmt="(I2)") day read(Infile(22:23), fmt="(I2)") start_hour read(Infile(24:25), fmt="(I2)") start_minute read(Infile(26:27), fmt="(I2)") start_sec read(Infile(31:32), fmt="(I2)") end_hour read(Infile(33:34), fmt="(I2)") end_minute read(Infile(35:36), fmt="(I2)") end_sec read(Infile(40:44), fmt="(I5)") orbit call time_0 % set_date ( & & year = year & & , month = month & & , day = day & & , hour = start_hour & & , minute = start_minute & & , second = start_sec ) call time_1 % set_date ( & & year = year & & , month = month & & , day = day & & , hour = end_hour & & , minute = end_minute & & , second = end_sec ) this % num_pix = 3200 this % num_elem = 768 this % num_seg = 4 this % num_elem_per_seg_general = 200 ! - later via option file ! call time_0 % print_data ! call time_1 % print_data end select end associate end subroutine set_dimensions ! =============================================== ! ! =============================================== subroutine read_l1b ( this , idx_seg) use viirs_read_mod , only : & viirs_data_config & , viirs_data_out & , get_viirs_data use planck_mod , only : & planck_temp use viewing_geometry_module, only : & glint_angle class ( sat_type) :: this integer , intent(in) :: idx_seg type ( viirs_data_config ) :: v_conf type ( viirs_data_out ) :: out integer :: i_mband ! - configure viirs interface this % num_elem_this_seg = this % num_elem_per_seg_general if ( idx_seg == this % num_seg) then this % num_elem_this_seg = this % num_elem - (idx_seg - 1 ) * this % num_elem_per_seg_general end if v_conf % chan_on_rfl_mband = this % chan_on (modis_chn_list) v_conf % file_gmtco_base = trim( this % file) v_conf % dir_1b = trim(this %l1b_path ) v_conf % offset (1) = 1 v_conf % offset (2) = (idx_seg -1 ) * this % num_elem_this_seg + 1 v_conf % count (1) = 3200 v_conf % count (2) = this % num_elem_this_seg call get_viirs_data ( v_conf, out ) allocate ( this % geo % lat (3200, this % num_elem_this_seg)) allocate ( this % geo % lon (3200, this % num_elem_this_seg)) allocate ( this % geo % sol_zen (3200, this % num_elem_this_seg)) allocate ( this % geo % sat_zen (3200, this % num_elem_this_seg)) allocate ( this % geo % rel_az (3200, this % num_elem_this_seg)) allocate ( this % geo % glint_zen (3200, this % num_elem_this_seg)) this % geo % lat = out % geo % lat this % geo % lon = out % geo % lon this % geo % sol_zen = out % geo % solzen this % geo % sat_zen = out % geo % satzen this % geo % rel_az = out % geo % relaz this % geo % glint_zen = glint_angle ( this % geo % sol_zen , this % geo % sat_zen , this % geo % rel_az ) do i_mband = 1, 16 if ( this % chan_on (modis_chn_list(i_mband)) .eqv. .false. ) cycle if ( i_mband < 12 ) then allocate ( this % data (modis_chn_list(i_mband)) % ref (3200,this % num_elem_this_seg)) this % data (modis_chn_list(i_mband)) % ref = out % mband (i_mband) % ref else allocate ( this % data (modis_chn_list(i_mband)) % rad (3200,this % num_elem_this_seg)) allocate ( this % data ( modis_chn_list(i_mband)) % bt (3200, this % num_elem_this_seg)) this % data (modis_chn_list(i_mband)) % rad = out % mband (i_mband) % rad call planck_temp (this % data ( modis_chn_list(i_mband)) % rad, 'VIIRS' & , modis_chn_list(i_mband) , this % data (modis_chn_list(i_mband) )% bt ) end if end do end subroutine ! ====================================== ! ! ====================================== subroutine deallocate_all ( this ) class ( sat_type) :: this integer :: i_mband if ( allocated ( this % geo % lon ) ) deallocate ( this % geo % lat ) if ( allocated ( this % geo % lon ) ) deallocate ( this % geo % lon ) if ( allocated ( this % geo % sol_zen ) ) deallocate ( this % geo % sol_zen ) if ( allocated ( this % geo % sat_zen ) ) deallocate ( this % geo % sat_zen ) if ( allocated ( this % geo % rel_az ) ) deallocate ( this % geo % rel_az ) if ( allocated ( this % geo % glint_zen ) ) deallocate ( this % geo % glint_zen ) do i_mband = 1 , 16 if ( allocated ( this % data (modis_chn_list(i_mband)) % ref ) ) deallocate ( this % data (modis_chn_list(i_mband)) % ref ) if ( allocated ( this % data (modis_chn_list(i_mband)) % rad ) ) deallocate ( this % data (modis_chn_list(i_mband)) % rad ) if ( allocated ( this % data (modis_chn_list(i_mband)) % bt ) ) deallocate ( this % data (modis_chn_list(i_mband)) % bt ) end do end subroutine deallocate_all end module sat_mod