subroutine process_HIRS_orbit(debug, month, sat, nscans, HIRS_2d, HIRS_1d, * cfsr_dir, year, s_jday, dom, along_track_idx, across_track_idx, * cloud_probability, land_cover, cloud_phase, csrbs, * shifted_FM_opt, output_prd_file, avhrr_ctp, outdata) c----------------------------------------------------------------------------------- c Description: Loop through HIRS L1b data and calculate cloud top pressures using c the CO2-slicing algorithm that follows the MODIS algorithm as closely c as possible. Use collocated PATMOS-X (AVHRR) data for cloudy vs. c clear sky discrimination. c c Input parameters: c debug Debug write flag c month HIRS data month c sat satellite c nscans Number scan lines in HIRS orbit c HIRS_2d 2D HIRS input data c HIRS_1d 1D HIRS input data c cfsr_dir directory that contains CFSR reanalysis profile data c year HIRS data year c s_jday HIRS starting day of year c dom HIRS starting day of month c along_track_idx AVHRR scan line numbers of collocated pixels c across_track_idx AVHRR element numbers of collocated pixels c cloud_probability collocated PATMOS-X (AVHRR GAC) cloud probabilities c land_cover collocated PATMOS-X (AVHRR GAC) land cover index c cloud_phase collocated PATMOS-X (AVHRR GAC) cloud phase index c csrbs monthly zonal mean clear sky radiance biases (CSRBs) c shifted_FM_opt shifted/unshifted forward model option c 1=unshifted (original), 2=shifted FM and Planck coefficients c output_prd_file output file name c avhrr_ctp PATMOS (AVHRR) cloud top pressure c c Output parameters: c outdata output data array c c Revision history: 12/10 R. Frey Original version c 01/11 R. Frey Changed dimensions of 'along_track_idx' c and 'across_track_idx' c 01/11 R. Frey Added AVHRR cloud fraction to calling c arguments of get_co2cld_nceprean1 c ('cloud_frc') c 01/11 R. Frey Added 'beg_elem', 'end_elem' for c +/- 32 deg. vza processing c 01/11 R. Frey Added 'rclr' from output of c get_co2cld_nceprean1 c 03/11 R. Frey Added 'land_cover' processing c subroutine get_land_fraction c 06/11 R. Frey Added array 'csrbs', subroutine c get_scene_type and associated c variables c 07/11 R. Frey Added 'shifted_FM_opt' c 10/13 R. Frey Renamed 'os_top_flag' to 'utls_flag' c 12/13 R. Frey Added 'cloud_phase' processing c 01/14 R. Frey Modified to use CFSR reanalysis files c 03/14 R. Frey Added 'avhrr_ctp' and call to get_avhrr_ctp c c Calls: integer function get_ancillary c subroutine get_satnode c subroutine get_cloud_fraction c subroutine get_land_fraction c subroutine get_co2cld_nceprean1 c subroutine get_scene_type c subroutine get_cloud_waterphase_fraction c subroutine get_avhrr_ctp c c----------------------------------------------------------------------------------- IMPLICIT NONE save c----------------------------------------------------------------------------------- c Parameters. integer ntbct parameter (ntbct = 8) integer num_2d_SDS parameter (num_2d_SDS = 13) integer num_1d_SDS parameter (num_1d_SDS = 3) real dtr parameter (dtr = 3.14159 / 180.0) c The following HIRS spots for +/- 32 deg. vza processing. integer beg_elem parameter (beg_elem = 1) integer end_elem parameter (end_elem = 56) c----------------------------------------------------------------------------------- c Calling argument variables. character*120 cfsr_dir, output_prd_file character*6 sat integer debug, nscans, year, s_jday, dom, along_track_idx(100,56,1100), * across_track_idx(100,56,1100), month, land_cover(409,15000), * shifted_FM_opt, cloud_phase(409,15000) real HIRS_2d(56,1100,num_2d_SDS), HIRS_1d(1100,num_1d_SDS), * cloud_probability(409,15000), outdata(1100, 56, 16), csrbs(180, 8, 3), * avhrr_ctp(409,15000) c Test c real tempdata(113,56,6) c----------------------------------------------------------------------------------- c Local variables. integer method, co2_flag, hc_flag, ci_flag, line, elem, s_time, jday, * nearnad_flag, landsea, dn_flag, out_flag, sat_node, * irt, andn, scene_type, j, nan_flag, bad_data_flag, utls_flag real nlat1, nlat2, node(1100), lat, lon, pres(26), temp(26), mixr(26), land, * sfctmp, prmsl, theta, tcold(ntbct), prsfc, sza_rad, cos_sza, tct, * pct, eca, ptrp, pwin, sza, cloud_frc, hirs_vis, rclr(ntbct), land_frc, * cth, wtrcld_frc, testbt, o3mr(6), tozone, mean_avhrr_ctp c----------------------------------------------------------------------------------- c External functions. integer get_ancillary c----------------------------------------------------------------------------------- c Loop through HIRS scan lines. c do line = 50, nscans do line = 1, nscans c----------------------------------------------------------------------------------- if(debug .gt. 0) then write(*,'(/,''Processing HIRS scan line '',i5)') line end if c----------------------------------------------------------------------------------- c Get satellite node (ascending=1 / descending=2). if( (line+1) .le. nscans) then nlat1 = HIRS_2d(28,line,1) nlat2 = HIRS_2d(28,line+1,1) if(nlat1 .eq. 0.0 .or. nlat2 .eq. 0.0) then node(line) = node(line-1) else call get_satnode(nlat1, nlat2, debug, sat_node) node(line) = sat_node*1.0 end if else node(line) = node(line-1) end if andn = int(node(line)) c Get scan day. jday = int(HIRS_1d(line,2)) c Get scan time. s_time = int(HIRS_1d(line,3)) c write(*,'(2f9.2,3i10)') nlat1, nlat2, line, jday, s_time c----------------------------------------------------------------------------------- c Loop through pixels on each scan line. do elem = beg_elem, end_elem c----------------------------------------------------------------------------------- c Get input HIRS data for current FOV. lat = HIRS_2d(elem,line,1) lon = HIRS_2d(elem,line,2) sza = HIRS_2d(elem,line,3) theta = HIRS_2d(elem,line,4) tcold(1) = HIRS_2d(elem,line,5) tcold(2) = HIRS_2d(elem,line,6) tcold(3) = HIRS_2d(elem,line,7) tcold(4) = HIRS_2d(elem,line,8) tcold(5) = HIRS_2d(elem,line,9) tcold(6) = HIRS_2d(elem,line,10) tcold(7) = HIRS_2d(elem,line,11) tcold(8) = HIRS_2d(elem,line,12) c Get day/night flag. Day=1, night=2. if(sza .le. 90.0) then dn_flag = 1 else dn_flag = 2 end if c Correct visible reflectance for solar zenith. sza_rad = dtr * sza cos_sza = cos(sza_rad) hirs_vis = HIRS_2d(elem,line,13) / cos_sza c----------------------------------------------------------------------------------- if(debug .gt. 0) then write(*,'(/,''Input data for line, element: '', 2i10)') line, elem write(*,'(7f15.3)') lat,lon,theta,sza,node(line),dn_flag*1.0,s_time*1.0 write(*,'(9f15.3,/)') tcold(1),tcold(2),tcold(3),tcold(4),tcold(5), & tcold(6),tcold(7),tcold(8),hirs_vis end if c----------------------------------------------------------------------------------- c Get ancillary data. c Use ancillary data from CFSR reanalysis files. irt = get_ancillary(debug, cfsr_dir, year, month, jday, dom, s_time, lat, & lon, pres, temp, mixr, land, sfctmp, prsfc, prmsl, & o3mr, tozone) landsea = int(land) c----------------------------------------------------------------------------------- c Get collocated AVHRR (from PATMOS-X) cloud data for the current HIRS FOV. c Collocated AVHRR GAC pixel cloud probabilities are combined to yield c a HIRS FOV cloud fraction. call get_cloud_fraction(debug, line, elem, along_track_idx, across_track_idx, * cloud_probability, cloud_frc) c----------------------------------------------------------------------------------- c Get AVHRR water cloud phase fraction (from PATMOS-X) for current HIRS FOV. call get_cloud_waterphase_fraction(debug, line, elem, along_track_idx, * across_track_idx, cloud_phase, wtrcld_frc) c----------------------------------------------------------------------------------- c Get collocated land cover fraction (portion of AVHRR GAC pixels in HIRS c IFOV that are classified as either land or coast by PATMOS-X). call get_land_fraction(debug, line, elem, along_track_idx, across_track_idx, * land_cover, land_frc) c----------------------------------------------------------------------------------- c Get collocated PATMOS-x (AVHRR) cloud top pressure statistics for the c current HIRS IFOV. call get_avhrr_ctp(debug, line, elem, along_track_idx, across_track_idx, * avhrr_ctp, mean_avhrr_ctp) c----------------------------------------------------------------------------------- c Get scene type. call get_scene_type(debug, landsea, andn, land_frc, scene_type) c----------------------------------------------------------------------------------- c Initialize CO2-slicing outputs. pct = -99.9 cth = -99.9 tct = -99.9 eca = -99.9 method = -99 co2_flag = 0 hc_flag = 0 ci_flag = 0 nearnad_flag = 0 utls_flag = 0 c----------------------------------------------------------------------------------- c Perform CO2-slicing algorithm if input data is valid. c Check here for calibration "stripe" or other bad input data. Assume all c bands are invalid if IR window data (band 8) is invalid. c Check for NaNs - probably temporary. nan_flag = 0 do j = 1, 6 testbt = tcold(j) if(testbt .ne. tcold(j)) nan_flag = 1 enddo if(HIRS_2d(elem,line,9) .gt. 0.0 .and. HIRS_2d(elem,line,9) .lt. 340.0 .and. * nan_flag .eq. 0) then bad_data_flag = 0 c Check collocated AVHRR cloud fraction corresponding to current FOV. if(cloud_frc .ge. 0.0) then call get_co2cld_cfsr( debug, sat, lat, lon, theta, year, month, & jday, pres, temp, mixr, o3mr, tozone, prsfc, prmsl, sfctmp, landsea, & tcold, cloud_frc, wtrcld_frc, csrbs, scene_type, shifted_FM_opt, & method, pct, eca, tct, ptrp, pwin, co2_flag, hc_flag, ci_flag, & utls_flag, nearnad_flag, rclr, cth ) end if else bad_data_flag = 1 end if c----------------------------------------------------------------------------------- c Results flag. if(cloud_frc .gt. 0.15 .and. pct .ne. -99.9) then c Valid cloud height from IRW, CO2-slicing, or MLC algorithms. out_flag = 3 else if(cloud_frc .ge. 0.0 .and. cloud_frc .le. 0.15 .and. pct .ne. -99.9) then c High cloud retrieval from CO2-slicing algorithm only. out_flag = 2 else if(cloud_frc .ge. 0.0 .and. cloud_frc .le. 0.15 .and. * pct .eq. -99.9 .and. bad_data_flag .eq. 0) then c No cloud height retrieval because HIRS FOV too clear. out_flag = 1 else c No valid cloud top pressure retrieval. out_flag = 0 end if c----------------------------------------------------------------------------------- c Fill output array. outdata(line, elem, 1) = lat outdata(line, elem, 2) = lon outdata(line, elem, 3) = pct outdata(line, elem, 4) = cth outdata(line, elem, 5) = tct outdata(line, elem, 6) = eca outdata(line, elem, 7) = method * 1.0 outdata(line, elem, 8) = dn_flag * 1.0 outdata(line, elem, 9) = node(line) outdata(line, elem, 10) = theta outdata(line, elem, 11) = s_time * 1.0 outdata(line, elem, 12) = cloud_frc outdata(line, elem, 13) = land_frc outdata(line, elem, 14) = out_flag * 1.0 outdata(line, elem, 15) = utls_flag * 1.0 outdata(line, elem, 16) = mean_avhrr_ctp c Test c outdata(line, elem, 16) = tcold(1) c outdata(line, elem, 17) = tcold(2) c outdata(line, elem, 18) = tcold(3) c outdata(line, elem, 19) = tcold(4) c outdata(line, elem, 20) = tcold(5) c----------------------------------------------------------------------------------- c End of CO2-slicing algorithm. c----------------------------------------------------------------------------------- enddo enddo c----------------------------------------------------------------------------------- return end