OB.DAAC Logo
NASA Logo
Ocean Color Science Software

ocssw V2022
find_v_viirs_nc4.f
Go to the documentation of this file.
1 ! @TODO: search() can't handle data that is at the last node of the table
2 ! TODO: Convert all of this to F95
3 
4 ! NOTES: Some "3.5" are sprinkled in here. The AOT tables only currently go up to AOT 3.5, but we
5 ! extrapolate out to a maximum of 5.0. However, we don't want to change the AE or SSA values, so
6 ! in some cases, the AOT values were limited to 3.5 or 3.5 is used in order to maintain the AE and SSA
7 ! calculations in here. At the end, the spectral AOT values are truncated to 3.5 if needed when returned
8 ! to modis.f90. Again, this is done in order to maintain the current AE and SSA values.
9  subroutine find_v_viirs(realbuf, tmpvg, outbuf, qa_flag, px_elev, lm, windsp, wv)
10 
11  use landcover, only: get_landcover
12  use calendars, only: gdatetime, season_from_doy, gregorian_from_doy
13  use modis_surface, only: get_brdfcorr_sr,
14  1 get_aot500,
15  1 terrain_flag,
16  1 terrain_flag_new,
17  1 get_geographic_zone,
18  1 get_background_aod, ! 9 January 2018 JLee
19  1 get_ler412,
20  1 get_ler470,
21  1 get_ler650,
22  1 get_modis_ler865,
23  1 get_swir_coeffs412,
24  1 get_swir_coeffs470,
25  1 get_swir_stderr412,
26  1 get_swir_stderr470,
27  1 get_swir_range
28 
29  use viirs_aerosol_luts, only: aero_412,
30  1 aero_470,
31  1 aero_650,
32  1 new_intep,
33  1 aero_412_abs,
34  1 aero_470_abs,
35  1 aero_412_dust,
36  1 aero_470_dust,
37  1 aero_650_dust,
38  1 aero_412_abs_dust,
39  1 aero_470_abs_dust
40 
41  include 'db.inc'
42  include 'newaottbl.inc'
43 
44  real realbuf(26), outbuf(21)
45  integer lm !land mask
46 c
47 c taking into account the terrain effect on rayleigh
48 c
49 c ---input parameters
50  real pdiff(3), Dstar1, tmpvg(7), bt11
51 c ---intermediate parameters
52  integer doy, ilat, ilon, ilat6, ilon6
53  real terrain_flag_new5
54  real trflg, pteran, x1, x2, xr470
55  real aot_mod(6),r470new, xday, r470ss2
56  real sza, scat_ang_2, amf
57  real outbufvg(21)
58  real model_frac, aod_frac
59 
60  real :: toa_ndvi, tmp412, tmp470, tmp670, ndvi_thold
61  real :: px_elev, bgaod
62  integer :: lc, gzflg, season, mod_sfc, gzflg_sav, regid
63  type(gdatetime) :: gdt1
64  integer :: status
65 
66  real :: r412_tbl, r470_tbl, r650_tbl, windsp, wv
67  real :: r412_135, r470_135, r650_135
68 
69  real :: tau_x412, tau_x412ss, tau_x412ss2, tau_x412_91
70  real :: tau_x412ss_995, tau_x412ss2_995, tau_x412ss_96
71  real :: tau_x412ss_97, tau_x412ss_94, tau_x412ss_95
72  real :: tau_x412ss_98, tau_x412ss_91, tau_smoke, tau_x412ss2_98
73 
74  real :: tau_x412_new_91, tau_x412_new_93, tau_x412_new_94
75  real :: tau_x412_new_96, tau_x412_new_995
76  real :: tau_x470_new_91, tau_x470_new_92, tau_x470_new_93
77  real :: tau_x470_new_94, tau_x470_new_95, tau_x470_new_96, tau_x470_new_995
78 
79  real :: tau_x412_dust_91, tau_x412_dust_93, tau_x412_dust_94
80  real :: tau_x412_dust_96, tau_x412_dust_995
81  real :: tau_x470_dust_91, tau_x470_dust_92, tau_x470_dust_93
82  real :: tau_x470_dust_94, tau_x470_dust_95, tau_x470_dust_96, tau_x470_dust_995
83 
84  real :: swir_coeffs412(3), swir_coeffs470(3)
85  real :: swir_stderr412, swir_stderr470
86  real :: swir_range(2)
87 
88  integer :: tau_x412_flag, tau_x412ss_flag, tau_x412ss2_flag, tau_x412_91_flag
89  integer :: tau_x412_flag_dust, tau_x470_flag_dust, tau_x470_flag_veg
90 
91  real :: tau_x470, tau_x470ss, tau_x470ss2, tau_x470_new
92  integer :: tau_x470_flag2, tau_x470ss_flag, tau_x470_new_flag, tau_x470_new_91_flag
93  real :: tau_x650
94  integer :: tau_x650_flag
95  integer :: tau_x412_flag2, tau_x470_flag, tau_x412_flag_91
96 
97  real :: refp_2100
98  real :: r470sv, r412sv, r470sv_veg, r412sv_veg
99  real :: tau_x412sv94, tau_x412sv96, tau_x412sv98
100  real :: tau_x470sv94, tau_x470sv96, tau_x470sv995
101  real :: tau_x470sv96_dust
102  integer :: tau_flagsv
103  logical :: abs_aero_flag
104 
105  logical :: sr_fail_flag, do_veg
106  integer :: alg_flag, brdf_flag
107 
108  logical dflag, dflag2
109  logical :: debug, use_alternate_brdf
110  integer :: lprint
111 
112 c ---output parameters
113  real xtau(3), alpha, ssa(3), tau550, sfc_typ
114  integer qa_flag(4)
115 c ---common parameter
116  common /xday/ doy
117 
118  integer handle_lut_out_of_bounds
119 
120  lprint = -999
121  debug = .false.
122 
123 c-----------------------------------------------------------------------
124 c initialization
125 c-----------------------------------------------------------------------
126 
127  mm = 10 ! solar zenith
128  nn = 46 ! satellite zenith
129  ll = 30 ! rel. azimuth
130  ma = 10 ! tau
131  mw = 8 ! ssa
132  tau_x470_flag_veg = -999 ! out-of-bound flag
133 
134 c-----------------------------------------------------------------------
135 c start processing the data
136 c-----------------------------------------------------------------------
137 
138  use_alternate_brdf = .false.
139 
140 c load the input data into local storage
141 
142  xlat = realbuf(1)
143  xlong = realbuf(2)
144  sza = realbuf(3)
145  xthet = realbuf(4)
146  xphi = realbuf(5)
147  ref650 = realbuf(6)
148  stdv = realbuf(13)
149  toa_ndvi = realbuf(14)
150  bt11 = realbuf(15)
151  dstar1 = realbuf(16)
152  pdiff(1) = realbuf(18)
153  pdiff(2) = realbuf(19)
154  pdiff(3) = realbuf(20)
155  pteran = realbuf(21)
156 
157 ! dump debug output
158 c if (abs(xlat-32.64) < 0.05 .AND. abs(xlong-(-114.5827)) < 0.05 .and. gzflg eq 31) then
159 c debug = .true.
160 c lprint = 1
161 c end if
162 ! if (abs(xlat-12.75) > 0.05 .OR. abs(xlong-23.92) > 0.07) go to 10
163 
164  if (xphi.gt.179.99) go to 10
165  if (xphi.lt.6.0) xphi = 6.
166 
167 c -- sun glint mask
168  cc = 3.14159/180.
169  psi = acos(cos(sza*cc)*cos(xthet*cc) +
170  1 sin(sza*cc)*sin(xthet*cc)*cos(xphi*cc))
171  glint_ang = psi/cc
172 c if (abs(psi/cc).lt.35.0) go to 10
173 
174 c -- scattering angle(scat_ang)
175  cc = 3.14159/180.
176  psi = acos(cos(sza*cc)*cos(xthet*cc) -
177  1 sin(sza*cc)*sin(xthet*cc)*cos(xphi*cc))
178  scat_ang = 180. - psi/cc
179 c if (scat_ang .gt. 175.) go to 10
180 
181 c -- air mass factor
182  amf = 1.0/cos(sza*cc)+1.0/cos(xthet*cc)
183 
184  ilat = floor(xlat*10.0) + 900 + 1
185  if (ilat > 1800) ilat = 1800
186  if (ilat < 1) ilat = 1
187 
188  ilon = floor(xlong*10.0) + 1800 + 1
189  if (ilon > 3600) ilon = 3600
190  if (ilon < 1) ilon = 1
191 
192  ilat6 = floor(xlat/0.06) + 1500 + 1
193  if (ilat > 3000) ilat = 3000
194  if (ilat < 1) ilat = 1
195 
196  ilon6 = floor(xlong/0.06) + 3000 + 1
197  if (ilon > 6000) ilon = 6000
198  if (ilon < 1) ilon = 1
199 
200  season = season_from_doy(2005, doy) ! year isn't currently available in here.
201  gdt1 = gregorian_from_doy(2005,doy) ! just use 2005 for now.
202 
203  trflg = terrain_flag(ilon,ilat)
204  terrain_flag_new5 = terrain_flag_new(ilon,ilat)
205  gzflg = get_geographic_zone(xlat, xlong, status)
206  if (status /= 0) then
207  print *, "ERROR: Failed to get geographic zone: ", xlat, xlong, status
208  return
209  end if
210 
211  gzflg_sav = gzflg
212 
213  regid = regid_2(ilon,ilat) ! 15 March 2018 JLee added
214 
215  ! 9 January 2018 JLee, read bg_aod
216  bgaod = get_background_aod(xlat, xlong, season, status)
217  if (status /= 0) then
218  print *, "ERROR: Failed to get background aod: ", xlat, xlong, status
219  return
220  end if
221 
222  xr470 = get_ler470(ilat, ilon, toa_ndvi, scat_ang, xphi)/100.0
223 
224  do i = 1, 3
225  if (pdiff(i).lt.0.0.and.pdiff(i).gt.-1.e-4)
226  1 pdiff(i) = 0.0
227 ! if (trflg.lt.0.0) pdiff(i) = 0.0
228  enddo
229  if (xr470.lt.0.) pdiff(2) = 0.
230 
231 ! -- disable pressure corrections
232 ! pdiff(:) = 0.0
233 
234  x1 = sza
235  x2 = xthet
236  refl1 = realbuf(11) ! 412 nm
237  refl3 = realbuf(9) ! 470 nm
238  refl6 = realbuf(7) ! 650 nm
239  refl11 = realbuf(12) ! 2250 nm
240  refl1_newgc = realbuf(10) ! 412 nm
241  refl3_newgc = realbuf(8) ! 470 nm
242 
243 
244 ! -- apply pressure correction to TOA reflectance over Morocco and West Sudan.
245  if (px_elev > 750.0 .AND. (xlat > 28.0 .AND. xlat < 37.0 .AND. xlong > -12.0 .AND.
246  1 xlong < 10.0)) then
247  refl3 = realbuf(9) + pdiff(2)/1.30 ! 470 nm
248  end if
249 
250  if (px_elev > 900.0 .AND. (xlat > 10.5 .AND. xlat < 19.5 .AND. xlong > 20.5 .AND.
251  1 xlong < 29.0)) then
252  refl3 = realbuf(9) + pdiff(2)/1.50 ! 470 nm
253  end if
254 
255  if (refl6 > 0.2) go to 10
256 
257  rr412_mod = realbuf(22)
258  rr470_mod = realbuf(23)
259  band26 = realbuf(17)
260  rat = refl6 / refl1
261  rat1 = rr470_mod / rr412_mod
262  rat_650_470 = ref650 / rr470_mod
263  rat_650_412 = ref650 / rr412_mod
264 
265  if (debug) then
266  print *, 'find_v, in: ', xlat, xlong, sza, xthet, xphi, realbuf(11), realbuf(9), realbuf(7), dstar1, lm
267  print *, 'ler, 412, 470, 650: ', realbuf(22), realbuf(23), realbuf(6)
268  endif
269 
270  5 continue
271 c---------------------------------------------------------------
272 c initialization
273 c---------------------------------------------------------------
274 c --intermediate parameters
275  dflag = .false.
276  dflag2 = .false.
277  w0_x = -999.
278  w0_int = -999.
279  w0_int_470 = -999.
280  w0_x412 = -999.
281  w0_x470 = -999.
282 
283  tau_x412 = -999.
284  tau_x412_91 = -999.
285  tau_x412ss = -999.
286  tau_x412ss_995 = -999.0
287  tau_x412ss2_995 = -999.0
288  tau_x412ss2_98 = -999.0
289  tau_x412ss_96 = -999.0
290  tau_x412ss_97 = -999.0
291  tau_x412ss_94 = -999.0
292  tau_x412ss_95 = -999.0
293  tau_x412ss_98 = -999.0
294  tau_x412ss_91 = -999.0
295  tau_smoke = -999.0
296  tau_x412_new_91 = -999.0
297  tau_x412_new_93 = -999.0
298  tau_x412_new_94 = -999.0
299  tau_x412_new_96 = -999.0
300  tau_x412_new_995 = -999.0
301  tau_x470_new = -999.0
302  tau_x470_new_91 = -999.0
303  tau_x470_new_92 = -999.0
304  tau_x470_new_93 = -999.0
305  tau_x470_new_94 = -999.0
306  tau_x470_new_95 = -999.0
307  tau_x470_new_96 = -999.0
308  tau_x470_new_995 = -999.0
309 
310  tau_x470 = -999.
311  tau_x470ss = -999.
312  tau_x470ss2 = -999.
313  tau_x650 = -999.
314  tau_x412_flag = -999
315  tau_x470_flag = -999
316  tau_x650_flag = -999
317  tau_x412_flag2 = -999
318  tau_x470_flag2 = -999
319  tau_x412_flag_91 = -999
320  xxrat = -999.
321  xxrat2 = -999.
322  aot = -999.
323  r412 = -999.0
324  r412new = -999.
325  r412ss = -999.0
326  r412ss2 = -999.0
327  r470 = -999.0
328  r470ss = -999.0
329  r470ss2 = -999.0
330  r470new = -999.
331  r470_sav = -999.0
332  r650 = -999.0
333  tau_x412_new = -999.
334  tau_x470_new = -999.
335 
336  tau_x412ss2 = -999.
337  tau_x470_new_91 = -999.
338  sr_fail_flag = .false.
339  abs_aero_flag = .false.
340  mod_sfc = -999
341 
342 c -- 2.2 um surface database aod parameters
343  r412sv = -999.0
344  r470sv = -999.0
345  r412sv_veg = -999.0
346  r470sv_veg = -999.0
347  tau_x412sv94 = -999.0
348  tau_x412sv96 = -999.0
349  tau_x412sv98 = -999.0
350  tau_x470sv94 = -999.0
351  tau_x470sv96 = -999.0
352  tau_x470sv995 = -999.0
353  tau_x470sv96_dust = -999.0
354  tau_flagsv = -999
355  swir_coeffs412(:) = -999.0
356  swir_coeffs470(:) = -999.0
357  swir_stderr412 = -999.0
358  swir_stderr470 = -999.0
359  swir_range(:) = -999.0
360 
361 c -- output parameters
362  tau550 = -999.
363  alpha = -999.
364  sfc_typ = -999.
365  alg_flag = 0
366 
367  do i = 1, 3
368  xtau(i) = -999.
369  ssa(i) = -999.
370  enddo
371 
372  do i = 1, 4
373  qa_flag(i) = 0
374  enddo
375 
376  do i = 1, 6
377  aot_mod(i) = -999.
378  enddo
379 
380  do i = 1,21
381  outbuf(i) = -999.
382  enddo
383  outbufvg(:) = -999.0
384 
385  xday = real(doy)
386 
387 c---------------------------------------------------------------
388 c screen for pixels outside reasonable ranges of reflectance
389 c---------------------------------------------------------------
390 c if (refl1.gt.0.0.and.refl1.lt.0.09.and.
391 c 1 refl6.gt.0.0.and.refl6.lt.0.14) go to 11
392 c if (refl1.gt.0.09.and.refl1.lt.0.50.and.
393 c 1 res.gt.6.0) go to 11
394 c go to 10
395 
396 11 continue
397 
398 ! -- over land, skip pixels with SZA > 72.0!
399  if (sza > 84.0) go to 10
400 
401 c--------------------------------------------------------
402 c load surface reflectance
403 c--------------------------------------------------------
404 c -- get base surf. reflc. values from surf. coeff. tables and save.
405  r412_tbl = get_ler412(ilat, ilon, toa_ndvi, scat_ang, xphi)
406  r470_tbl = get_ler470(ilat, ilon, toa_ndvi, scat_ang, xphi)
407  r650_tbl = get_ler650(ilat, ilon, toa_ndvi, scat_ang, xphi)
408  r865_tbl = get_modis_ler865(ilat, ilon) ! No 865 yet for VIIRS. Just use MODIS value.
409 
410  r412_135 = get_ler412(ilat, ilon, toa_ndvi, 135.0, xphi)
411  r470_135 = get_ler470(ilat, ilon, toa_ndvi, 135.0, xphi)
412  r650_135 = get_ler650(ilat, ilon, toa_ndvi, 135.0, xphi)
413  r865_135 = get_modis_ler865(ilat, ilon)
414 
415 c -- set to surface reflectance values to default table values.
416  r412 = r412_tbl
417  r470 = r470_tbl
418  r650 = r650_tbl
419  r865 = r865_tbl
420 
421  if (r865.lt.12.0.and.glint_ang.lt.30.0) go to 10 ! sun glint mask
422 
423  if (debug) print *, 'glint, scat. ang, r865: ', glint_ang, scat_ang, r865
424  if (debug) print *, 'ilat, ilon, toa_ndvi, xphi: ', ilat, ilon, toa_ndvi, xphi
425  if (debug) print *, "r412, r470, r650: ", r412, r470, r650
426 
427 c -- out of scope, skip.
428 c if (r412.gt.30.0) go to 10
429 c if (r470.gt.50.0) go to 10
430 c if (r650.gt.60.0) go to 10
431 
432 ! -- initialize other surface variables
433  r470ss = r470_tbl
434  r470new = r470_tbl
435 
436  r412ss = r412_tbl
437  r412ss2 = r412_tbl
438  r412new = r412_tbl
439 
440 !---------------------------------------------------------------------------------------------------
441 ! START NEW AERONET-DERIVED SURFACE REFLECTIVITY
442 !---------------------------------------------------------------------------------------------------
443 ! -- get landcover, then surface reflectivity for pixel..
444  lc = get_landcover(xlat, xlong, status)
445  if (status .ne. 0) then
446  print *, "ERROR: Failed to get land cover value: ", i, j, xlat, xlong, status
447  return
448  end if
449 
450 ! -- turn off tropical Sahel transition zone (zone 27) if it's not winter.
451  if (gzflg == 27 .AND. (xday >= 60 .AND. xday < 335)) gzflg = 5
452  if (gzflg == 27 .AND. realbuf(22) > 8.0) gzflg = 26
453 
454  brdf_flag = -1
455  tmp412=-999.0 ; tmp470=-999.0 ; tmp670=-999.0
456  brdf_flag = get_brdfcorr_sr(xlat, xlong, xphi, scat_ang, xthet, amf, px_elev, gdt1%month,
457  & toa_ndvi, stdv, gzflg, lc, bgaod, tmp412, tmp470, tmp670, use_alternate_brdf=use_alternate_brdf,
458  & debug=debug)
459  if (brdf_flag == 0 .OR. brdf_flag == 1) then
460  r412 = tmp412
461  r470 = tmp470
462  if (tmp670 > -900.0) then
463  r650 = tmp670
464  end if
465  end if
466  if (debug) then
467  print *, 'surface reflc, 412, 490, 670: ', r412, r470, r650
468  end if
469 
470  if (gzflg == 31 .and. px_elev >= 750.0) gzflg = 13
471 
472  if (gzflg == 18 .and. toa_ndvi < 0.3 .and. lc /= 6 .and. lc /= 4 .and.
473  1 realbuf(22) < 12.0 .and. realbuf(22)/realbuf(23) < 0.8) then ! modified by CH 2/3/2017
474  r412 = r412 * 1.6
475  r470 = r470 * 1.3
476  endif
477 
478 
479 ! -- insert special function to calculate surface reflectivity over N.Africa coast
480 ! -- e.g. Blida, Saada -- geozone 2. Overriding values from BRDF in modis_surface.f95.
481 ! if (gzflg == 2) then
482 ! mod_sfc = 8
483 ! call rsfc470(dflag2,mod_sfc,xday,xthet,xphi,scat_ang,r470)
484 ! if (dflag2) go to 10
485 ! endif
486 
487 ! -- special function surface reflectivity for Hamim and all of Arabian Peninsula.
488  if (gzflg >= 6 .AND. gzflg <= 11 .AND. gzflg /= 10) then
489  if (r650_135 > 32.0 .AND. (r650_135/r412_135) > 3.7) then
490 
491 ! -- r412_135+0.5 == value from table was a bit too low here. Adjust upwards to
492 ! -- bring the AOT back down.
493  mod_sfc = 9
494  call newsfc412_arab(dflag2,mod_sfc,xday,ilat,ilon,xthet,xphi,
495  & scat_ang,terrain_flag_new5,r412_135,r412new)
496  if ((r412_135 < 7.25) .OR. r650_135 > 42.0 .OR. dstar1 >1.03 .OR. use_alternate_brdf) then
497  continue
498  else
499  r412 = r412new
500  end if
501  end if
502  end if
503 
504 
505 ! -- use surface table value as last resort. Skip if undefined.
506  if (r412 < -900.0) r412 = r412_tbl
507  if (r470 < -900.0) r470 = r470_tbl
508  if (r650 < -900.0) r650 = r650_tbl
509 
510 ! -- check if we have any surface information.
511 ! -- if not, just skip the pixel if we have an AERONET BRDF but still
512 ! -- failed (brdf_flag == -2).
513 ! -- Otherwise, set sr_fail_flag and try the veg retrieval.
514  if (r412 < -900.0 .OR. r470 < -900.0 .OR. r650 < -900.0) then
515  if (brdf_flag == -2) then
516  go to 10
517  else
518  sr_fail_flag = .true.
519  go to 637 ! go to vege. retrieval, no need for SR database values.
520  end if
521  end if
522 
523 ! -- 2.2 um surface database approach first. Otherwise, no retrieval
524 ! results for certain conditions due to a lot of go to statement
525 ! below
526 
527 ! jlee added to test swir vs. vis surface table for 488 nm surface
528 ! reflectance
529  ! output: r470sv, 488 nm surface reflectance
530 
531  !if(xlong.gt.-20.0.and.xlong.lt.70.0.and.xlat.gt.5.0.and.xlat.lt.45.0) then
532  swir_coeffs412 = get_swir_coeffs412(ilat6,ilon6)
533  swir_stderr412 = get_swir_stderr412(ilat6,ilon6)
534  swir_coeffs470 = get_swir_coeffs470(ilat6,ilon6)
535  swir_stderr470 = get_swir_stderr470(ilat6,ilon6)
536  swir_range = get_swir_range(ilat6,ilon6)
537 
538 ! print *, swir_coeffs(:), swir_stderr, swir_range(:)
539 
540  refp_2100 = refl11*3.14159/cos(sza*cc)*100 !reflectance in percentage, not I over F
541  if(refp_2100.ge.(swir_range(1)-2).and.refp_2100.le.(swir_range(2)+2)) then
542  if(swir_coeffs470(1).gt.-900.0.and.swir_coeffs470(2).gt.-900.0.and.
543  1 swir_coeffs470(3).gt.-900.0.and.swir_stderr470.lt.1.0) then
544 
545  r470sv = swir_coeffs470(1)
546  1 + swir_coeffs470(2)*refp_2100
547  1 + swir_coeffs470(3)*refp_2100*refp_2100
548  r470sv_veg = r470sv
549  endif
550  !==========
551  if(swir_coeffs412(1).gt.-900.0.and.swir_coeffs412(2).gt.-900.0.and.
552  1 swir_coeffs412(3).gt.-900.0.and.swir_stderr412.lt.1.0) then
553 
554  r412sv = swir_coeffs412(1)
555  1 + swir_coeffs412(2)*refp_2100
556  1 + swir_coeffs412(3)*refp_2100*refp_2100
557  r412sv_veg = r412sv
558  endif
559  endif
560  !endif
561 
562  ! AOD using 490 nm band
563  if (r470sv > 0.0 .and. r470sv < 24.0) then
564  if (r470sv < 1.0) r470sv = 1.0
565  !========
566  imod = 2 ! w0 = 0.94
567  call aero_470(dflag,refl3_newgc,x1,x2,x3,mm,nn,ll,ma,
568  1 imod,r470sv,tau_x470sv94,tau_flagsv,trflg,0.0,debug)
569  if (dflag) go to 10
570  status = handle_lut_out_of_bounds(gzflg,tau_flagsv,tau_x470sv94)
571  if (status /= 0) then
572  print *,
573  1 "ERROR: Failed to check/reset AOT out of bounds condition: ",status
574  return
575  endif
576  !========
577  imod = 3 ! w0 = 0.96
578  call aero_470(dflag,refl3_newgc,x1,x2,x3,mm,nn,ll,ma,
579  1 imod,r470sv,tau_x470sv96,tau_flagsv,trflg,0.0,debug)
580  if (dflag) go to 10
581  status = handle_lut_out_of_bounds(gzflg,tau_flagsv,tau_x470sv96)
582  if (status /= 0) then
583  print *,
584  1 "ERROR: Failed to check/reset AOT out of bounds condition: ",status
585  return
586  endif
587  !========
588  imod = 4 ! w0 = 0.995
589  call aero_470(dflag,refl3_newgc,x1,x2,x3,mm,nn,ll,ma,
590  1 imod,r470sv,tau_x470sv995,tau_flagsv,trflg,0.0,debug)
591  if (dflag) go to 10
592  status = handle_lut_out_of_bounds(gzflg,tau_flagsv,tau_x470sv995)
593  if (status /= 0) then
594  print *,
595  1 "ERROR: Failed to check/reset AOT out of bounds condition: ",status
596  return
597  endif
598  !========
599  endif
600 
601  ! AOD using 412 nm band
602  if (r412sv > 0.0 .and. r412sv < 20.0) then
603  if (r412sv < 1.0) r412sv = 1.0
604  imod = 5 ! w0 = 0.94
605  call aero_412(dflag,refl1_newgc,x1,x2,x3,mm,nn,ll,ma,
606  1 imod,r412sv,tau_x412sv94,tau_flagsv,trflg,0.0,debug)
607  if (dflag) go to 10
608  status = handle_lut_out_of_bounds(gzflg,tau_flagsv,tau_x412sv94)
609  if (status /= 0) then
610  print *,
611  1 "ERROR: Failed to check/reset AOT out of bounds condition: ",status
612  return
613  endif
614  !========
615  imod = 6 ! w0 = 0.96
616  call aero_412(dflag,refl1_newgc,x1,x2,x3,mm,nn,ll,ma,
617  1 imod,r412sv,tau_x412sv96,tau_flagsv,trflg,0.0,debug)
618  if (dflag) go to 10
619  status = handle_lut_out_of_bounds(gzflg,tau_flagsv,tau_x412sv96)
620  if (status /= 0) then
621  print *,
622  1 "ERROR: Failed to check/reset AOT out of bounds condition: ",status
623  return
624  endif
625  !========
626  imod = 7 ! w0 = 0.98
627  call aero_412(dflag,refl1_newgc,x1,x2,x3,mm,nn,ll,ma,
628  1 imod,r412sv,tau_x412sv98,tau_flagsv,trflg,0.0,debug)
629  if (dflag) go to 10
630  status = handle_lut_out_of_bounds(gzflg,tau_flagsv,tau_x412sv98)
631  if (status /= 0) then
632  print *,
633  1 "ERROR: Failed to check/reset AOT out of bounds condition: ",status
634  return
635  endif
636  endif
637 
638  !r470sv adjustment for tau_x470sv96_dust
639  ddx = 0.
640  if (px_elev < 500.0) then
641  if (xphi < 90.0 .and. r470sv > 0.0 .and. r470sv < 24.0) then
642  ddx = xthet*1.5/65.0
643  if (r412_tbl > 12.0) ddx = xthet*0.5/65.0
644  if (r412_tbl > 14.0) ddx = 0.0
645  r470sv = r470sv - ddx
646  endif
647  if (xphi >= 90.0 .and. r412_tbl > 12.0 .and. r470sv > 0.0 .and. r470sv < 24.0) then
648  ddx = xthet*0.8/65.0
649  r470sv = r470sv + ddx
650  endif
651  endif
652 
653  ! output: tau_x470sv, 470 nm AOD from default aerosol model (SSA = 0.96)
654  if (r470sv > 0.0 .and. r470sv < 24.0) then
655  if (r470sv < 1.0) r470sv = 1.0
656  imod = 3 ! w0 = 0.96, dust
657  call aero_470_dust(dflag,refl3_newgc,x1,x2,x3,mm,nn,ll,ma,
658  1 imod,r470sv,tau_x470sv96_dust,tau_flagsv,trflg,0.0,debug)
659  if (dflag) go to 10
660  status = handle_lut_out_of_bounds(gzflg,tau_flagsv,tau_x470sv96_dust)
661  if (status /= 0) then
662  print *,
663  1 "ERROR: Failed to check/reset AOT out of bounds condition: ", status
664  return
665  endif
666  endif
667 ! end jlee added
668 
669  r650_sav = r650
670  r470_sav = r470
671  r412_sav = r412
672 
673 c --- direct copy/paste from seawifs code
674 !-- 470 nm
675 
676  r470new = r470_sav
677 
678 ! -- start tweaking r412ss, r412ss2, r470ss---------------------------------------
679 ! -- new bi-directional surface 2
680  ddx = 1.8
681  if (xday > 150.0 .and. xday < 258.0) ddx = 1.0 ! Jun, Jul
682 
683  ddx2 = 1.0
684  if (xday > 150.0 .and. xday < 258.0) ddx2 = 1.5 ! Jun, Jul
685 
686  r470ss = r470_tbl
687  r470ss2 = r470_tbl
688 
689  dda1 = 2.0
690  dda2 = 1.0
691  if (scat_ang >= 100.0 .and. scat_ang < 155.0)
692  1 r470ss = r470_tbl + dda1*(scat_ang-100.)/55.
693  if (scat_ang >= 155.0)
694  1 r470ss = r470_tbl+dda1
695  r470ss = r470ss - ddx - ddx2
696 
697  dda1 = 2.0
698  dda2 = 1.0
699  if (r470_tbl > 9.0) then
700  if (scat_ang >= 100.0 .and. scat_ang < 155.0)
701  1 r470ss = r470_tbl + dda1*(scat_ang-100.)/55.
702  if (scat_ang >= 155.0)
703  1 r470ss = r470_tbl+dda1
704  r470ss = r470ss - ddx - ddx2
705  endif
706 
707 ! Libya - Egypt 1
708  if (xlat >21. .and. xlong > 12.0) then
709  dda1 = 2.5
710  if (r412_tbl > 12.0) dda1 = 3.0
711  dda2 = 1.0
712  if (r412_tbl > 10.0) then
713  if (scat_ang >= 100.0 .and. scat_ang < 155.0)
714  1 r470ss = r470_tbl + dda1*(scat_ang-100.)/55.
715  if (scat_ang >= 155.0)
716  1 r470ss = r470_tbl+dda1
717  r470ss = r470ss - ddx
718  endif
719  endif
720 
721 ! Libya - Egypt 2
722  if (xlat >15. .and. xlong > 22.0) then
723  dda1 = 2.5
724  if (r412_tbl > 12.0) dda1 = 3.0
725  dda2 = 1.0
726  if (r412_tbl > 10.0) then
727  if (scat_ang >= 100.0 .and. scat_ang < 155.0)
728  1 r470ss = r470_tbl + dda1*(scat_ang-100.)/55.
729  if (scat_ang >= 155.0)
730  1 r470ss = r470_tbl+dda1
731  r470ss = r470ss - ddx
732  endif
733  endif
734 
735 ! Central Algeria
736  if (xday > 243.0 .and. xday < 274.0) then ! Sept
737  if (xlat >22.5 .and. xlat <30. .and. xlong > -4.0.and. xlong < 1.) then
738  dda1 = 2.5
739  if (r412_tbl > 10.) then
740  if (scat_ang >= 100.0 .and. scat_ang < 155.0)
741  1 r470ss = r470_tbl + dda1*(scat_ang-100.)/55.
742  if (scat_ang >= 155.0)
743  1 r470ss = r470_tbl+dda1
744  r470ss = r470ss - ddx - ddx2
745  endif
746  endif
747  endif
748 
749 ! N. Algeria 1
750  if (xlat >31.5 .and. xlong > 4. .and. xlong < 10.) then
751  dda1 = 1.5
752  dda2 = 1.0
753  if (r412_tbl > 9.0) then
754  if (scat_ang >= 100.0 .and. scat_ang < 155.0)
755  1 r470ss = r470_tbl + dda1*(scat_ang-100.)/55.
756  if (scat_ang >= 155.0)
757  1 r470ss = r470_tbl+dda1
758  r470ss = r470ss - ddx
759  endif
760  endif
761 
762 ! NE Mauritania 1
763  if (xlat>22.5 .and. xlat<30.0.and. xlong>-11.0.and. xlong< -5.) then
764  if (r412_tbl > 9.) then
765  r470ss = r470ss + 0.0
766  endif
767  endif
768 
769 ! NE Mauritania 2
770  if (xlat>22.5 .and. xlat<26.0.and. xlong>-13.0.and. xlong< -11.001) then
771  if (r412_tbl > 9.) then
772  r470ss = r470ss + 0.0
773  endif
774  endif
775 
776 ! NE Mauritania 3
777  if (xlat>20.0 .and. xlat<30.0.and. xlong>-20.0.and. xlong< -12.0) then
778  if (r412_tbl > 9.) then
779  r470ss = r470ss + 0.0
780  endif
781  endif
782 
783 ! NE Mauritania 3
784  if (xlat>15.0 .and. xlat<20.0 .and. xlong< -14.9) then
785  if (r412_tbl > 10.5) then
786  if (xday >= 244.0 .and. xday < 274.0) r470ss = r470_tbl - ddx - ddx2
787  endif
788  endif
789 
790 
791  if (xphi > 90.0) then
792  dda1 = 0.0
793  dda2 = 1.4
794  if (xday > 150.0 .and. xday < 258.0) dda2 = 2.0
795  if (xthet >= 30.0 .and. xthet < 70.0) dda1 = (xthet-30.0) * dda2/40.
796  if (xthet >= 70.0) dda1 = dda2
797  r470ss = r470ss + dda1
798  endif
799 
800  if (dstar1 > 1.1) r470ss = r470_tbl - ddx - ddx2
801  if (dstar1 > 1.06.and.xday > 150.0 .and. xday < 258.0) r470ss = r470_tbl - ddx - ddx2
802  if (dstar1 > 1.04.and.wv>1.5.and.xday > 150.0 .and. xday < 258.0) r470ss = r470_tbl - ddx - ddx2
803 
804  if (gzflg >= 6 .and. gzflg <= 11) then
805 
806  ddx = 0.5 +0.8
807  if (xday > 59.0 .and. xday < 151.0) ddx = 0.5+0.5 ! Mar, Apr, May
808  if (xday > 150.0 .and. xday < 258.0) ddx = 1.0+1.5 ! Jun, Jul, Aug
809  if (xday > 243.0 .and. xday < 335.0) ddx = 0.7+1.5 ! Sep, Oct, Nov
810  if (xday > 243.0 .and. xday < 335.0.and.r412_tbl >11.0) ddx = 0.7+1.0 ! Sep, Oct, Nov
811 
812  ddx2 = 0.0
813  if (xday > 150.0 .and. xday < 258.0.and.r412_tbl <11.0) ddx2 = 0.5 ! Jun, Jul, Aug
814 
815  r470ss = r470 - ddx
816  r470ss2 = r470
817 
818  dda1 = r470_tbl*0.1 +0.5
819  dda2 = r470_tbl*0.0
820  if (scat_ang >= 90.0 .and. scat_ang < 170.0)
821  1 r470ss = r470_tbl + dda1*(scat_ang-90.)/80.
822  if (scat_ang >= 170.0)
823  1 r470ss = r470_tbl+dda1+dda2*(scat_ang-170.)/10.
824  r470ss = r470ss - ddx - ddx2
825 
826  dda1 = r470_tbl*0.1 +1.0
827  dda2 = r470_tbl*0.0
828  if (r412_tbl > 11.0) dda1 = r470_tbl*0.1 +1.0
829  if (r412_tbl > 12.0) dda1 = r470_tbl*0.1 +1.0
830  if (xday > 120.0 .and. xday < 152.0 .and. r412_tbl > 11.0) dda1 = dda1 - 1.5 ! May
831  if (xday > 59.0 .and. xday < 244.0.and.r412_tbl > 10.3) dda1 = r470_tbl*0.1 +1.5 ! spring, summer
832  if (xday > 243.0 .and. xday < 335.0.and.r412_tbl > 11.0) dda1 = r470_tbl*0.1 +1.0 ! fall
833 
834  if (r412_tbl > 9.0) then
835  scat_ang_2 = scat_ang
836  if (scat_ang_2 >= 90.0 .and. scat_ang_2 < 160.0)
837  1 r470ss = r470_tbl + dda1*(scat_ang_2-90.)/70.
838  if (scat_ang_2 >= 160.0)
839  1 r470ss = r470_tbl+dda1+dda2*(scat_ang_2-160.)/20.
840  r470ss = r470ss - ddx - ddx2
841  endif
842 
843  if (xday > 150.0 .and. xday < 244.0) then ! Jun,Jul, Aug
844  if (dstar1 > 1.1) r470ss = r470_tbl
845  endif
846 
847  if (xday > 243.0 .and. xday < 274.0.and.dstar1 > 1.02) then ! Sept
848  r470ss = r470_tbl-1.0
849  endif
850  if (xday > 120.0 .and. xday < 152.0.and.dstar1 > 1.08) then ! May
851  r470ss = r470_tbl-1.0
852  endif
853 
854 ! Southern Arabian Pen
855  if (xday > 152.0 .and. xday < 244.0) then ! Jun, Jul, Aug
856 
857  ddx = 0.0
858  if (xday > 152.0 .and. xday < 182.) ddx = 0.5 ! Jun
859  if (xday > 181.0 .and. xday < 244.0) ddx = 1.5 ! Jul, Aug
860  dda1 = 1.5
861  dda2 = 0.0
862  if (r412_tbl > 11.0) dda1 = 1.8
863  if (r412_tbl > 12.0) dda1 = 2.0
864  if (xday > 181.0 .and. xday < 244.0) then ! Jul, Aug
865  if (r412_tbl > 11.0) dda1 = 1.8
866  if (r412_tbl > 12.0) dda1 = 2.0
867  endif
868 
869  if (xlat < 18.0.and. xlong <= 49.0) then
870  if (r412_tbl > 11.) then
871  if (scat_ang >= 100.0 .and. scat_ang < 175.0)
872  1 r470ss = r470_tbl-ddx + dda1*(scat_ang-100.)/75.
873  if (scat_ang >= 175.0)
874  1 r470ss = r470_tbl-ddx +dda1+dda2*(scat_ang-175.)/5.
875  endif
876  endif
877 
878  if (xlat < 18.5 .and. xlong > 49.0.and. xlong <= 52.5) then
879  if (r412_tbl > 11.) then
880  if (scat_ang >= 100.0 .and. scat_ang < 175.0)
881  1 r470ss = r470_tbl-ddx + dda1*(scat_ang-100.)/75.
882  if (scat_ang >= 175.0)
883  1 r470ss = r470_tbl-ddx +dda1+dda2*(scat_ang-175.)/5.
884  endif
885  endif
886 
887  if (xlat < 23.0 .and. xlong > 52.5) then
888  if (r412_tbl > 11.) then
889  if (scat_ang >= 100.0 .and. scat_ang < 175.0)
890  1 r470ss = r470_tbl-ddx + dda1*(scat_ang-100.)/75.
891  if (scat_ang >= 175.0)
892  1 r470ss = r470_tbl-ddx +dda1+dda2*(scat_ang-175.)/5.
893  endif
894  endif
895 
896  endif ! end of southern Arabian Pen
897 
898  dda2 = 1.01
899  if (xday > 151.0 .and. xday < 244.0) dda2 = 1.1
900  if (xphi > 90.0 .and. dstar1 <dda2) then
901  ddx = 0.0
902  dda1 = 1.4
903  if (xday > 151.0 .and. xday < 244.0) dda1 = 1.8
904  if (xday > 151.0 .and. xday < 244.0.and.r412_tbl > 10.3.and.(rr470_mod-r470ss) > 0.0) dda1 = 2.2
905  if (xday > 243.0 .and. xday < 335.0) dda1 = 1.8
906  ! if (xday > 243.0 .and. xday < 335.0) dda1 = 2.5
907  if (xday > 243.0 .and. xday < 335.0.and.r412_tbl > 10.3.and.(rr470_mod-r470ss) > 0.0) dda1 = 2.5
908  if (xday > 243.0 .and. xday < 335.0.and.r412_tbl > 11.) dda1 = 1.8
909  if (xthet >= 30.0 .and. xthet < 70.0) ddx = (xthet-30.0) * dda1/40.
910  if (xthet >= 70.0) ddx = dda1
911  r470ss = r470ss + ddx
912  endif
913 
914  endif
915 
916  if (r470 >= 24.0) r470 = 23.9
917  if (r470 < 1.0) r470 = 1.0
918  if (r470ss >= 24.0) r470ss = 23.9
919  if (r470ss < 1.0 .AND. (gzflg > 0 .AND. (gzflg <= 11 .OR. gzflg == 27))) go to 10
920 
921  if (lprint > 0) print *,'final r470,r470ss =', r470,r470ss
922 
923 !-- 412 nm
924 
925  ddx = 0.0
926  if (xday > 32.0 .and. xday < 60.0) ddx = 0.5 ! Feb
927  if (xday > 59.0 .and. xday < 152.0) ddx = 1.0 + 0.4+0.3 ! Mar,Apr, May
928  if (xday > 151.0 .and. xday < 213.0) ddx = 1.0 + 0.8+0.5 ! Jun, Jul
929  if (xday > 212.0 .and. xday < 244.0) ddx = 0.5+0.5 ! Aug
930  if (xday > 243.0 .and. xday < 335.0) ddx = 0.5+0.5 ! Sep, Oct, Nov
931  if (xday > 273.0 .and. xday < 305.0) ddx = 0.7+0.5 ! Oct
932 
933  ddx1 = 0.5
934  if (xday < 32.0) ddx1 = 0.5 ! Jan
935  if (xday > 59.0 .and. xday < 152.0) ddx1 = 0.7 ! Mar,Apr, May
936  if (xday > 151.0 .and. xday < 244.0) ddx1 = 0.5 ! Jun,Jul,Aug
937  if (xday > 243.0 .and. xday < 305.0) ddx1 = 0.5 ! Sep, Oct
938  if (xday > 304.0) ddx1 = 0.5 ! Nov, Dec
939 
940  ddx2 = 0.7
941  if (xday < 32.0) ddx2 = 0.7 ! Jan
942  if (xday > 59.0 .and. xday < 152.0) ddx2 = 0.7 ! Mar,Apr,May
943  if (xday > 243.0 .and. xday < 305.0) ddx2 = 0.7 ! Sep, Oct
944  if (xday > 304.0) ddx2 = 0.7 ! Nov, Dec
945 
946  ddx3 = 0.3
947  if (xday > 0.0 .and. xday < 152.0) ddx3 = 1.0 ! Jan,Feb,Mar,Apr,May
948  ddx33 = 0.4
949  if (xday > 32.0 .and. xday < 60.0) ddx33 = 0.7 ! Feb
950  if (xday > 59.0 .and. xday < 152.0) ddx33 = 1.0 ! Mar,Apr, May
951  if (xday > 151.0 .and. xday < 244.0) ddx33 = 0.7 ! Jun,Jul,Aug
952  if (xday > 243.0 .and. xday < 305.0) ddx33 = 0.4 ! Sep, Oct
953  if (xday > 304.0) ddx33 = 0.6 ! Nov, Dec
954 
955  r412ss2 = r412_tbl - ddx
956 
957  if (scat_ang >= 120.0 .and. scat_ang < 170.0)
958  1 r412ss2 = r412_tbl- ddx + ddx1*(scat_ang-120.)/50.
959  if (scat_ang >= 170.0)
960  1 r412ss2 = r412_tbl- ddx+ddx1+0.0*(scat_ang-170.)/10.
961 
962  if (r412_tbl > 9.0) then
963  if (scat_ang >= 120.0 .and. scat_ang < 170.0)
964  1 r412ss2 = r412_tbl- ddx + ddx2*(scat_ang-120.)/50.
965  if (scat_ang >= 170.0)
966  1 r412ss2 = r412_tbl- ddx+ddx2+0.0*(scat_ang-170.)/10.
967  endif
968 
969  if (r412_tbl > 11.0) then
970  if (scat_ang >= 110.0 .and. scat_ang < 160.0)
971  1 r412ss2 = r412_tbl - ddx+ ddx3*(scat_ang-110.)/50.
972  if (scat_ang >= 160.0)
973  1 r412ss2 = r412_tbl- ddx+ddx3+0.0*(scat_ang-160.)/20.
974  endif
975 
976  if (r412_tbl > 12.0) then
977  if (scat_ang >= 110.0 .and. scat_ang < 160.0)
978  1 r412ss2 = r412_tbl - ddx+ ddx33*(scat_ang-110.)/50.
979  if (scat_ang >= 160.0)
980  1 r412ss2 = r412_tbl- ddx+ddx33+0.0*(scat_ang-160.)/20.
981  endif
982 
983  if (xphi > 90.0) then
984  ddx = 0.0
985  ddx33 = 1.4
986  if (xday > 243.0 .and. xday < 305.0) ddx33 = 0.0
987  if (r412_tbl > 12.0) then
988  ddx33 = 2.0
989  if (xday > 243.0 .and. xday < 305.0) ddx33 = 1.5
990  endif
991  if (xthet < 70.0) ddx = xthet * ddx33/70.
992  if (xthet >= 70.0) ddx = ddx33
993  r412ss2 = r412ss2 + ddx
994  endif
995 
996 ! -- new bi-directional surface 2
997 
998  ddx = 0.0
999  if (xday > 0.0 .and. xday < 32.0) ddx = 0.5 + 0.7 ! Jan
1000  if (xday > 31.0 .and. xday < 60.0) ddx = 0.5 + 0.4 ! Feb
1001  if (xday > 59.0 .and. xday < 152.0) ddx = 0.5 + 0.4+0.8 ! Mar, Apr, May
1002  if (xday > 151.0 .and. xday < 244.0) ddx = 0.5 + 0.4+1.0 ! Jun, Jul, Aug
1003  if (xday > 243.0 .and. xday < 305.0) ddx = 0.5 + 0.4+1.0 ! Sep, Oct
1004  if (xday > 304.0 .and. xday < 335.0) ddx = 0.5 + 0.4 ! Nov
1005  if (xday > 334.0) ddx = 0.5 + 0.4 ! Dec
1006  ddx2 = 1.0
1007  if (xday > 59.0 .and. xday < 152.0) ddx2 = 1.0 ! Mar, Apr,May
1008  if (xday > 151.0 .and. xday < 244.0) ddx2 = 1.2 ! Jun,Jul,Aug
1009  if (xday > 243.0 .and. xday < 305.0) ddx2 = 0.7 ! Sep, Oct
1010  if (xday > 304.0 .and. xday < 335.0) ddx2 = 1.0 ! Nov
1011  ddx3 = 1.5
1012  if (xday < 32.0) ddx3 = 1.0
1013  if (xday > 151.0 .and. xday < 244.0) ddx3 = 1.5 ! Jun,Jul,Aug
1014  if (xday > 243.0 .and. xday < 305.0) ddx3 = 0.8 ! Sep, Oct
1015  ddx33 = 0.0
1016  if (xday > 151.0 .and. xday < 244.0) ddx33 = 0.0 ! Jun,Jul,Aug
1017 
1018  if (xday > 151.0 .and. xday < 244.0) then
1019  dda1 = rr412_mod - r412_tbl
1020  if (dda1 <-1.5 .and.dstar1 >1.0 ) ddx2 = 0. ! Jun,Jul,Aug
1021  endif
1022 
1023  r412ss = r412_tbl - ddx
1024 
1025  if (xday > 151.0 .and. xday < 274.0) then
1026  dd = (xlong - 8.0) /4.0
1027  if (xlong > 12.0) dd = 1.
1028  if (xlong <= 8.0) dd = 0.
1029  dda1 = ddx2
1030  ddx2 = dda1*dd
1031  endif
1032 
1033  if (scat_ang >= 100.0 .and. scat_ang < 160.0)
1034  1 r412ss = r412_tbl- ddx + ddx2*(scat_ang-100.)/60.
1035  if (scat_ang >= 160.0)
1036  1 r412ss = r412_tbl- ddx+ddx2+ddx33*(scat_ang-160.)/20.
1037 
1038  if (xday > 243.0 .and. xday < 305.0) then
1039  if (scat_ang >= 140.0 .and. scat_ang < 160.0)
1040  1 r412ss = r412_tbl- ddx + ddx2*(scat_ang-140.)/20.
1041  if (scat_ang >= 160.0)
1042  1 r412ss = r412_tbl- ddx+ddx2
1043  endif
1044 
1045  if (r412_tbl > 12.0) then
1046  if (scat_ang >= 120.0 .and. scat_ang < 170.0)
1047  1 r412ss = r412_tbl- ddx + ddx3*(scat_ang-120.)/50.
1048  if (scat_ang >= 170.0)
1049  1 r412ss = r412_tbl- ddx+ddx3+0.0*(scat_ang-170.)/10.
1050  endif
1051 
1052  if (xday > 151.0 .and. xday < 244.0) then
1053  if (r412_tbl > 10.0) r412ss = r412ss + 0.5
1054  endif
1055 
1056 ! E. Central Algeria, long narrow desert
1057  dda1 = 3.2
1058  dda2 = 9.5
1059  dda3 = 100.0
1060  dda4 = 0.5
1061  dda5 = 70.
1062  if (xday > 181.0 .and. xday < 305.0) then ! summer, Sep, Oct
1063  dda1 = 2.5
1064  dda2 = 9.5
1065  dda3 = 100.0
1066  dda4 = 0.0
1067  dda5 = 70.
1068  endif
1069  if (xday > 151.0 .and. xday < 182.0) dda1 = 1.5
1070  if (xday > 243.0 .and. xday < 274.0) dda1 = 0.7
1071  if (xday > 59.0 .and. xday < 152.0) dda1 = 1.0
1072 
1073  if (xlat>27. .and. xlat<36. .and. xlong>2.5 .and. xlong<11.5) then
1074  if (r412_135 > dda2) then
1075  if (scat_ang >= dda3 .and. scat_ang < 170.0)
1076  1 r412ss = r412_tbl- ddx + dda1*(scat_ang-dda3)/dda5
1077  if (scat_ang >= 170.0)
1078  1 r412ss = r412_tbl- ddx+dda1+dda4*(scat_ang-170.)/10.
1079  endif
1080  endif
1081 
1082 ! Libya - Egypt 1
1083  if (xlat >20. .and. xlong > 14.9) then
1084  dda1 = r412_tbl*0.08 +0.5
1085  if (r412_tbl > 12.0) dda1 = r412_tbl*0.08 + 0.8
1086  if (xday > 59.0 .and. xday < 91.0 .and.r412_tbl > 12.0) dda1 = r412_tbl*0.08 + 1.5 ! Mar
1087  if (xday > 151.0 .and. xday < 244.0 .and.r412_tbl > 12.0) dda1 = r412_tbl*0.08 + 0.4 ! summer
1088  dda2 = r412_tbl*0.0
1089  if (r412_tbl > 9.6) then
1090  if (scat_ang >= 110.0 .and. scat_ang < 160.0)
1091  1 r412ss = r412_tbl- ddx + dda1*(scat_ang-110.)/50.
1092  if (scat_ang >= 160.0)
1093  1 r412ss = r412_tbl- ddx+dda1+dda2*(scat_ang-160.)/20.
1094  if (xday > 151.0 .and. xday < 244.0) then
1095  if (xday > 151.0 .and. xday < 182.0) then
1096  if (r412_tbl >10.0.and.xlat>20.0 .and. xlat<24.8.and. xlong>24.7.and.xlong<30.0) go to 211
1097  if (r412_tbl >10.0.and.xlat>17.5 .and. xlat<=20.0.and. xlong>24.0.and.xlong<30.0) go to 211
1098  endif
1099  dda1 = r412_tbl- ddx+r412_tbl*0.08 +0.5
1100  if (scat_ang >= 110.0 .and. scat_ang < 160.0) r412ss = dda1+0.8*(scat_ang-110.)/50.
1101  if (scat_ang >= 160.0) r412ss = dda1+0.8
1102  if (r412_tbl > 14.0) then
1103  dda1 = r412_tbl- ddx+r412_tbl*0.08 +0.5
1104  if (scat_ang >= 110.0 .and. scat_ang < 160.0) r412ss = dda1+1.0*(scat_ang-110.)/50.
1105  if (scat_ang >= 160.0) r412ss = dda1+1.0
1106  endif
1107  endif
1108  endif
1109  endif
1110 211 continue
1111 
1112 ! Libya - Egypt 2
1113  if (xlat >15. .and. xlong > 22.0) then
1114  dda1 = r412_tbl*0.08 +0.5
1115  if (r412_tbl > 12.0) dda1 = r412_tbl*0.08 + 0.8
1116  if (xday > 59.0 .and. xday < 91.0 .and.r412_tbl > 12.0) dda1 = r412_tbl*0.08 + 1.5 ! Mar
1117  if (xday > 151.0 .and. xday < 244.0 .and.r412_tbl > 12.0) dda1 = r412_tbl*0.08 + 0.4 ! summer
1118  dda2 = r412_tbl*0.0
1119  if (r412_tbl > 9.6) then
1120  if (scat_ang >= 110.0 .and. scat_ang < 160.0)
1121  1 r412ss = r412_tbl- ddx + dda1*(scat_ang-110.)/50.
1122  if (scat_ang >= 160.0)
1123  1 r412ss = r412_tbl- ddx+dda1+dda2*(scat_ang-160.)/20.
1124  if (xday > 151.0 .and. xday < 244.0) then
1125  if (xday > 151.0 .and. xday < 182.0) then
1126  if (r412_tbl >10.0.and.xlat>20.0 .and. xlat<24.8.and. xlong>24.7.and.xlong<30.0) go to 212
1127  if (r412_tbl >10.0.and.xlat>17.5 .and. xlat<=20.0.and. xlong>24.0.and.xlong<30.0) go to 212
1128  endif
1129  dda1 = r412_tbl- ddx+r412_tbl*0.08 +0.5
1130  if (scat_ang >= 110.0 .and. scat_ang < 160.0) r412ss = dda1+0.8*(scat_ang-110.)/50.
1131  if (scat_ang >= 160.0) r412ss = dda1+0.8
1132  if (r412_tbl > 14.0) then
1133  dda1 = r412_tbl- ddx+r412_tbl*0.08 +0.5
1134  if (scat_ang >= 110.0 .and. scat_ang < 160.0) r412ss = dda1+1.0*(scat_ang-110.)/50.
1135  if (scat_ang >= 160.0) r412ss = dda1+1.0
1136  endif
1137  endif
1138  endif
1139  endif
1140 212 continue
1141 
1142 ! N. Algeria 1
1143  if (xlat >31.5 .and. xlong > 4. .and. xlong < 10.) then
1144  dda1 = r412_tbl*0.08 + 0.5
1145  dda2 = r412_tbl*0.05 + 0.5
1146  if (xday > 151.0 .and. xday < 244.0) then
1147  dda1 = r412_tbl*0.08 -0.3
1148  dda2 = 0.0
1149  endif
1150  if (xday > 243.0 .and. xday < 274.0) then
1151  dda1 = r412_tbl*0.08 -0.3
1152  dda2 = 0.0
1153  endif
1154  if (r412_tbl > 9.4) then
1155  if (scat_ang >= 120.0 .and. scat_ang < 170.0)
1156  1 r412ss = r412_tbl - ddx+0.3+ dda1*(scat_ang-120.)/50.
1157  if (scat_ang >= 170.0)
1158  1 r412ss = r412_tbl- ddx+0.3+dda1+dda2*(scat_ang-170.)/10.
1159  endif
1160  endif
1161 33 continue
1162 
1163 ! N. Algeria 2
1164  if (xlat >30.0 .and. xlat <32.0 .and. xlong > 4.8 .and. xlong < 7.) then
1165  dda1 = r412_tbl*0.08 + 0.5
1166  dda2 = 0.0
1167  if (xday > 151.0 .and. xday < 244.0) then
1168  dda1 = r412_tbl*0.08 -0.3
1169  dda2 = 0.0
1170  endif
1171  if (xday > 243.0 .and. xday < 274.0) then
1172  dda1 = r412_tbl*0.08 -0.3
1173  dda2 = 0.0
1174  endif
1175  if (r412_tbl > 9.4) then
1176  if (scat_ang >= 100.0 .and. scat_ang < 160.0)
1177  1 r412ss = r412_tbl- ddx+0.3 + dda1*(scat_ang-100.)/60.
1178  if (scat_ang >= 160.0)
1179  1 r412ss = r412_tbl- ddx+0.3+dda1+dda2*(scat_ang-160.)/20.
1180  endif
1181  endif
1182 
1183 ! Chad - Libya border
1184  if (xday > 31.0 .and. xday < 60.0) go to 36 ! use for all months except for Feb
1185  if (xlat >20. .and. xlat <25. .and. xlong >15.0.and. xlong <17.5) then
1186  if (r412_tbl > 12.0) then
1187  r412ss = r412ss + 0.5
1188  endif
1189  endif
1190 36 continue
1191 
1192  dda2 = 0.5
1193  if (xday > 59.0 .and. xday < 121.0) dda2 = 1.0
1194 ! NE Mauritania 1
1195  if (xlat>22.5 .and. xlat<30.0.and. xlong>-11.0.and. xlong< -5.) then
1196  if (r412_tbl > 9.) then
1197  r412ss = r412ss + dda2
1198  if (xday > 151.0 .and. xday < 274.0) r412ss = r412_tbl-1.5 ! summer
1199  endif
1200  endif
1201 
1202  dda2 = 0.5
1203  if (xday > 59.0 .and. xday < 121.0) dda2 = 1.0
1204 ! NE Mauritania 2
1205  if (xlat>22.5 .and. xlat<26.0.and. xlong>-13.0.and. xlong< -11.001) then
1206  if (r412_tbl > 9.) then
1207  r412ss = r412ss + dda2
1208  if (xday > 151.0 .and. xday < 274.0) r412ss = r412_tbl-1.5 ! summer
1209  endif
1210  endif
1211 
1212  dda2 = 0.0
1213  if (xday > 59.0 .and. xday < 121.0) dda2 = 1.0
1214 ! NE Mauritania 3
1215  if (xlat>=20.0 .and. xlat<29.0.and. xlong< -12.5) then
1216  if (r412_tbl > 10.5) then
1217  r412ss = r412ss + dda2
1218  if (xday > 151.0 .and. xday < 274.0) r412ss = r412_tbl-1.5 ! summer
1219  endif
1220  endif
1221 
1222  dda2 = 0.0
1223  if (xday > 59.0 .and. xday < 121.0) dda2 = 1.0
1224 ! NE Mauritania 3
1225  if (xlat>15.0 .and. xlat<20.0 .and. xlong< -14.9) then
1226  if (r412_tbl > 10.5) then
1227  r412ss = r412ss + dda2
1228  if (xday > 151.0 .and. xday < 244.0) r412ss = r412_tbl-1.5 ! summer
1229  if (xday >= 244.0 .and. xday < 274.0) r412ss = r412ss-1.5
1230  endif
1231  endif
1232 
1233 ! Lake Chad
1234  if (xlat >10.0 .and. xlat <21.0 .and. xlong > 10.0 .and. xlong < 20.0) then
1235  if (r412_tbl > 12.) then
1236  if (xday > 181.0 .and. xday < 244.0) r412ss = r412_tbl ! summer
1237  endif
1238  endif
1239 
1240 ! Morocco 1
1241  if (xlat>30.0 .and. xlat<36.0 .and. xlong<= -7.5) then
1242  if (r412_tbl > 5.8) then
1243  if (scat_ang >= 100.0 .and. scat_ang < 175.0)
1244  1 r412ss = r412_tbl + 4.0*(scat_ang-100.)/75.
1245  if (scat_ang >= 175.0)
1246  1 r412ss = r412_tbl+4.0+0.0*(scat_ang-175.)/5.
1247  endif
1248  endif
1249 
1250 ! Morocco 2
1251  if (xlat>31.5 .and. xlat<36.0 .and. xlong<= -5.0.and. xlong> -7.5) then
1252  if (r412_tbl > 5.8) then
1253  if (scat_ang >= 100.0 .and. scat_ang < 175.0)
1254  1 r412ss = r412_tbl + 4.0*(scat_ang-100.)/75.
1255  if (scat_ang >= 175.0)
1256  1 r412ss = r412_tbl+4.0+0.0*(scat_ang-175.)/5.
1257  endif
1258  endif
1259 
1260 
1261  if (xday > 151.0 .and. xday < 244.0) go to 37 ! winter, spring
1262 
1263  dda2 = 0.5
1264  if (xday > 243.0 .and. xday < 274.0) dda2 = 0.0 ! Sept
1265 ! Central Algeria
1266  if (xlat >22.5 .and. xlat <30. .and. xlong > -4.0.and. xlong < 1.) then
1267  if (r412_tbl > 10.) then
1268  r412ss = r412ss + dda2
1269  endif
1270  endif
1271 
1272 ! Mali - Algeria border
1273  if (xlat>20. .and. xlat<22.501.and. xlong>-2.0.and. xlong< 2.) then
1274  if (r412_tbl > 10.) then
1275  r412ss = r412ss + dda2-1.0
1276  endif
1277  endif
1278 
1279  if (xday > 243.0 .and. xday < 274.0) then
1280  if (xlat>16. .and. xlat<22.501.and. xlong>-12.5.and. xlong< 2.) then
1281  if (r412_tbl > 12.) then
1282  r412ss = r412ss - 1.0
1283  endif
1284  endif
1285  endif
1286 
1287 37 continue
1288 
1289  if (xday > 59.0 .and. xday < 244.0) then ! summer
1290  if (xphi > 90.0) then
1291  ddx = 0.0
1292  ddx33 = 0.6
1293  if (r412_tbl > 12.0) ddx33 = 0.9
1294  if (xthet < 70.0) ddx = xthet * ddx33/70.
1295  if (xthet >= 70.0) ddx = ddx33
1296  r412ss = r412ss + ddx
1297  endif
1298 
1299  if (xday > 59.0 .and. xday < 152.0) go to 38
1300  if (r412_tbl > 9.6) then
1301  if (xlat >20. .and. xlong > 14.9) go to 38
1302  if (xlat >15. .and. xlong > 22.0) go to 38
1303  endif
1304  endif
1305 38 continue
1306 
1307 ! -- special treatment for gzone 21 to alleviate a hot spot.
1308  if (gzflg == 21) then
1309  if (xphi > 90.0) then
1310  ddx = 0.0
1311  ddx33 = 3.0
1312  if (xthet < 70.0) ddx = xthet * ddx33/70.
1313  if (xthet >= 70.0) ddx = ddx33
1314  r412 = r412 + ddx
1315  endif
1316  endif
1317 
1318  if (lprint > 0) print *,'scat_ang,b,f r412ss =',r412_tbl,scat_ang,r412ss
1319 
1320  if (r412ss2 > 20.0) r412ss2 = 19.9
1321  if (r412ss2 < 1.0) r412ss2 = 1.0
1322  if (r412 > 20.0 .and. r412 < 40.0) r412 = 19.9
1323  if (r412 > 40.0) go to 10
1324  if (r412 < 1.0) r412 = 1.0
1325  if (r412ss >= 20.0 .and. r412ss < 40.0) r412ss = 19.9
1326  if (r412ss < 1.0) r412ss = 1.0
1327 
1328  if (r470 >= 24.0) r470 = 23.9
1329  if (r470 < 1.0) r470 = 1.0
1330  if (r470ss >= 24.0) r470ss = 23.9
1331  if (r470ss < 1.0) r470ss = 1.0
1332 
1333  if (lprint > 0) print *,'final r412,r412ss =',r412,r412ss
1334 ! write (888,'(3(F10.5,1X))') xlat, xlong, r412ss
1335 
1336  if (r650 >= 47.0) r650 = 46.9
1337  if (r650 < 1.0) r650 = 1.0
1338 
1339 c -- END DIRECT COPY SEAWIFS code
1340 
1341  if (r412.gt.12.0.and.r412.lt.80.) qa_flag(4) = 2
1342 
1343 c--------------------------------------------------------
1344 c cloud screening
1345 c--------------------------------------------------------
1346 c-- band26 values should be fill value(-999) so most of these
1347 c-- test are non-functional. confirmed w/ dr. hsu 2013-03-01 --cb.
1348  if (debug) print *,'--- start cloud screening ----'
1349  if (debug) print *, 'band26, realbuf(7), realbuf(11): ', band26, realbuf(7), realbuf(11)
1350  if (band26.lt.0.0.and.realbuf(11).lt.0.0.and.
1351  1 realbuf(7).gt.0.0) go to 620
1352  if (band26.gt.0.0.and.realbuf(11).lt.0.0.and.
1353  1 realbuf(7).gt.0.0) go to 10
1354  rat = ref650/ rr412_mod
1355  rat1 = rr470_mod / rr412_mod
1356  if (debug) print *, 'rat, rat1, ref650: ', rat, rat1, ref650
1357  if (ref650.gt.45.0.and.rat.lt.1.4) go to 10
1358  if (ref650.gt.56.0.and.rat.lt.1.3) go to 10
1359 
1360  if (debug) print *, 'trflg, rr412_mod, rat1, r412: ', trflg, rr412_mod, rat1, r412
1361  if (trflg.gt.0.0.and.rr412_mod.gt.10.0.and.
1362  1 rat1.gt.1.25) go to 50
1363 
1364  if (r412.gt.6.0) then
1365  if (band26.gt.0.0.and.rr412_mod.gt.6.0) go to 10
1366  else
1367  if (band26.gt.0.0.and.rr412_mod.gt.10.0.and.rat1.lt.1.25)
1368  1 go to 10
1369  endif
1370  if (debug) print *, '--- end cloud screening ---'
1371 
1372 50 continue
1373 
1374 c--------------------------------------------------------
1375 c special case: for thick strong-absorbing dust plume,
1376 c go to 3-channels
1377 c--------------------------------------------------------
1378 
1379  if (ref650 > 32.0 .and. (gzflg >= 6 .and. gzflg /= 21 .and. gzflg /= 24
1380  1 .and. gzflg /= 13 .and. gzflg /= 31) .and. dstar1 > 1.1) then
1381  xxrat = 0.8
1382 
1383  abs_aero_flag = .true. ! -- set flag for use below.
1384 
1385  go to 620
1386  endif
1387  80 continue
1388 
1389 ! -- vegetated area in zone 1 next to DMN_Maine_Soroa site in summer. Tweak.
1390  if (gzflg == 1 .AND. (gdt1%month >= 6 .AND. gdt1%month <= 8)) then
1391  if (r650_135 < 16.0) then
1392  if (toa_ndvi < 0.18) then
1393  r412 = min(r412 + 2.0, 20.0)
1394  r470new = min(r470new + 2.0, 24.0)
1395  r650 = min(r650 + 2.0, 47.0)
1396  else
1397  r412 = max(r412 - 2.0, 1.0)
1398  r470new = max(r470new - 2.0, 1.0)
1399  r650 = max(r650 - 2.0, 1.0)
1400  end if
1401  end if
1402  end if
1403 
1404 c--------------------------------------------------------
1405 c preliminary retrieval on aot
1406 c--------------------------------------------------------
1407 
1408 c
1409 c for moderate aot, Use 412 - 470 nm pair
1410 c
1411 c retrieving 470 nm aot
1412 c
1413  if (debug) print *, '--- starting aot470 --- '
1414  refl = refl3
1415  x3 = xphi
1416 
1417  tau_x470_flag = -999
1418  tau_x470_flag2 = -999
1419 
1420 ! -- retrieve using most aerosol models -- in preparation for call to get_aot500.
1421  tau_x470_new_91 = -999.0 ; tau_x470_new_92 = -999.0 ; tau_x470_new_93 = -999.0
1422  tau_x470_new_94 = -999.0 ; tau_x470_new_95 = -999.0 ; tau_x470_new_96 = -999.0
1423  tau_x470_new_995 = -999.0
1424 
1425  imod = 3 ! w0 = 0.96
1426  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1427  1 imod,r470,tau_x470,tau_x470_flag,trflg,0.0,debug)
1428  if (dflag) go to 10
1429  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag, tau_x470)
1430  if (status /= 0) then
1431  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1432  return
1433  end if
1434 
1435 ! if (gzflg < 0 .or. gzflg > 11) go to 81 ! gzflg > 11 everywhere except N. Africa/SAP.
1436 
1437  imod = 3 ! w0 = 0.96
1438  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1439  1 imod,r470ss,tau_x470ss,tau_x470_flag2,trflg,0.0,debug)
1440  if (dflag) go to 10
1441  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470ss)
1442  if (status /= 0) then
1443  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1444  return
1445  end if
1446 
1447  if (tau_x470.lt.0.0601.and.dstar1.gt.1.1.and.rat1.gt.1.6)
1448  1 then
1449  imod = 3 ! w0 = 0.96
1450  r470ss = r470ss - 1.0
1451  if (r470ss .lt. 1.0) r470ss = 1.0
1452  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1453  1 imod,r470ss,tau_x470ss,tau_x470_flag2,trflg,0.0,debug)
1454  if (dflag) go to 10
1455  go to 313
1456  endif
1457 
1458  if (tau_x470.lt.0.5.and.dstar1.gt.0.98.and.rat1.gt.1.46)
1459  1 then
1460  imod = 3 ! w0 = 0.96
1461  r470ss = r470ss - 1.0
1462  if (dstar1.gt.1.04) r470ss = r470ss - 0.5
1463  if (r470ss .lt. 1.0) r470ss = 1.0
1464  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1465  1 imod,r470ss,tau_x470ss2,tau_x470_flag2,trflg,0.0,debug)
1466  if (dflag) go to 10
1467  endif
1468 313 continue
1469 
1470  if (r470new >= 24.0) go to 81
1471  if (r470new .lt. 1.0) r470new = 1.0
1472 
1473  imod = 3 ! w0 = 0.96
1474  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1475  1 imod,r470new,tau_x470_new,tau_x470_flag2,trflg,0.0,debug)
1476  if (dflag) go to 10
1477  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new)
1478  if (status /= 0) then
1479  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1480  return
1481  end if
1482 
1483  tau_x470_new_96 = tau_x470_new
1484 
1485  imod = 4 ! w0=0.995
1486  call aero_470(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1487  1 imod,r470new,tau_x470_new_995,tau_x470_flag2,trflg,0.0,debug)
1488  if (dflag) go to 10
1489  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_995)
1490  if (status /= 0) then
1491  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1492  return
1493  end if
1494 
1495 
1496  imod = 2
1497  model_frac = 0.5 ! w0=0.95
1498  call aero_470(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1499  1 imod,r470new,tau_x470_new_95,tau_x470_flag2,trflg,model_frac,debug)
1500  if (dflag) go to 10
1501  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_95)
1502  if (status /= 0) then
1503  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1504  return
1505  end if
1506 
1507  imod = 2 ! w0=0.94
1508  call aero_470(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1509  1 imod,r470new,tau_x470_new_94,tau_x470_flag2,trflg,0.0,debug)
1510  if (dflag) go to 10
1511  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_94)
1512  if (status /= 0) then
1513  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1514  return
1515  end if
1516 
1517  imod = 1
1518  model_frac = 2.0/3.0 ! w0=0.93
1519  call aero_470(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1520  1 imod,r470new,tau_x470_new_93,tau_x470_flag2,trflg,model_frac,debug)
1521  if (dflag) go to 10
1522  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_93)
1523  if (status /= 0) then
1524  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1525  return
1526  end if
1527 
1528  imod = 1
1529  model_frac = 1.0/3.0 ! w0=0.92
1530  call aero_470(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1531  1 imod,r470new,tau_x470_new_92,tau_x470_flag2,trflg,model_frac,debug)
1532  if (dflag) go to 10
1533  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_92)
1534  if (status /= 0) then
1535  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1536  return
1537  end if
1538 
1539  imod = 1 ! w0 = 0.91
1540  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1541  1 imod,r470new,tau_x470_new_91,tau_x470_flag2,trflg,0.0,debug)
1542  if (dflag) go to 10
1543  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_91)
1544  if (status /= 0) then
1545  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1546  return
1547  end if
1548 
1549  if (r470 < 12.0 .and. rr470_mod > 11.0) then
1550  rat_470_412 = rr470_mod / rr412_mod
1551  if (rat_470_412 > 1.85) then
1552  imod = 1 ! w0 = 0.91
1553  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1554  1 imod,r470new,tau_x470_new_91,tau_x470_flag2,trflg,0.0,debug)
1555  if (dflag) go to 10
1556  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_91)
1557  if (status /= 0) then
1558  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1559  return
1560  end if
1561  endif
1562  endif
1563 
1564  if (xthet > 62.0) then
1565  rat_470_412 = rr470_mod / rr412_mod
1566  if (rat_470_412 > 1.7) then
1567  imod = 1 ! w0 = 0.91
1568  call aero_470(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1569  1 imod,r470new,tau_x470_new_91,tau_x470_flag2,trflg,0.0,debug)
1570  if (dflag) go to 10
1571  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag2, tau_x470_new_91)
1572  if (status /= 0) then
1573  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1574  return
1575  end if
1576  end if
1577  end if
1578 
1579 !--------------------------------------------------------
1580 ! 470 nm aod using dust table
1581 
1582  tau_x470_flag_dust = -999
1583 
1584 ! -- retrieve using most aerosol models -- in preparation for call to get_aot500.
1585  tau_x470_dust_91 = -999.0 ; tau_x470_dust_92 = -999.0 ; tau_x470_dust_93 = -999.0
1586  tau_x470_dust_94 = -999.0 ; tau_x470_dust_95 = -999.0 ; tau_x470_dust_96 = -999.0
1587  tau_x470_dust_995 = -999.0
1588 
1589  imod = 3 ! w0 = 0.96
1590  call aero_470_dust(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1591  1 imod,r470new,tau_x470_dust_96,tau_x470_flag_dust,trflg,0.0,debug)
1592  if (dflag) go to 10
1593  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_96)
1594  if (status /= 0) then
1595  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1596  return
1597  end if
1598 
1599  imod = 4 ! w0=0.995
1600  call aero_470_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1601  1 imod,r470new,tau_x470_dust_995,tau_x470_flag_dust,trflg,0.0,debug)
1602  if (dflag) go to 10
1603  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_995)
1604  if (status /= 0) then
1605  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1606  return
1607  end if
1608 
1609  imod = 2
1610  model_frac = 0.5 ! w0=0.95
1611  call aero_470_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1612  1 imod,r470new,tau_x470_dust_95,tau_x470_flag_dust,trflg,model_frac,debug)
1613  if (dflag) go to 10
1614  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_95)
1615  if (status /= 0) then
1616  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1617  return
1618  end if
1619 
1620  imod = 2 ! w0=0.94
1621  call aero_470_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1622  1 imod,r470new,tau_x470_dust_94,tau_x470_flag_dust,trflg,0.0,debug)
1623  if (dflag) go to 10
1624  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_94)
1625  if (status /= 0) then
1626  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1627  return
1628  end if
1629 
1630  imod = 1
1631  model_frac = 2.0/3.0 ! w0=0.93
1632  call aero_470_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1633  1 imod,r470new,tau_x470_dust_93,tau_x470_flag_dust,trflg,model_frac,debug)
1634  if (dflag) go to 10
1635  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_93)
1636  if (status /= 0) then
1637  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1638  return
1639  end if
1640 
1641  imod = 1
1642  model_frac = 1.0/3.0 ! w0=0.92
1643  call aero_470_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1644  1 imod,r470new,tau_x470_dust_92,tau_x470_flag_dust,trflg,model_frac,debug)
1645  if (dflag) go to 10
1646  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_92)
1647  if (status /= 0) then
1648  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1649  return
1650  end if
1651 
1652  imod = 1 ! w0 = 0.91
1653  call aero_470_dust(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1654  1 imod,r470new,tau_x470_dust_91,tau_x470_flag_dust,trflg,0.0,debug)
1655  if (dflag) go to 10
1656  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_91)
1657  if (status /= 0) then
1658  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1659  return
1660  end if
1661 
1662  if (r470 < 12.0 .and. rr470_mod > 11.0) then
1663  rat_470_412 = rr470_mod / rr412_mod
1664  if (rat_470_412 > 1.85) then
1665  imod = 1 ! w0 = 0.91
1666  call aero_470_dust(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1667  1 imod,r470new,tau_x470_dust_91,tau_x470_flag_dust,trflg,0.0,debug)
1668  if (dflag) go to 10
1669  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_91)
1670  if (status /= 0) then
1671  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1672  return
1673  end if
1674  endif
1675  endif
1676 
1677  if (xthet > 62.0) then
1678  rat_470_412 = rr470_mod / rr412_mod
1679  if (rat_470_412 > 1.7) then
1680  imod = 1 ! w0 = 0.91
1681  call aero_470_dust(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1682  1 imod,r470new,tau_x470_dust_91,tau_x470_flag_dust,trflg,0.0,debug)
1683  if (dflag) go to 10
1684  status = handle_lut_out_of_bounds(gzflg, tau_x470_flag_dust, tau_x470_dust_91)
1685  if (status /= 0) then
1686  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1687  return
1688  end if
1689  end if
1690  end if
1691 ! end 470 nm aod using dust table
1692 
1693  if (debug) print *, '--- end aot470 ---'
1694 81 continue
1695 
1696 c if (lprint > 0) print *,'tau_x470,tau_x470ss,tau_x470_new=',
1697 c 1 tau_x470,tau_x470ss,tau_x470_new
1698 
1699 
1700 c--------------------------------------------------------
1701 c retrieving 412 nm aot
1702 c--------------------------------------------------------
1703 
1704  refl = refl1
1705  x3 = xphi
1706 
1707  tau_x412_flag = -999
1708  tau_x412_flag2 = -999
1709 
1710  r412new = r412 ! transitional zone
1711 
1712  imod = 5 ! w0 = 0.94
1713  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1714  1 imod,r412,tau_x412,tau_x412_flag2,trflg,0.0,debug)
1715  if (dflag) go to 10
1716  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag2, tau_x412)
1717  if (status /= 0) then
1718  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1719  return
1720  end if
1721 !---
1722  if (r412 > 11.0 .and. tau_x412 < 0.4) then
1723  tau_x412_91 = tau_x412 * 2.
1724  go to 630
1725  else if (r412 > 12.0) then
1726  tau_x412_91 = tau_x412 * 2.
1727  go to 630
1728  endif
1729 
1730  refl = refl1
1731  x3 = xphi
1732 
1733  tau_x412_flag_91 = -999
1734 
1735  imod = 4 ! w0 = 0.91
1736  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1737  1 imod,r412,tau_x412_91,tau_x412_flag_91,trflg,0.0,debug)
1738  if (dflag) go to 10
1739  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag_91, tau_x412_91)
1740  if (status /= 0) then
1741  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1742  return
1743  end if
1744 
1745 630 continue
1746 
1747  tau_x412_flag = -999
1748  tau_x412_new_91 = -999.0 ; tau_x412_new_93 = -999.0 ; tau_x412_new_94 = -999.0
1749  tau_x412_new_96 = -999.0 ; tau_x412_new_995 = -999.0
1750 
1751  tau_x412ss_995 = -999.0 ; tau_x412ss2_995 = -999.0 ; tau_x412ss_96 = -999.0
1752  tau_x412ss_97 = -999.0 ; tau_x412ss_94 = -999.0 ; tau_x412ss_95 = -999.0
1753  tau_x412ss_98 = -999.0
1754 
1755  if (r412new < 1.0) go to 631
1756  if (r412new >= 20.0) go to 631
1757 
1758  imod = 8 ! w0=0.995
1759  call aero_412(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1760  1 imod,r412new,tau_x412_new_995,tau_x412_flag,trflg,0.0,debug)
1761  if (dflag) go to 10
1762  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412_new_995)
1763  if (status /= 0) then
1764  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1765  return
1766  end if
1767 ! print *
1768 ! print *, 'aero_412, tau_x412_new_995: ', refl, sza, xthet, xphi, imod, r412new, tau_x412_new_995, tau_x412_flag, trflg, dflag
1769 ! print *
1770 
1771  imod = 6 ! w0=0.96
1772  call aero_412(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1773  1 imod,r412new,tau_x412_new_96,tau_x412_flag,trflg,0.0,debug)
1774  if (dflag) go to 10
1775  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412_new_96)
1776  if (status /= 0) then
1777  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1778  return
1779  end if
1780 
1781  imod = 5 ! w0=0.94
1782  call aero_412(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1783  1 imod,r412new,tau_x412_new_94,tau_x412_flag,trflg,0.0,debug)
1784  if (dflag) go to 10
1785  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412_new_94)
1786  if (status /= 0) then
1787  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1788  return
1789  end if
1790 
1791  imod = 4 ! w0=0.93
1792  model_frac = 2.0/3.0
1793  if (dflag) print *,'calling aero_412, model_frac: imod, model_frac: ', imod, model_frac
1794  call aero_412(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1795  1 imod,r412new,tau_x412_new_93,tau_x412_flag,trflg,model_frac,debug)
1796  if (dflag) go to 10
1797  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412_new_93)
1798  if (status /= 0) then
1799  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1800  return
1801  end if
1802 
1803  imod = 4 ! w0=0.91
1804  call aero_412(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1805  1 imod,r412new,tau_x412_new_91,tau_x412_flag,trflg,0.0,debug)
1806  if (dflag) go to 10
1807  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412_new_91)
1808  if (status /= 0) then
1809  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1810  return
1811  end if
1812 
1813 ! -- arabian peninsula
1814  if (gzflg >= 6 .and. gzflg <= 11) then
1815  tau_x412_new = tau_x412_new_94
1816  end if
1817 
1818  if (gzflg == 6) then
1819  if (xday >= 60.0 .and. xday < 274.0) go to 602
1820  tau_x470_new = tau_x412_new
1821  end if
1822 
1823 602 continue
1824 
1825  if (gzflg == 8) then
1826  if (xday >= 182.0 .and. xday < 274.0) go to 603
1827  tau_x470_new = tau_x412_new
1828  end if
1829 603 continue
1830 
1831  if (gzflg == 7) tau_x470_new = tau_x412_new
1832 
1833  if (gzflg == 9) tau_x470_new = tau_x412_new
1834 
1835  if (gzflg == 10) then
1836  if (xday >= 60.0 .and. xday < 274.0) go to 605
1837  tau_x470_new = tau_x412_new
1838  end if
1839 605 continue
1840 
1841 !---------------------------------------------------------------
1842 ! 412 nm aod using dust table
1843 
1844  tau_x412_flag_dust = -999
1845  tau_x412_dust_91 = -999.0 ; tau_x412_dust_93 = -999.0 ; tau_x412_dust_94 = -999.0
1846  tau_x412_dust_96 = -999.0 ; tau_x412_dust_995 = -999.0
1847 
1848  imod = 8 ! w0=0.995
1849  call aero_412_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1850  1 imod,r412new,tau_x412_dust_995,tau_x412_flag_dust,trflg,0.0,debug)
1851  if (dflag) go to 10
1852  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag_dust, tau_x412_dust_995)
1853  if (status /= 0) then
1854  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1855  return
1856  end if
1857 
1858  imod = 6 ! w0=0.96
1859  call aero_412_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1860  1 imod,r412new,tau_x412_dust_96,tau_x412_flag_dust,trflg,0.0,debug)
1861  if (dflag) go to 10
1862  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag_dust, tau_x412_dust_96)
1863  if (status /= 0) then
1864  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1865  return
1866  end if
1867 
1868  imod = 5 ! w0=0.94
1869  call aero_412_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1870  1 imod,r412new,tau_x412_dust_94,tau_x412_flag_dust,trflg,0.0,debug)
1871  if (dflag) go to 10
1872  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag_dust, tau_x412_dust_94)
1873  if (status /= 0) then
1874  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1875  return
1876  end if
1877 
1878  imod = 4 ! w0=0.93
1879  model_frac = 2.0/3.0
1880  call aero_412_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1881  1 imod,r412new,tau_x412_dust_93,tau_x412_flag_dust,trflg,model_frac,debug)
1882  if (dflag) go to 10
1883  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag_dust, tau_x412_dust_93)
1884  if (status /= 0) then
1885  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1886  return
1887  end if
1888 
1889  imod = 4 ! w0=0.91
1890  call aero_412_dust(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1891  1 imod,r412new,tau_x412_dust_91,tau_x412_flag_dust,trflg,0.0,debug)
1892  if (dflag) go to 10
1893  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag_dust, tau_x412_dust_91)
1894  if (status /= 0) then
1895  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1896  return
1897  end if
1898 ! end 412 nm aod using dust table
1899 
1900 631 continue
1901 
1902  if ((gzflg <= 5 .and. gzflg > 0) .OR. gzflg == 27) then
1903 
1904  tau_x412_flag = -999
1905 ! if (r412ss > 19.89) print *, 'out, r412ss: ', r412ss,xlat,xlong
1906 
1907  imod = 4 ! w0 = 0.91
1908  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1909  1 imod,r412ss,tau_x412ss_91,tau_x412_flag,trflg,0.0,debug)
1910  if (dflag) go to 10
1911  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss_91)
1912  if (status /= 0) then
1913  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1914  return
1915  end if
1916 
1917  imod = 5 ! w0 = 0.94
1918  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1919  1 imod,r412ss,tau_x412ss,tau_x412_flag,trflg,0.0,debug)
1920  if (dflag) go to 10
1921  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss)
1922  if (status /= 0) then
1923  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1924  return
1925  end if
1926 
1927  tau_x412ss_94 = tau_x412ss
1928  imod = 8 ! w0=0.995
1929  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1930  1 imod,r412ss,tau_x412ss_995,tau_x412_flag,trflg,0.0,debug)
1931  if (dflag) go to 10
1932  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss_995)
1933  if (status /= 0) then
1934  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1935  return
1936  end if
1937 
1938  imod = 7 ! w0=0.98
1939  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1940  1 imod,r412ss,tau_x412ss_98,tau_x412_flag,trflg,0.0,debug)
1941  if (dflag) go to 10
1942  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss_98)
1943  if (status /= 0) then
1944  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1945  return
1946  end if
1947 
1948 
1949  imod = 6 ! w0=0.96
1950  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
1951  1 imod,r412ss,tau_x412ss_96,tau_x412_flag,trflg,0.0,debug)
1952  if (dflag) go to 10
1953  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss_96)
1954  if (status /= 0) then
1955  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1956  return
1957  end if
1958 
1959  imod = 5 ! w0=0.95
1960  model_frac = 0.5
1961  call aero_412(dflag,refl,sza,xthet,xphi,mm,nn,ll,ma,
1962  1 imod,r412ss,tau_x412ss_95,tau_x412_flag,trflg,model_frac,debug)
1963  if (dflag) go to 10
1964  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss_95)
1965  if (status /= 0) then
1966  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
1967  return
1968  end if
1969 
1970 ! General surfaces, tau_x412ss = default (0.94)
1971  if (xday < 60.0 .or. xday > 334.0) then ! Dec, Jan, Feb
1972  if (tau_x412ss_97 <0.6) tau_x412ss = tau_x412ss_96
1973  endif
1974  if (xday > 243.0 .and. xday < 274.0) then ! Sep
1975  if (tau_x412ss_97 <0.6) tau_x412ss = tau_x412ss_96
1976  endif
1977  if (xday > 273.0 .and. xday < 335.0) then ! Oct, Nov
1978  if (tau_x412ss < 0.5) tau_x412ss = tau_x412ss_96
1979  endif
1980  if (xday > 59.0 .and. xday < 121.0) then ! Mar,Apr
1981  if (tau_x412ss < 0.5) tau_x412ss = tau_x412ss_96
1982  endif
1983 
1984 ! Bright Surfaces
1985  ddx = 10.5
1986  if (xday > 151.0 .and. xday < 244.0) ddx = 12.0 ! Jun, Jul, Aug
1987  if (xday > 243.0 .and. xday < 274.0) ddx = 12.0 ! Sep
1988  if (r412_tbl > ddx ) then
1989  if (tau_x412ss_97 < 0.6) tau_x412ss = tau_x412ss_96
1990  if (xday > 243.0 .or. xday < 60.0) then ! Sep,Nov,Dec,Jan,Feb
1991  if (tau_x412ss_995 < 0.5) tau_x412ss = tau_x412ss_995
1992  if (xday > 273.0 .and. xday < 305.0) then ! Oct
1993  if (tau_x412ss_97 < 0.6) tau_x412ss = tau_x412ss_96
1994  endif
1995  if (xlat >10.0 .and. xlat <21.0 .and. xlong > 10.0 .and. xlong < 20.0)
1996  1 tau_x412ss = tau_x412ss_96
1997  endif
1998  if (xday > 59.0 .and. xday < 121.0) then ! Mar,Apr
1999  tau_x412ss = tau_x412ss_96
2000  if (xlat >10.0 .and. xlat <21.0 .and. xlong > 10.0 .and. xlong < 20.0)
2001  1 tau_x412ss = tau_x412ss_96
2002  if (dstar1 < 1.01 .and. tau_x412ss>0.6) tau_x412ss = tau_x412ss_995
2003  if (dstar1 > 1.01 .and. dstar1 < 1.04 .and. tau_x412ss>0.6) tau_x412ss = tau_x412ss_98
2004  endif
2005  endif
2006 
2007 ! Deserts near Egypt and Sudan
2008  if (r412_tbl > 10.5 .and. dstar1 > 1.04) then
2009  if (xlat >20.0 .and. xlat <25.0 .and. xlong > 25.0 .and. xlong < 30.0)
2010  1 tau_x412ss = tau_x412ss_94
2011  endif
2012 
2013 ! Deserts near Libya and Egypt
2014  ddx = 1.01
2015  if (xday > 151.0 .and. xday < 305.0) ddx = 1.04
2016 
2017  dda1 = tau_x412ss
2018 
2019  if (xlat >20. .and. xlong > 12.9) then
2020  if (r412_tbl > 8.5) then
2021  if (xlat <22. .and. xlong < 17.0) go to 632
2022  if (dstar1 .lt. ddx) then
2023  dd = (xlong - 11.9) /3.0
2024  if (xlong > 14.9) dd = 1.
2025  tau_x412ss = tau_x412ss-(tau_x412ss-tau_x412ss_995)*dd
2026  dda2 = tau_x412ss
2027  dd = (xlat - 18.0) /4.0
2028  if (xlat > 22.0) dd = 1.
2029  tau_x412ss = dda1-(dda1-dda2)*dd
2030  endif
2031  endif
2032  endif
2033  if (xlat >15. .and. xlat <=20. .and. xlong > 22.0) then
2034  if (r412_tbl > 8.5) then
2035  if (dstar1 .lt. ddx) then
2036  dda2 = tau_x412ss_995
2037  dd = (xlat - 18.0) /4.0
2038  if (xlat > 22.0) dd = 1.
2039  if (xlat < 18.0) dd = 0.
2040  tau_x412ss = dda1-(dda1-dda2)*dd
2041  endif
2042  endif
2043  endif
2044  if (xlat >19.0 .and. xlat <=20.0 .and. xlong > 19.8) then
2045  if (r412_tbl > 8.5) then
2046  if (dstar1 .lt. ddx) then
2047  dda2 = tau_x412ss_995
2048  dd = (xlat - 18.0) /4.0
2049  if (xlat > 22.0) dd = 1.
2050  tau_x412ss = dda1-(dda1-dda2)*dd
2051  endif
2052  endif
2053  endif
2054 
2055  if (xday > 243.0 .and. xday < 274.0) then
2056  if (xlat >19.0 .and. xlat <=21.0 .and. xlong > 19.8.and. xlong <23.5) then
2057  if (r412_tbl > 8.5) then
2058  if (dstar1 .lt. ddx) tau_x412ss = tau_x412ss_995
2059  endif
2060  endif
2061  endif
2062 
2063 ! NE Mauritania 3
2064  if (xday > 243.0 .and. xday < 274.0) then
2065  if (xlat>22.5 .and. xlat<30.0.and. xlong>-11.0.and. xlong< -5.) then
2066  if (r412_tbl > 9.0) tau_x412ss = tau_x412ss_94
2067  endif
2068  if (xlat>22.5 .and. xlat<26.0.and. xlong>-13.0.and. xlong< -11.001) then
2069  if (r412_tbl > 9.0) tau_x412ss = tau_x412ss_94
2070  endif
2071  if (xlat>=20.0 .and. xlat<29.0.and. xlong< -12.5) then
2072  if (r412_tbl > 10.5) tau_x412ss = tau_x412ss_94
2073  endif
2074  if (xlat>15.0 .and. xlat<20.0 .and. xlong< -14.9) then
2075  if (r412_tbl > 10.5) tau_x412ss = tau_x412ss_94
2076  endif
2077  endif
2078 
2079  if (xday > 151.0 .and. xday < 274.0) go to 632 ! Jun, Jul, Aug
2080 
2081 ! NE Mauritania 1
2082  if (xlat>20.0 .and. xlat<30.0 .and. xlong< -12.5) then
2083  if (r412_tbl > 9.) then
2084  if (dstar1 .lt. 1.01) tau_x412ss = tau_x412ss_995
2085  endif
2086  endif
2087 
2088 ! NE Mauritania 2
2089  if (xlat>26.0 .and. xlat<30.0.and. xlong>=-12.5.and. xlong< -11.001) then
2090  if (r412_tbl > 9.) then
2091  if (dstar1 .lt. 1.01) tau_x412ss = tau_x412ss_995
2092  endif
2093  endif
2094 
2095 ! Morocco 1
2096  if (xlat>30.0 .and. xlat<36.0 .and. xlong<= -7.5) then
2097  if (r412_tbl > 5.8) then
2098  if (dstar1 .lt. 1.01) tau_x412ss = tau_x412ss_995
2099  endif
2100  endif
2101 
2102 ! Morocco 2
2103  if (xlat>31.5 .and. xlat<36.0 .and. xlong<= -5.0.and. xlong> -7.5) then
2104  if (r412_tbl > 5.8) then
2105  if (dstar1 .lt. 1.01) tau_x412ss = tau_x412ss_995
2106  endif
2107  endif
2108 
2109 632 continue
2110 
2111  go to 635
2112 633 continue
2113  if (tau_x412ss_97 < 0.6) tau_x412ss = tau_x412ss_96
2114  if (xday > 181.0 .and. xday < 274.0) then ! Jul, Aug, Sep
2115  if (tau_x412ss_995 <0.5) tau_x412ss = tau_x412ss_995
2116  endif
2117  if (xday > 273.0 .and. xday < 335.0) then ! Oct, Nov
2118  if (tau_x412ss_995 <0.5) tau_x412ss = tau_x412ss_995
2119  endif
2120  if (xday > 59.0 .and. xday < 121.0) then ! Mar,Apr
2121  if (tau_x412ss_995 <0.5) tau_x412ss = tau_x412ss_995
2122  endif
2123 
2124 635 continue
2125  endif
2126 
2127  if (gzflg >= 6 .and. gzflg <= 11) then
2128 
2129  tau_x412_flag = -999
2130  if (r412ss2 > 20.0) print *, 'out, r412ss2: ', r412ss2
2131  if (r412ss2 < 1.0) print *, 'less, r412ss2: ', r412ss2
2132 
2133  imod = 5 ! w0 = 0.94
2134  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
2135  1 imod,r412ss2,tau_x412ss2,tau_x412_flag,trflg,0.0,debug)
2136  if (dflag) go to 10
2137  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss2)
2138  if (status /= 0) then
2139  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
2140  return
2141  end if
2142 
2143  imod = 7 ! w0=0.98
2144  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
2145  1 imod,r412ss2,tau_x412ss2_98,tau_x412_flag,trflg,0.0,debug)
2146  if (dflag) go to 10
2147  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss2_98)
2148  if (status /= 0) then
2149  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
2150  return
2151  end if
2152 
2153  imod = 8 ! w0=0.995
2154  call aero_412(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
2155  1 imod,r412ss2,tau_x412ss2_995,tau_x412_flag,trflg,0.0,debug)
2156  if (dflag) go to 10
2157  status = handle_lut_out_of_bounds(gzflg, tau_x412_flag, tau_x412ss2_98)
2158  if (status /= 0) then
2159  print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
2160  return
2161  end if
2162  end if
2163 
2164 c If surface is vegetated, derive the surface reflectivity at 490 nm and 670 nm
2165 c from 865 nm surface ler and ndvi. specifically exclude sahara and arabia via gzflg
2166 c flag. exclude areas with aeronet-based brdf info via brdf_flag == 0.
2167 c--------------------------------------------------------------------------------------------------
2168  637 continue
2169 
2170 ! -- set up our NDVI threshold over which the veg code can run.
2171 ! -- if we failed to get surface info above (sr_fail_flag == .true.),
2172 ! -- let veg try no matter the NDVI (unless superseded by another test
2173 ! -- later on.
2174  ndvi_thold = 0.1
2175  if (sr_fail_flag .eqv. .true.) then
2176  ndvi_thold = -999.0
2177  end if
2178 
2179  if (lc == 6) then
2180  ndvi_thold = 0.2
2181  end if
2182 
2183 ! -- increase NDVI threshold over high elevation or barren areas
2184 ! -- in the United States (13)
2185  if (gzflg == 13) then
2186  if (px_elev >= 500.0 .OR. lc == 6) then
2187  ndvi_thold = 0.3
2188  end if
2189  endif
2190 
2191 ! -- increase NDVI threshold over high elevation areas in North and South America.
2192 ! -- this covers Central America which does not currently have a region assigned to it.
2193  if (gzflg < -900 .AND. (xlat < 30.0 .AND. xlong < -30.0)) then
2194  if (px_elev >= 500) then
2195  ndvi_thold = 0.3
2196  end if
2197  end if
2198 
2199 ! -- decide if we need to run the veg retrieval.
2200  do_veg = .false.
2201  if (toa_ndvi >= ndvi_thold .AND. brdf_flag /= 0) then
2202  do_veg = .true.
2203  end if
2204 
2205 ! -- override do_veg from above to exclude vegetated algorithm from the following
2206 ! -- zones.
2207  if (gzflg > 0 .AND. gzflg <= 11) do_veg = .false. ! N. Africa, Arabia
2208  if (gzflg == 15) do_veg = .false. ! India, Kanpur region
2209  if (gzflg == 20) do_veg = .false. ! India, Thar desert
2210  if (gzflg == 19) do_veg = .false. ! India, Pune region
2211  if (gzflg == 34) do_veg = .false. ! India, high elevation
2212  if (gzflg == 21) do_veg = .false. ! Asia dust transport region
2213  if (gzflg == 23) do_veg = .false. ! Iraq, Iran region
2214  if (gzflg == 24) do_veg = .false. ! Taklimakan desert
2215  if (gzflg == 27) do_veg = .false. ! Sahel transition region
2216 
2217 ! -- some special logic over China in Winter.
2218  if (gzflg == 16 .AND. (xday >= 335 .OR. xday < 60)) then
2219  if (toa_ndvi > 0.2) then
2220  if (lc == 1 .OR. lc == 2 .OR. lc == 4) then
2221  do_veg = .true.
2222  end if
2223  if (lc == 3 .AND. xlat < 25.0) then
2224  do_veg = .true.
2225  end if
2226  end if
2227  end if
2228 
2229 ! -- some special logic over India, only test NDVI regardless of
2230 ! brdf_flag, adjustment of ndvi_thold above won't trigger veg
2231 ! retrievals due to the brdf_flag condition, 3 January 2018 JLee
2232  if ((gzflg == 15 .or. gzflg == 19) .and. toa_ndvi > 0.4) then
2233  do_veg = .true.
2234  end if
2235 
2236 ! -- special logic over N. America urban areas
2237  if (gzflg == 13 .and. toa_ndvi > 0.4) then
2238  do_veg = .true.
2239  end if
2240 
2241 !! -- special logic over Sahel region
2242 ! if ((gzflg == 26 .or. gzflg == 27) .and. toa_ndvi > 0.4) then
2243 ! do_veg = .true.
2244 ! end if
2245 
2246  if (debug) print *, "ndvi, gzflg, brdf_flag, sr_fail_flag, do_veg: ", toa_ndvi, gzflg, brdf_flag, sr_fail_flag, do_veg
2247 
2248 c -- skip pixel if no surf. reflc. value and not suitable for vege. retrieval.
2249  if ((do_veg .eqv. .false.) .AND. (sr_fail_flag .eqv. .true.)) go to 10
2250 
2251  if (do_veg) then ! .OR. oob_samer) then
2252  alg_flag = 1
2253 c if (ndvi670 >= 0.1 .AND. (gzflg < -900 .OR. gzflg > 11) .AND. brdf_flag /= 0) then
2254 c print *, 'doing veg retrieval: ', itmp, jtmp, alg_flag
2255 c get current season and tweak for swf_aero_veg input.
2256 c iopss = season_from_doy(yr, doy)
2257 c iopss = iopss - 1
2258 c if (iopss == 0) iopss = 4
2259 
2260 c... *******************************************************************
2261 c... *******************************************************************
2262 c... *******************************************************************
2263 c-----------------------------------------------------------------------
2264 c... do retrieval over vegetated surfaces(ndvi=>0.1)
2265 c-----------------------------------------------------------------------
2266 c... imod=2 ! ssa_490=0.94 --> You may open this line and remove
2267 c... do loop "do 2500 imod=1,4", or use
2268 c... --> "call swf_aero_veg(nvalx470,nvalx650,iopss,2,sza,xthet,"
2269 c... ------------------------------------------------------
2270 c gdt1 = gregorian_from_doy(yr,doy)
2271 c ioprg=0 ! region index initialization @@new@@
2272 c idx=int((xlong-(-180.0))/0.10)+1 ! @@new@@
2273 c idy=int((xlat-(-90.0))/0.10)+1 ! @@new@@
2274 c if(idx.ge.1.and.idx.le.3600.and.idy.ge.1.and.idy.le.1800) then ! @@new@@
2275 c ioprg=veg_regions(idx,idy) ! @@new@@
2276 c else ! @@new@@
2277 c ioprg=0 ! @@new@@
2278 c endif ! @@new@@
2279 c imod = 2
2280 c select case (ioprg)
2281 c case (6) ! S. Africa
2282 c select case (gdt1%month)
2283 c case (6:11) ! Jun-Nov use more absorbing model (0.89).
2284 c imod = 1
2285 c case default
2286 c imod = 2
2287 c end select
2288 c case default
2289 c imod = 2
2290 c end select
2291 
2292 c tau_x470_flag = -999
2293 c tau_x650_flag = -999
2294 
2295  tau_x470_flag_veg = -999
2296  call find_v_veg(gdt1%month,season,realbuf,tmpvg,
2297  1 r412sv_veg,r470sv_veg,gzflg,outbufvg,tau_x470_flag_veg)
2298  if (outbufvg(7) < 0.0 .and. tau_x470_flag_veg == 1) go to 805
2299 
2300 c translate outbufvg back to local variables.
2301  do i=1,3
2302  xtau(i) = outbufvg(i)
2303  ssa(i) = outbufvg(i+3)
2304  outbuf(i+3) = ssa(i)
2305  enddo
2306  tau550 = outbufvg(7)
2307  alpha = outbufvg(8)
2308  r412 = -999.0
2309  r470 = outbufvg(11)
2310  r650 = outbufvg(12)
2311  tau_x650_flag = outbufvg(10)
2312  xthet = outbufvg(13)
2313  scat_ang = outbufvg(14)
2314  sfc_typ = outbufvg(15)
2315 
2316  if (debug) print *, "veg final: aot550, aot: ", tau550, xtau(1),
2317  & xtau(2), xtau(3)
2318 
2319  ! jlee test 20170726, use 2.2 um surface reflectance for this
2320  ! condition, Dr. Hsu found a high bias over bright vegetated
2321  ! surfaces in CONUS, Africa, etc.
2322  if ((outbufvg(3).gt.0.2.and.alpha.lt.-0.5).or.
2323  & (outbufvg(3).gt.0.2.and.(realbuf(23)/realbuf(6)).lt.0.7).or.
2324  & (outbufvg(2).gt.0.4.and.(realbuf(23)/realbuf(6)).lt.0.7))
2325  & tau550 = tau_x470sv96
2326 
2327 c Return fill for 412 nm, no dark target retrieval yet.
2328 c xtau(1) = -999.0
2329 c tau_x412_flag2 = -999
2330 c if (lprint > 0) print *, 'after veg: ', i, j, xtau(2), xtau(4), tau550, &
2331 c & alpha, tau_x470_flag, tau_x650_flag
2332 c
2333 c! -- reset AOT to zone-specific value if we ran off the LUT.
2334 c status = handle_lut_out_of_bounds(gzflg, tau_x470_flag, xtau(2))
2335 c if (status /= 0) then
2336 c print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
2337 c return
2338 c end if
2339 c status = handle_lut_out_of_bounds(gzflg, tau_x650_flag, xtau(4))
2340 c if (status /= 0) then
2341 c print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
2342 c return
2343 c end if
2344 c status = handle_lut_out_of_bounds(gzflg, tau_x470_flag, tau550)
2345 c if (status /= 0) then
2346 c print *, "ERROR: Failed to check/reset AOT out of bounds condition: ", status
2347 c return
2348 c end if
2349 
2350 
2351 c write (500,'(7(F10.5,1X))') xlat, xlong, r470, r650, xtau(2), xtau(4), tau550
2352 
2353 c if (lprint > 0) print *, 'after LUT check, veg: ', i, j, xlat, xlong, ndvi670, &
2354 c & r470, r650, xtau(2), xtau(4), tau550, alpha
2355 
2356 c! Check consistency
2357  if (xtau(2) < 0.0 .or. xtau(3) < 0.0 .or. tau550 < 0.0) then
2358  xtau(1) = -999.0
2359  xtau(2) = -999.0
2360  xtau(3) = -999.0
2361  ssa(1) = -999.0
2362  ssa(2) = -999.0
2363  ssa(3) = -999.0
2364  tau550 = -999.0
2365  alpha = -999.0
2366  r412 = -999.0
2367  r470 = -999.0
2368  r650 = -999.0
2369  tau_x470_flag = -999
2370  tau_x650_flag = -999
2371  alg_flag = -999
2372  outbufvg(:) = -999.0
2373  end if
2374 
2375  if (debug) print *, 'final: ', xlat, xlong, xtau(2), xtau(3), tau550, alpha
2376 c write(888,'(7(F12.5), I3)') xlat, xlong, xtau(2), xtau(4), r470, r650, z_ndvi, proc_flag(i,j)
2377 c-----------------------------------------------------------------------
2378 c... end of retrieval over vegetated surfaces
2379 c-----------------------------------------------------------------------
2380 c... *******************************************************************
2381 c... *******************************************************************
2382 c... *******************************************************************
2383 
2384 c r470 = get_ssr_490(yr, doy, ler650(i,j)*100.0, ler865(i,j)*100.0, status)
2385 c if (status < 0) then
2386 c print *, "ERROR: Unable to derive vegetated surface "//
2387 c 1 "reflectivity for 490 nm: ", status
2388 c stop
2389 c cycle
2390 c end if
2391 c
2392 c r650 = get_ssr_670(yr, doy, ler650(i,j)*100.0, ler865(i,j)*100.0, status)
2393 c if (status < 0) then
2394 c print *, "ERROR: Unable to derive vegetated surface "//
2395 c 1 "reflectivity for 670 nm: ", status
2396 c stop
2397 c cycle
2398 c end if
2399  go to 865
2400  end if
2401 
2402  alg_flag = 0
2403 
2404 c if (lprint > 0) print *,'tau_x412,tau_x412_new,tau_x412ss,tau_x412ss2=',
2405 c 1 tau_x412,tau_x412_new,tau_x412ss,tau_x412ss2
2406 c--------------------------------------------------------
2407 c retrieving 412 nm ssa(412 - 470 nm)
2408 c--------------------------------------------------------
2409  if (dstar1.gt.1.05.and.rat1.gt.1.6) go to 620
2410  if (tau_x470.lt.0.2.and.tau_x470.gt.0.0) go to 805
2411  if (rr412_mod.gt.20.0.and.tau_x470.lt.0.0) go to 620
2412  if (tau_x470.lt.0.0) go to 620
2413 
2414  refl = refl1
2415  x3 = xphi
2416  tau_x = tau_x470
2417  w0_x = -999.
2418 
2419  if (tau_x >= 3.5) tau_x = 3.5-0.0001 ! search2 dies if tau_x >= tau(10), 3.5
2420  call aero_412_abs(dflag,refl,x1,x2,x3,mm,nn,ll,
2421  1 r412,tau_x,w0_x)
2422  if (dflag) go to 10
2423 
2424  w0_int = w0_x
2425 
2426  if (w0_x.lt.0.0) go to 10
2427 
2428  w0_int_470 = w0_x +(0.976 -w0_x)*(470.-412.)/(650.-412.)
2429  if (w0_int_470.lt.0.0) w0_int_470 = -999.
2430 
2431 c--------------------------------------------------------
2432 c retrieving 650 nm aot
2433 c--------------------------------------------------------
2434  if (dstar1 > 1.1 .and. xlong < -30.0) go to 620
2435  if (xthet.gt.60.0.and.xphi.lt.90.0.and.tau_x470.gt.0.5)
2436  1 go to 620
2437  if (xlong < -30.0) then ! mostly N.America
2438  if (r650.lt.30.0.and.tau_x470.lt.0.7.and.tau_x470.gt.0.0)
2439  1 go to 805
2440  else
2441  if (r650.lt.30.0.and.r650.gt.15.0.and.tau_x470.lt.0.7.and.tau_x470.gt.0.0)
2442  1 go to 805
2443  end if
2444 c if (r650.lt.30.0.and.tau_x470.lt.0.7.and.tau_x470.gt.0.0)
2445 c 1 go to 805
2446  if (r650.ge.30.0.and.tau_x470.lt.1.0.and.tau_x470.gt.0.0)
2447  1 go to 805
2448  if (scat_ang.gt.165.0.and.tau_x470.lt.1.8) go to 805
2449 
2450 620 continue
2451  refl = refl6
2452  x3 = xphi
2453 
2454  tau_x650_flag = -999
2455 
2456  call aero_650(dflag,refl,x1,x2,x3,mm,nn,ll,ma,
2457  1 r650,tau_x650,tau_x650_flag,tau_x470_flag2,
2458  1 tau_x412, tau_x470,tau_x412_flag_91,trflg)
2459  if (dflag) go to 10
2460  status = handle_lut_out_of_bounds(gzflg, tau_x650_flag, tau_x650)
2461  if (status /= 0) then
2462  print *, "ERROR: Failed to check/reset AOT out of bounds
2463  1 condition: ", status
2464  return
2465  end if
2466 
2467 !--------------------------------------------------------
2468 ! Retrieving 412 and 470 nm SSA
2469 !--------------------------------------------------------
2470  if (tau_x650 < 0.0) go to 805
2471 c if (tau_x650 > 3.5) go to 805
2472 
2473 !--- 412 nm
2474 
2475  refl = refl1
2476  x3 = xphi
2477  tau_x = tau_x650
2478  if (rat1.lt.1.8.and.tau_x470ss.gt.0.8) tau_x = tau_x470ss
2479  w0_x = -999.
2480  if (tau_x > 3.5) tau_x = 3.5
2481  if (tau_x > 3.5) then
2482  goto 805
2483  endif
2484  if (tau_x >= 3.5) tau_x = 3.5-0.0001
2485 
2486  call aero_412_abs(dflag,refl,x1,x2,x3,mm,nn,ll,
2487  1 r412,tau_x,w0_x)
2488 
2489  if (dflag) go to 10
2490  w0_x412 = w0_x
2491  if (w0_x < 0.0) go to 10
2492 
2493 !--- 470 nm
2494 
2495  refl = refl3
2496  x3 = xphi
2497  tau_x = tau_x650
2498  w0_x = -999.
2499  if (tau_x >= 3.5) tau_x = 3.5-0.0001
2500 ! ADDED BY COREY
2501  if (r470 > 24.0) go to 805
2502 ! END ADDED BY COREY
2503  call aero_470_abs(dflag2,refl,x1,x2,x3,mm,nn,ll,
2504  1 r470,tau_x,w0_x470)
2505  if (dflag2) go to 10
2506 
2507  if (w0_x470 < 0.0) w0_x470 = -999.
2508 
2509  805 continue
2510 !
2511 ! -- Selecting Models
2512 !
2513 
2514  call aero_mod (tau_x412,tau_x470,tau_x650,tau_x412_91,aot_mod)
2515 
2516  if (xlong < -30.0 .and. xlat >= 10.0 .and. dstar1 > 1.1 .and.
2517  1 realbuf(23) > 12.0 ) then
2518  if ((gzflg == 13 .or. gzflg == 31) .and. realbuf(22)/realbuf(6) < 0.3) go to 807
2519  tau_x412 = tau_x650
2520  tau_x470 = tau_x650
2521  aot = tau_x650
2522  go to 860
2523  endif
2524  807 continue
2525 
2526  if (tau_x412 > 0.0 .and. tau_x470 > 0.0) xxrat = tau_x412 / tau_x470
2527  if (tau_x650 > 0.0 .and. tau_x470 > 0.0) xxrat2 = tau_x470 / tau_x650
2528 
2529 ! -- just a hack to avoid AE changes now that AOT values can be extrapolated
2530 ! -- past 3.5 whereas before they couldn't and would be truncated to 3.5.
2531  if (tau_x412 > 0.0 .AND. tau_x470 > 0.0) xxrat = min(tau_x412, 3.5)/min(tau_x470, 3.5)
2532  if (tau_x650 > 0.0 .and. tau_x470 > 0.0) xxrat2 = min(tau_x470, 3.5)/min(tau_x650, 3.5)
2533 
2534  if (xxrat < 0.0) then
2535  alpha = -999.
2536  go to 806
2537  end if
2538 
2539  dd = alog(412./490.)
2540  alpha = alog(xxrat)
2541  alpha = -1.*alpha/dd
2542 
2543  806 continue
2544 
2545  dd = ref650 / rr470_mod
2546  dd1 = ref650 / rr412_mod
2547  dd2 = rr470_mod / rr412_mod
2548  sfcdd = r650 / r412
2549  view = xthet
2550 
2551  if (r650 < 30.0 .and. r650 > 9.0 .and. sfcdd < 2.4) go to 850
2552  if (r650 < 30.0 .and. r650 > 9.0 .and. sfcdd > 2.4 .and.
2553  1 sfcdd < 3.2 .and. trflg > 0.) go to 850
2554 
2555  aot = aot_mod(3)
2556 
2557  if (px_elev > 500.0 .and. r650 > 10.0 .and. xlong < -30.0) aot = tau_x412
2558  if (px_elev > 500.0 .and. 1./rat_650_470 < 0.6 .and. tau_x650 > 0.6) then
2559  if (xlong < -30.0 .and. xlat < 10.0) go to 10
2560  if (xlong < -30.0 .and. xlat >= 10.0 .and. lc /= 6) go to 10
2561  endif
2562 
2563  if (tau_x412_91 < 0.0 .and. tau_x650 > 0.0 .and. r650 < 30.0) aot = aot_mod(1)
2564  if (tau_x412_91 < 0.0 .and. r650 >= 30.0 .and. tau_x412 > 0.0) aot = aot_mod(4)
2565  if (aot < 0.3 .and. aot > 0.0) go to 860
2566  if (w0_int < 0.92 .and. aot > 0.0) aot = aot_mod(2)
2567  if (r650 < 30.0 .and. tau_x412_91 > 1.2
2568  1 .and. tau_x412_91 > tau_x650 .and. tau_x650 > tau_x470
2569  1 .and. w0_int >= 0.939) aot = aot_mod(1)
2570  if (tau_x412_91 < 0.0 .and. tau_x470 < 0.0 .and. tau_x650 > 0.0) aot = aot_mod(1)
2571  if (aot < 0.0 .and. tau_x650 > 0.0) aot = aot_mod(1)
2572  if (tau_x412 < 0.0 .and. tau_x470 < 0.0 .and. dd1 > 1.35) alpha = -0.4
2573  if (tau_x412 < 0.0 .and. tau_x470 < 0.0 .and. dd1 <= 1.35) alpha = 1.8
2574  if (tau_x412 < 0.0 .and. tau_x470 > 0.0 .and. dd1 <= 1.4) alpha = 1.8
2575  if (aot < 0.0) go to 10
2576  go to 860
2577 
2578 850 continue
2579 
2580  aot = aot_mod(2)
2581 
2582  if (tau_x650 > aot .and. xlong > -30.0) aot = aot_mod(1)
2583  if (px_elev > 500.0 .and. r650 > 10.0 .and. xlong < -30.0) aot = tau_x412
2584  if (px_elev > 500.0 .and. 1./rat_650_470 < 0.6 .and. tau_x650 > 0.6) then
2585  if (xlong < -30.0 .and. xlat < 10.0) go to 10
2586  if (xlong < -30.0 .and. xlat >= 10.0 .and. lc /= 6) go to 10
2587  endif
2588 
2589  if (aot < 0.0) go to 10
2590  if (gzflg /= 16 .and. tau_x650 > 1.0 .and. dd2 < 1.1 .and. r650 > 15.) go to 10
2591 
2592  if (tau_x412 > 0.0 .and. tau_x650 > 1.9 .and. dd2 < 1.2) go to 10
2593  if (tau_x412 > 0.0 .and. tau_x470 > 1.9 .and. tau_x650 > 0.4
2594  1 .and. dd2 < 1.2 .and. gzflg /= 31 .and. gzflg /= 13) go to 10
2595  if (gzflg /= 16 .and. tau_x412 > 0.0 .and. aot > 1.0 .and. dd2 < 1.1 .and.
2596  1 tau_x650 > 0.4 .and. r650 > 13.0 .and. r650 <= 15.0 .and. gzflg /= 31
2597  1 .and. gzflg /= 13) go to 10
2598 
2599  if (tau_x412 > 0.0 .and. r412 > 18.0 .and. tau_x650 > 0.4
2600  1 .and. dd < 1.4 .and. r650 <= 15.) go to 10
2601 
2602  if (gzflg /= 16 .and. tau_x412 > 1.0 .and. xxrat > 1.2 .and. r650 > 15.0
2603  1 .and. tau_x650 > 0.4 .and. gzflg /= 31 .and. gzflg /= 13) go to 10
2604 
2605  if (tau_x650 > 1.5 .and. xxrat > 1.05 .and. tau_x412 < tau_x650 .and. tau_x412 > 0.0) go to 10
2606  if (tau_x650 > 1.5 .and. tau_x412 < 0.0 .and. tau_x470 < 0.0 .and.
2607  1 dd2 < 1.6 .and. w0_x412 > 0.96) go to 10
2608  if (tau_x650 > 1.2 .and. tau_x412 < 0.0 .and. xxrat2 > 1.2 .and. w0_x412 > 0.97) go to 10
2609 
2610 
2611  if (tau_x412 < 0.0 .and. tau_x470 < 0.0 .and. dd1 > 1.35) alpha = -0.4
2612  if (tau_x412 < 0.0 .and. tau_x470 < 0.0 .and. dd1 <= 1.35) alpha = 1.8
2613  if (tau_x412 < 0.0 .and. tau_x470 > 0.0 .and. dd1 <= 1.4) alpha = 1.8
2614 
2615  if (tau_x650 > 2.0 .and. tau_x412 < 0.0 .and. tau_x470 > 2.0
2616  1 .and. tau_x470 < 3.0 .and. xxrat2 < 1.2 .and. xxrat2 > 1.0)
2617  1 go to 10
2618 
2619  if (tau_x650 > 2.0 .and. tau_x412 < 0.0 .and. tau_x470 > 3.0
2620  1 .and. xxrat2 < 1.45 .and. xxrat2 > 1.0) go to 10
2621 
2622  if (tau_x412 < 0.0 .and. tau_x470 > 0.0 .and. xxrat2 > 2.)
2623  1 aot = aot_mod(5)
2624  if (tau_x412 > 1.5 .and. tau_x470 > 0.0 .and. tau_x650 < 0.3)
2625  1 aot = aot_mod(5)
2626  if (alpha > 1.0 .and. tau_x470 > 0.2) aot = aot_mod(5)
2627  if (alpha > 1.0 .and. tau_x470 <= 0.2) aot = aot_mod(5)*0.75
2628  if (tau_x412 < 0.0 .and. tau_x470 < 0.0 .and. dd1 <= 1.35) aot = aot_mod(6)
2629 
2630 860 continue
2631  tau550 = aot
2632 c --- additional cloud screening
2633 ! if ((gzflg <= 6 .OR. gzflg == 27) .and. (gzflg > 0 .and. gzflg /= 2)) then
2634 ! if (tau_x412ss > 0.6.and.rat1 <1.4.and.Dstar1 < 1.0.and.r412_tbl < 10.) go to 10
2635 ! endif
2636 
2637  ! -- skip all of the stuff below, just use aot.
2638  if (abs_aero_flag .eqv. .true.) then
2639  goto 864
2640  end if
2641 
2642  ! Overwrite DB AOD (NDVI < NDVI_thold and no AERONET BRDF)
2643  ! Default, not using default because we do not know the impact on other regions
2644 ! aot = tau_x470_new_96
2645 ! if (aot < 0.7) then
2646 ! model_frac = 1.0-aot/0.7
2647 ! aot = tau_x470_new_96*model_frac + tau_x470_new_94*(1.0-model_frac)
2648 ! else
2649 ! aot = tau_x470_new_94
2650 ! endif
2651 
2652  ! North America, test needed
2653 ! if (tau_x470sv_94 > 0.0 .and. tau_x470sv_96 > 0.0 .and. tau_x470sv_995 > 0.0) then
2654 ! if (regid == 1) then
2655 ! aot = tau_x470sv_995
2656 ! if (aot < 0.7) then
2657 ! model_frac = 1.0-aot/0.7
2658 ! aot = tau_x470sv_995*model_frac + tau_x470sv_96*(1.0-model_frac)
2659 ! else
2660 ! aot = tau_x470sv_96
2661 ! endif
2662 ! endif
2663 ! endif
2664 
2665  ! Europe
2666  if (regid == 8 .and. gzflg < 0) then
2667  aot = tau_x470_new_995
2668  if (aot < 0.7) then
2669  model_frac = 1.0-aot/0.7
2670  aot = tau_x470_new_995*model_frac + tau_x470_new_96*(1.0-model_frac)
2671  else
2672  aot = tau_x470_new_96
2673  endif
2674  endif
2675 
2676  ! Korea and Japan
2677  if (regid == 9 .and. gzflg < 0) then
2678  aot = tau_x470_new_96
2679  if (aot < 0.7) then
2680  model_frac = 1.0-aot/0.7
2681  aot = tau_x470_new_96*model_frac + tau_x470_new_94*(1.0-model_frac)
2682  else
2683  aot = tau_x470_new_94
2684  endif
2685  endif
2686  tau550 = aot
2687 
2688  ! -- force Tinga_Tingana scheme for AOT @ 500nm over barren surfaces even when gzflag is undefined.
2689  if (lc == 6 .AND. (gzflg /= 16 .AND. gzflg /= 2) .AND. gzflg /= 14 .AND. gzflg /= 21 .AND.
2690  1 gzflg /= 10 .AND. gzflg /= 20 .AND. gzflg /= 30 .and. gzflg/=31) then
2691  tau550 = get_aot500(xlat, xlong, 0.0, scat_ang, season, toa_ndvi, 12, lc, stdv,
2692  1 tau_x412_new_91, tau_x412_new_93, tau_x412_new_94,
2693  1 tau_x412_new_96, tau_x412_new_995, tau_x470_new_91,
2694  1 tau_x470_new_92, tau_x470_new_93, tau_x470_new_94,
2695  1 tau_x470_new_95, tau_x470_new_96, tau_x470_new_995,
2696  1 tau_x412_dust_91, tau_x412_dust_93, tau_x412_dust_94,
2697  1 tau_x412_dust_96, tau_x412_dust_995,tau_x470_dust_91,
2698  1 tau_x470_dust_92, tau_x470_dust_93, tau_x470_dust_94,
2699  1 tau_x470_dust_95, tau_x470_dust_96, tau_x470_dust_995,
2700  1 alpha, status, (lprint > 0))
2701  if (status /= 0) then
2702  tau550 = aot
2703  endif
2704  endif
2705 
2706  if (gzflg > 0) then
2707 ! -- get AOT @ 500nm based on AERONET regions and case studies
2708  tau550 = get_aot500(xlat, xlong, px_elev, scat_ang, season, toa_ndvi, gzflg, lc, stdv,
2709  1 tau_x412_new_91, tau_x412_new_93, tau_x412_new_94,
2710  1 tau_x412_new_96, tau_x412_new_995, tau_x470_new_91,
2711  1 tau_x470_new_92, tau_x470_new_93, tau_x470_new_94,
2712  1 tau_x470_new_95, tau_x470_new_96, tau_x470_new_995,
2713  1 tau_x412_dust_91, tau_x412_dust_93, tau_x412_dust_94,
2714  1 tau_x412_dust_96, tau_x412_dust_995,tau_x470_dust_91,
2715  1 tau_x470_dust_92, tau_x470_dust_93, tau_x470_dust_94,
2716  1 tau_x470_dust_95, tau_x470_dust_96, tau_x470_dust_995,
2717  1 alpha, status, (lprint > 0))
2718  if (status /= 0) then! revert to manual assignment
2719  tau550 = aot
2720  endif
2721 
2722 ! -- force Kanpur (gzflg=15, lc=3) AOT models over gzflg 20 (Thar Desert)
2723 c if (gzflg == 20) then
2724 c tau550 = get_aot500(xlat, xlong, 0.0, scat_ang, season, toa_ndvi, 15, 3, stdv,
2725 c 1 tau_x412_new_91, tau_x412_new_93, tau_x412_new_94,
2726 c 1 tau_x412_new_96, tau_x412_new_995, tau_x470_new_91,
2727 c 1 tau_x470_new_92, tau_x470_new_93, tau_x470_new_94,
2728 c 1 tau_x470_new_95, tau_x470_new_96, tau_x470_new_995,
2729 c 1 alpha, status, (lprint > 0))
2730 c if (status /= 0) then
2731 c tau550 = aot
2732 c endif
2733 c endif
2734 
2735 ! -- force Beijing (gzflg=16, low elev) AOT models over gzflg 21 (Asian desert)
2736  if (gzflg == 21) then
2737  tau550 = get_aot500(xlat, xlong, 0.0, scat_ang, season, toa_ndvi, 16, 4, stdv,
2738  1 tau_x412_new_91, tau_x412_new_93, tau_x412_new_94,
2739  1 tau_x412_new_96, tau_x412_new_995, tau_x470_new_91,
2740  1 tau_x470_new_92, tau_x470_new_93, tau_x470_new_94,
2741  1 tau_x470_new_95, tau_x470_new_96, tau_x470_new_995,
2742  1 tau_x412_dust_91, tau_x412_dust_93, tau_x412_dust_94,
2743  1 tau_x412_dust_96, tau_x412_dust_995,tau_x470_dust_91,
2744  1 tau_x470_dust_92, tau_x470_dust_93, tau_x470_dust_94,
2745  1 tau_x470_dust_95, tau_x470_dust_96, tau_x470_dust_995,
2746  1 alpha, status, (lprint > 0))
2747  if (status /= 0) then
2748  tau550 = aot
2749  endif
2750  endif
2751 
2752 ! -- force SW Asia (Pakistan, Iraq, etc) (gzflg=23) to use Tinga_Tingana
2753  if (gzflg == 23) then
2754  tau550 = get_aot500(xlat, xlong, 0.0, scat_ang, season, toa_ndvi, 12, 6, stdv,
2755  1 tau_x412_new_91, tau_x412_new_93, tau_x412_new_94,
2756  1 tau_x412_new_96, tau_x412_new_995, tau_x470_new_91,
2757  1 tau_x470_new_92, tau_x470_new_93, tau_x470_new_94,
2758  1 tau_x470_new_95, tau_x470_new_96, tau_x470_new_995,
2759  1 tau_x412_dust_91, tau_x412_dust_93, tau_x412_dust_94,
2760  1 tau_x412_dust_96, tau_x412_dust_995,tau_x470_dust_91,
2761  1 tau_x470_dust_92, tau_x470_dust_93, tau_x470_dust_94,
2762  1 tau_x470_dust_95, tau_x470_dust_96, tau_x470_dust_995,
2763  1 alpha, status, (lprint > 0))
2764  if (status /= 0) then
2765  tau550 = aot
2766  endif
2767  endif
2768 
2769 ! -- force Tinga_Tingana scheme for AOT @ 500nm over barren surfaces
2770  if (lc == 6 .AND. (gzflg /= 16 .AND. gzflg /= 2) .AND. gzflg /= 14 .AND. gzflg /= 20 .AND.
2771  1 gzflg /= 21 .AND. gzflg /= 23 .AND. gzflg /= 24 .AND. gzflg /= 10 .AND. gzflg /= 30 .and. gzflg/=31) then
2772  tau550 = get_aot500(xlat, xlong, 0.0, scat_ang, season, toa_ndvi, 12, lc, stdv,
2773  1 tau_x412_new_91, tau_x412_new_93, tau_x412_new_94,
2774  1 tau_x412_new_96, tau_x412_new_995, tau_x470_new_91,
2775  1 tau_x470_new_92, tau_x470_new_93, tau_x470_new_94,
2776  1 tau_x470_new_95, tau_x470_new_96, tau_x470_new_995,
2777  1 tau_x412_dust_91, tau_x412_dust_93, tau_x412_dust_94,
2778  1 tau_x412_dust_96, tau_x412_dust_995,tau_x470_dust_91,
2779  1 tau_x470_dust_92, tau_x470_dust_93, tau_x470_dust_94,
2780  1 tau_x470_dust_95, tau_x470_dust_96, tau_x470_dust_995,
2781  1 alpha, status, (lprint > 0))
2782  if (status /= 0) then
2783  tau550 = aot
2784  endif
2785 
2786  if ((gzflg <= 6 .OR. gzflg == 27) .and. (gzflg > 0 .and. gzflg /= 2)) then
2787  tau550 = tau_x412ss
2788  if (dstar1 > 1.08 .and. tau_x470ss > tau550) tau550 = tau_x470ss
2789  if (xday > 151.0.and.xday < 244.0) then ! Jun, Jul, Aug
2790  tau550 = tau_x412ss
2791  if (dstar1 > 1.09.and. wv <= 1.5 .and. tau_x470ss > tau550) tau550 = tau_x470ss
2792  if (dstar1 > 1.07 .and. tau_x470ss/tau_x412ss > 3.5) tau550 = tau_x470ss
2793  endif
2794  ddx = 0.0
2795  if (xday > 243.0.and.xday < 274.0) ddx = 1.0 ! Sep
2796  if (xday > 120.0.and.xday < 152.0) ddx = 1.0 ! May
2797  if (xday > 151.0.and.xday < 244.0) ddx = 1.0 ! Jun, Jul, Aug
2798  if (ddx > 0.0) then
2799  if (xday > 243.0.and.xday < 274.0) then
2800  if (dstar1 > 1.01 .and. r412_tbl <12.0 .and. tau_x470ss > tau550) tau550 = tau_x470ss
2801  if (dstar1 > 0.98 .and. wv > 1.7 .and. tau_x470ss > tau550) tau550 = tau_x470ss
2802  endif
2803  if (tau_x470ss > tau550) then
2804  if (xlong <= -5.0) tau550 = tau_x470ss
2805  if (xlong > -5.0 .and. xlong < 0.0) then
2806  dd = (5.+xlong) / 5.
2807  tau550 = (1.-dd)* tau_x470ss + dd * tau_x412ss
2808  endif
2809  if (xday > 120.0.and.xday < 152.0) then
2810  if (xlong <= -5.0) tau550 = tau_x470ss
2811  if (xlong > -5.0 .and. xlong < 10.0) then
2812  dd = (5.+xlong) / 15.
2813  tau550 = (1.-dd)* tau_x470ss + dd * tau_x412ss
2814  endif
2815  endif
2816  endif
2817  endif
2818  if (tau_x470ss > tau550) then
2819  if (xlong >= 35.0) tau550 = tau_x470ss
2820  if (xlong > 27.0 .and. xlong < 35.0) then
2821  dd = (xlong-27.0) / 8.
2822  tau550 = dd* tau_x470ss + (1.-dd) * tau_x412ss
2823  endif
2824  endif
2825 
2826 
2827  if (tau550 < 0.1) tau550 = tau550 + 0.05 ! check for geo zone
2828 
2829  if (dstar1 > 1.1 .and. tau_x470ss > tau_x412ss) tau550 = tau_x470ss !new
2830 
2831  if (xday > 151.0.and.xday < 258.0.and.
2832  1 dstar1 > 1.2.and.tau_x650 < 0.0.and.tau550 < 0.8)
2833  1 tau550 = tau_x470ss * 2.
2834 
2835  if (dstar1 > 1.2.and.tau_x412ss < 0.5 .and.tau550< 0.7)
2836  1 tau550 = dstar1
2837 
2838  if (wv < 0.45.and.tau_x470ss >0.9.and.tau_x412ss < 0.7) tau550 = tau_x412ss !new
2839  if (dstar1 > 1.08 .and. tau_x470ss > tau550) tau550 = tau_x470ss
2840 
2841  end if
2842 
2843  if (gzflg >= 6 .and. gzflg <= 11 .AND. gzflg /= 10) then
2844  tau550 = tau_x470ss
2845  if (tau550 < 0.45) tau550 = (tau_x470ss+tau_x412ss2)/2.
2846  if (xday > 243.0 .and. xday < 274.0) then
2847  tau550 = tau_x470ss
2848  if (tau550 < 0.45.and. tau_x412ss2 > tau550) tau550 = (tau_x470ss+tau_x412ss2)/2.
2849  endif
2850  if (xday > 120.0 .and. xday < 151.0) then ! May
2851  tau550 = tau_x470ss
2852  if (tau550 < 0.45.and. tau_x412ss2 > tau550) tau550 = (tau_x470ss+tau_x412ss2)/2.
2853  end if
2854  if (xday > 90.0 .and. xday < 121.0) then ! Apr
2855  tau550 = tau_x412ss2
2856  if (tau550 < 0.5 .and. tau_x470ss > tau550) tau550 = tau_x470ss
2857  end if
2858  if (xlat > 29.5 .and. dstar1 < 0.97) then
2859  if (tau550 > 0.6) tau550 = tau_x412ss2_995
2860  end if
2861 
2862  if (xday < 182.0) then
2863  if (xlat > 20.0.and. xlat < 28.0.and.xlong > 42.0.and. xlong < 50.0) then
2864  if (dstar1 < 1.06.and.r412_tbl > 14.0.and.tau550 > 0.8) go to 10
2865  end if
2866  end if
2867 
2868  if (dstar1 > 1.1 .and. tau_x470ss <0.45) tau550 = dstar1
2869 
2870  if (xday > 243.0 .and. xday < 274.0) then ! Sept
2871  if (dstar1 < 0.98.and.r412_tbl > 12.0.and.tau550 > 0.8) go to 10
2872  endif
2873  if (xday > 273.0 .and. xday < 305.0) then ! Oct
2874  if (dstar1 < 1.04.and.r412_tbl > 12.0.and.tau550 > 0.75) go to 10
2875  endif
2876 
2877  if (xday > 150.0 .and. xday < 258.0) then ! Jun, Jul, Aug
2878  if (xlat > 20.0.and. xlat < 28.0.and.xlong > 42.0.and. xlong < 50.0) then
2879  if (dstar1 < 1.06.and.r412_tbl > 12.8.and.tau550 > 0.78) go to 10
2880  end if
2881  if (dstar1 < 1.04.and.r412_tbl > 12.0.and.tau550 > 0.78.and. xlat>21.0.and.xlong<53.0) go to 10
2882  if (dstar1 < 1.04.and.r412_tbl > 12.0.and.tau550 > 0.78.and. xlat>23.0.and.xlong>=53.0) go to 10
2883  end if
2884 
2885  end if
2886 
2887  end if ! end of barren
2888  if (gzflg < 6 .and. gzflg >0 .and. dstar1 > 1.08 .and. tau_x470ss > tau550) tau550 = tau_x470ss
2889 
2890  end if
2891 
2892 ! -- for Arabian peninsula, just use AOT412 w/ w0=0.96. Override all of that above.
2893  if (gzflg >= 6 .and. gzflg <= 11 .AND. gzflg /= 10) then
2894  if ((r650_135 > 32.0) .AND. (r650_135/r412_135) > 3.7) then
2895  if (xlat > 20.5 .and. xlong < 50.0) go to 861
2896  if (xlat > 19.0.and. xlat < 20.5 .and. xlong > 43.0 .and. xlong < 48.0) go to 861
2897  if ((r412_135 < 7.25) .OR. r650_135 > 42.0 .OR. dstar1 >1.03 .OR. use_alternate_brdf) then
2898  continue
2899 
2900  else if (xday >= 152.0 .AND. xday <= 243.0) then ! summer
2901  if (dstar1 < 0.99) tau550 = (tau_x412_new_94+tau_x412_new_96)/2.0 ! orig = 0.96
2902 
2903  else if (xday >= 244.0 .AND. xday < 335.0) then
2904  tau550 = (tau_x412_new_94+tau_x412_new_96)/2.0 ! orig = 0.995
2905 
2906  else
2907  tau550 = tau_x412_new_94
2908  if (xday >= 120.0 .AND. xday < 152.0.AND.tau550<0.3.AND.tau_x470ss>tau550) tau550 = tau_x470ss
2909  end if
2910 
2911  end if
2912  end if
2913  861 continue
2914 
2915  if (lprint > 0) print *, 'after get_aot500: xlat, xlong, tau550: ', gzflg, xlat, xlong, tau550
2916 ! Weakly absorbing dust
2917 ! if (Dstar1 < 1.01.and.tau550 > 0.8.and.rat1 < 1.6.and.
2918 ! 1 gzflg <6 .and. gzflg >0)! .and. lc .eq. 6)
2919 ! 1 tau550 = tau_x412ss_98
2920 ! if (lprint > 0) print *, 'after absorbing dust check: tau550, Dstar1, rat1: ',
2921 ! 1 tau550, Dstar1, rat1, rr470_mod, rr412_mod
2922 
2923 
2924 ! Smoke pixels
2925  if (((gzflg <6 .and. gzflg >0) .OR. gzflg == 27) .and. xlat <20.0 .and. lc <5) then
2926  if (xday > 335.0 .or. xday < 32.0) then
2927  if (tau550 < 0.065 .or. tau_x470 <0.05) then
2928  dd412 = rr412_mod - r412_tbl
2929  dd470 = rr470_mod - r470_tbl
2930  dd650 = ref650 - r650_tbl
2931  if (debug) print *, 'smoke check, dd412, dd470, dd650: ', dd412, dd470, dd650
2932  !if (dd412 < 0.0 .AND. dd470 < 0.0 .AND. dd650 < 0.0) then
2933  if (dd412 < 1.0) then
2934  call smoke_mod(xthet,scat_ang,dd412,dd470,dd650,
2935  1 tau_x412ss,tau_x412ss_91,tau_smoke)
2936  if (debug) print *, 'smoke_mod, dd412, dd470, dd650,t412ss, t412_91, tsmoke: ',
2937  1 dd412, dd470, dd650, tau_x412ss, tau_x412ss_91, tau_smoke
2938  tau550 = tau_smoke
2939  if (tau_smoke < 0.25) then
2940  tau550 = tau_smoke *2.0
2941  if (debug) print *, 'smoke detected, aot550 reset: ', tau550
2942  end if
2943  end if
2944  end if
2945  end if
2946  end if
2947  if (lprint > 0) print *, 'after smoke check: tau550: ', tau550
2948 
2949  if (lprint > 0) print *,
2950  1 'gzflg,lc, tau_x412ss,tau_x470ss,tau550 =',
2951  1 gzflg,lc,tau_x412ss,tau_x470ss,tau550
2952 
2953 864 if (lprint > 0) then
2954  print *, 'lat,lon,gzflg,r650,tau_x650,tau_x470,tau550: ', xlat, xlong, gzflg, r650,
2955  1 tau_x650, tau_x470, tau550
2956  print *, 'abs_aero_flag: ', abs_aero_flag
2957  end if
2958 
2959 c -- try to detect and remove smoke over clouds. for example, see rgb:
2960 c -- myd021km.a2007081.0640*.hdf
2961  if (debug) print *, 'smoke detection: ', xlat, xlong, r650, tau_x650, refl6, w0_x470
2962  if (r650 <8.0 .and. tau_x650 > 3.49 .and. refl6 > 0.1 .and.w0_x470 > 0.999) go to 10
2963 ! if (xday > 0.0 .and. xday < 121.0) then
2964 ! if (gzflg==16 .and. (w0_x470 > 0.96 .and. alpha < 1.0) .and. Dstar1 <1.0 .and. tau550 >1.0) then
2965 ! go to 10
2966 ! end if
2967 ! endif
2968 
2969 ! -- set tau550 to tau_x470ss over high elevation to minimize the effects of
2970 ! -- Rayleigh scattering.
2971 ! -- Morocco region
2972  if (px_elev > 750.0 .AND. (xlat > 28.0 .AND. xlat < 37.0 .AND. xlong > -12.0 .AND.
2973  1 xlong < 10.0)) then
2974  tau550 = tau_x470
2975  endif
2976 
2977 ! -- West Sudan region
2978  if (px_elev > 750.0 .AND. (xlat > 10.5 .AND. xlat < 19.5 .AND. xlong > 20.5 .AND.
2979  1 xlong < 29.0)) then
2980  tau550 = tau_x470
2981  endif
2982 
2983 ! -- Somalia
2984  if (lc == 6 .AND. gzflg < 6 .AND. (xlat > -5.0 .AND. xlat < 18.0 .AND. xlong > 35.0 .AND.
2985  1 xlong < 52.0)) then
2986  tau550 = tau_x470ss
2987  endif
2988 
2989 ! -- filter out dark Iraq oil fire plumes
2990 ! 1/31/2017
2991  if ((rr412_mod - r412_135) < -1.0.and.px_elev < 1000.0 .and.
2992  1 gzflg >= 6 .and. gzflg <= 11) go to 10
2993 
2994 ! -- filter out noisy high AOD over US desert
2995  if (gzflg == 13 .and. lc == 6 .and. tau550 > 1.2 .and. dstar1 < 1.0 .and.
2996  1 realbuf(6)/realbuf(22) > 1.65) go to 10
2997 
2998  if (gzflg_sav == 31 .and. tau550 > 0.6 .and. dstar1 < 1.1 .and.
2999  1 realbuf(22)/realbuf(6) < 0.5) go to 10
3000 
3001 c -- check d* for dust plumes globally except for n.africa(1-5, 26, 27),
3002 c -- arabian peninsula(6-11), taklimakan desert(24), and beijing(east china, 16)
3003 c -- zones. limit to locations at less than 4500 m elevation. otherwise,
3004 c -- abnormally high aot's are seen over Tibetan Plateau.
3005 c ----------------------------------------------------------------------------
3006 c ----------------------------------------------------------------------------
3007 c -- Increased threshold to 1.15 from 1.05 over globe. Created exception for
3008 c -- ConUS (1.05). Needs to be tested and results analysed after next SIPS
3009 c -- large-scale test!
3010 c ----------------------------------------------------------------------------
3011 c ----------------------------------------------------------------------------
3012 .OR..AND..AND..AND. if ((gzflg < 0 gzflg > 11) (gzflg /= 24 gzflg /= 26
3013 .AND. & gzflg /= 27 gzflg /= 16)) then
3014  if (px_elev < 4750.0) then
3015  dda1 = 1.15
3016 .and. if (gzflg_sav == 13 realbuf(22)/realbuf(6) > 0.3) dda1 = 1.05
3017 .and. if (gzflg_sav == 31 realbuf(22)/realbuf(6) > 0.5) dda1 = 1.05
3018  if (Dstar1 > dda1) tau550 = tau_x650
3019 .and. if (Dstar1 > dda1 tau_x470 > tau550) tau550 = tau_x470
3020 .and..and. if (gzflg_sav == 31 gzflg == 13 realbuf(22)/realbuf(6) > 0.3
3021 .and..and. 1 Dstar1 > 1.05 tau_x470 > tau550) tau550 = tau_x470
3022  end if
3023  end if
3024 
3025  xtau(1) = tau_x412
3026  xtau(2) = tau_x470
3027  xtau(3) = tau_x650
3028 
3029  if (w0_x412 > 0.0) w0_x = w0_x412
3030  if (w0_x470 > 0.0) w0_int_470 = w0_x470
3031 
3032  if (alpha < 1.0) ssa(1) = w0_x
3033  if (alpha < 1.0) ssa(2) = w0_int_470
3034  if (ssa(1) > 0.0) ssa(3) = 0.976
3035 c-- Set Surface Type
3036 
3037  sfc_typ = -999.
3038 
3039 .gt. if (terrain_flag_new50.0) sfc_typ = 7.
3040 .lt..and..gt. if (terrain_flag_new55.5terrain_flag_new50.0)
3041  1 sfc_typ = 10.
3042 
3043 C dd = xsfc650(ilon,ilat) - xsfc650_bk2(ilon,ilat)
3044 C dd1 = xsfc470b_bk(ilon,ilat) - xsfc470_bk(ilon,ilat)
3045 
3046 .lt..and..gt..and. if (terrain_flag_new55.5terrain_flag_new50.0
3047 .le..and..lt. 1 xlat15.0xlong26.0) then
3048  sfc_typ = 1.
3049 .gt..and..ge.C if (xlong10.0abs(dd1)0.2) then
3050 .gt.C if (xsfc470_bk(ilon,ilat)10.0) sfc_typ = 2.
3051 .le..and.C if (xsfc470_bk(ilon,ilat)10.0
3052 .gt..and.C 1 xsfc470_bk(ilon,ilat)9.0
3053 .lt.C 1 xsfc650_bk(ilon,ilat)22.0) sfc_typ = 2.
3054 C endif
3055 .gt.C if (xlong10.0) then
3056 .gt..and.C if (xsfc470_bk(ilon,ilat)10.0
3057 .lt.C 1 xsfc650_bk(ilon,ilat)23.0) sfc_typ = 2.
3058 C endif
3059  endif
3060 .lt..and..gt..and.C if (terrain_flag_new55.5terrain_flag_new50.0
3061 .gt..and..le..and..ge.C 1 xlat15.0xlat20.0dd1.6) then
3062 C sfc_typ = 1.
3063 .gt..and..ge.C if (xlong10.0abs(dd1)0.2) then
3064 .gt.C if (xsfc470_bk(ilon,ilat)10.0) sfc_typ = 2.
3065 .gt..and.C if (xsfc470_bk(ilon,ilat)10.0
3066 .lt.C 1 xsfc650_bk(ilon,ilat)23.0) sfc_typ = 2.
3067 .le..and.C if (xsfc470_bk(ilon,ilat)10.0
3068 .gt..and.C 1 xsfc470_bk(ilon,ilat)9.0
3069 .lt.C 1 xsfc650_bk(ilon,ilat)22.0) sfc_typ = 2.
3070 C endif
3071 C endif
3072 .lt..and..gt..and.C if (terrain_flag_new55.5terrain_flag_new50.0
3073 .ge..and..le..and..ge.C 1 xlong26.0xlat20.0dd1.6) then
3074 C sfc_typ = 1.
3075 .gt..and..ge.C if (xlong10.0abs(dd1)0.2) then
3076 .gt.C if (xsfc470_bk(ilon,ilat)10.0) sfc_typ = 2.
3077 .gt..and.C if (xsfc470_bk(ilon,ilat)10.0
3078 .lt.C 1 xsfc650_bk(ilon,ilat)23.0) sfc_typ = 2.
3079 .le..and.C if (xsfc470_bk(ilon,ilat)10.0
3080 .gt..and.C 1 xsfc470_bk(ilon,ilat)9.0
3081 .lt.C 1 xsfc650_bk(ilon,ilat)22.0) sfc_typ = 2.
3082 C endif
3083 C endif
3084 
3085 .gt..and..lt. if (terrain_flag_new52.5terrain_flag_new53.5)
3086  1 sfc_typ = 3.
3087 
3088 .gt..and..lt. if (terrain_flag_new51.5terrain_flag_new52.5) then
3089  sfc_typ = 5.
3090 .gt..and..lt.C if (xsfc470_bk(ilon,ilat)6.7xlong-7.9)
3091 C 1 sfc_typ = 4.
3092  endif
3093 
3094 .gt..and..lt..and.C if (terrain_flag_new54.5terrain_flag_new55.5
3095 .gt..and..lt..and.C 1 xlat31.3xlong-7.9xsfc470_bk(ilon,ilat)
3096 .gt..and..lt.C 1 6.7xsfc470_bk(ilon,ilat)9.5)
3097 C 1 sfc_typ = 4.
3098 
3099 .gt..and..lt. if (terrain_flag_new53.5terrain_flag_new54.5)
3100  1 sfc_typ = 6.
3101 
3102 .gt..and..lt. if (terrain_flag_new58.5terrain_flag_new59.5)
3103  1 sfc_typ = 8.
3104 
3105 .gt..and..lt. if (terrain_flag_new59.5terrain_flag_new510.5)
3106  1 sfc_typ = 9.
3107 
3108 .gt..and..lt..and. if (terrain_flag_new55.5terrain_flag_new56.5
3109 .gt..and..lt..and..gt..and. 1 xlat24.6xlat25.2xlong46.1
3110 .lt..and..ge..and. 1 xlong46.8r412_tbl7.0
3111 .le. 1 r412_tbl9.9)
3112  1 sfc_typ = 9.
3113 
3114 c -- Set Flags
3115 .gt..or..gt. if (tau_x650_flag0tau_x470_flag0) qa_flag(4)= 3
3116 .lt..and..lt..and. if (tau_x4120.0tau_x4700.0
3117 .lt. 1 tau_x6500.0) then
3118  qa_flag(1)= 0
3119  qa_flag(2)= 0
3120  qa_flag(3)= 0
3121  qa_flag(4)= 2
3122  endif
3123 .gt..or..gt..or. if (tau_x4120.0tau_x4700.0
3124 .gt. 1 tau_x6500.0) qa_flag(1)= 1
3125 .lt..or..gt. if (tau_x4120.05tau_x4700.05)
3126  1 qa_flag(2)= 1
3127 .ge..and..lt. if (tau_x4120.05tau_x4120.3)
3128  1 qa_flag(2)= 2
3129 .ge. if (tau_x4120.3)
3130  1 qa_flag(2)= 3
3131 .lt..and..gt. if (alpha0.5alpha0.0)
3132  1 qa_flag(3)= 1
3133 .ge..and..lt. if (alpha0.5alpha1.4)
3134  1 qa_flag(3)= 2
3135 
3136  865 continue
3137 
3138 ! -- in zone 1, in summer, sometimes the BRDF for low NDVI pixels is too high
3139 ! -- resulting in an abnormally low AOT. In these cases, perform a second retrieval
3140 ! -- using a constant fit (an alternate fit) for 470nm that will bring up
3141 ! -- those low areas.
3142 .AND..AND..AND. if (tau550 < 0.1 (gdt1%month >= 6 gdt1%month <= 8) gzflg == 1) then
3143 .NOT. if ( use_alternate_brdf) then
3144  use_alternate_brdf = .true.
3145  goto 5
3146  end if
3147  end if
3148 
3149 ! -- over Arabian Peninsula, trouble with some very low AOT's in certain areas. redo
3150 ! -- the retrieval if in summer.
3151  if (.NOT. use_alternate_brdf .AND. (tau550 < 0.25 .AND.
3152  & (gdt1%month >= 6 .AND. gdt1%month< = 8) .AND. gzflg >= 6 .AND. gzflg <= 11 .AND. gzflg /= 10)) then
3153  if (r650_135 > 32.0 .AND. (r650_135/r412_135) > 3.7) then
3154  use_alternate_brdf = .true.
3155  goto 5
3156  end if
3157  end if
3158 
3159 ! from MODIS code for swir-vis surface database
3160 c -- use tau_x470sv96_dust over high elevation
3161 c
3162  if (gzflg /= 15 .and. gzflg /= 19 .and. xlong> -20.0 .and. xlong< 70.0 .and.
3163  & xlat> 5.0 .and. xlat< 45.0) then
3164 
3165  if (gzflg < 6 .and. gzflg > 0) then
3166  if (dstar1 < 1.1) then
3167  if (r412_tbl > 12.0 .and. dstar1 < 0.95 .and. tau_x470sv96_dust > 0.6) tau_x470sv96_dust = -999.
3168  dd = 999.
3169  if (xday >= 91.0 .and. xday < 244.0) dd = 10.0
3170  if (xlat >= 20.0 .and. xlong <= dd .and. tau_x470sv96_dust > 0.0) then
3171  tau550 = tau_x470sv96_dust
3172  endif
3173 
3174  if (xday < 91.0 .or. xday >= 244.0) go to 871
3175  if (xlong >= 10.0 .and. xlong <15.0 .and. xlat > 15.0) then
3176  dd = (xlong - 10.0) / 5.
3177  if (tau_x470sv96_dust > 0.0 .and. tau550 > 0.0) then
3178  tau550 = tau_x470sv96_dust + (tau550 - tau_x470sv96_dust) * dd
3179  endif
3180  endif
3181 871 continue
3182 
3183  dd = 999.
3184  if (xday >= 91.0 .and. xday < 244.0) dd = 10.0
3185  if (xlat >= 15.0 .and. xlat <20.0 .and. xlong < dd) then
3186  dd = (xlat - 15.0) / 5.
3187  if (tau_x470sv96_dust > 0.0 .and. tau550 > 0.0) then
3188  tau550 = tau550 + (tau_x470sv96_dust - tau550) * dd
3189  endif
3190  endif
3191 
3192  tau550_sav = tau550
3193 
3194  if (tau_x470sv96_dust > tau550) then
3195  tau550 = tau_x470sv96_dust
3196  endif
3197 
3198  if (xlat >= 25.0 .and. xlong >= 20) then
3199  tau550 = tau550_sav
3200  endif
3201 
3202  if (xlong >= 10.0 .and. xlong <20.0 .and. xlat > 20.0) then
3203  dd = (xlong - 10.0) / 10.
3204  if (tau550_sav > 0.0 .and. tau550 > 0.0) then
3205  tau550 = tau550 + (tau550_sav - tau550) * dd
3206  endif
3207  endif
3208 
3209  if (xlat >= 20.0 .and. xlat <25.0 .and. xlong > 10) then
3210  dd = (xlat - 20.0) / 5.
3211  if (tau550_sav > 0.0 .and. tau550 > 0.0) then
3212  tau550 = tau550 + (tau550_sav - tau550) * dd
3213  endif
3214  endif
3215 
3216 
3217  if (r412_tbl > 12.0 .and. tau_x470sv96_dust < 0.0) then
3218  if (realbuf(23) > 25.0 .and. realbuf(23)/realbuf(6) > 0.7
3219  & .and. realbuf(23)/realbuf(6) < 0.85) go to 867
3220  tau550 = -999.0
3221 867 continue
3222  endif
3223  if (xday >= 152.0 .and. xday < 244.0 .and. dstar1 >= 1.04 .and. tau_x470sv96_dust > tau550) then
3224  tau550 = tau_x470sv96_dust
3225  endif
3226  endif
3227  if (dstar1 >= 1.1 .and. tau_x470sv96_dust > tau550) then
3228  tau550 = tau_x470sv96_dust
3229  endif
3230  if (dstar1 >= 1.1 .and. tau_x470sv96_dust <0.0 .and. tau550 < 0.3) then
3231  tau550 = tau_x470sv96_dust
3232  endif
3233  endif
3234 
3235  if (gzflg == 2 .or. (gzflg > 5 .and. gzflg < 12)) then
3236  tau550 = tau_x470sv96_dust
3237  endif
3238 
3239  if (gzflg == 30) then
3240  if (xlong > 65.0 .and. xlong< 70.0) then
3241  dd = (xlong - 65.0) / 5.
3242  if (tau_x470sv96_dust > 0.0 .and. tau550 > 0.0) then
3243  tau550 = tau550 + (tau_x470sv96_dust - tau550) * (1.-dd)
3244  endif
3245  if (tau_x470sv96_dust < 0.0) then
3246  tau550 = tau_x470sv96_dust
3247  endif
3248  else
3249  tau550 = tau_x470sv96_dust
3250  endif
3251  endif
3252 
3253  if (gzflg /= 30 .and. tau_x470sv96_dust > tau550 .and.
3254  & (px_elev > 750.0 .and. (xlat> 14.0. or. xlong< 5.5 .or. xlong> 19.0))
3255  & .and. (px_elev > 750.0 .and. (xlat> 14.0 .or. xlong > -6.37))) then
3256  tau550 = tau_x470sv96_dust
3257  endif
3258 
3259  if (gzflg == 32) then! Portugal and Spain
3260  if (season == 1) then
3261  if (tau_x470sv96 < 0.3) then
3262  tau550 = tau_x470sv96
3263  elseif (tau_x470sv96 < 0.6) then
3264  aod_frac = (tau_x470sv96-0.3)/(0.6-0.3)
3265  tau550 = (1.0-aod_frac)*tau_x470sv96+aod_frac*tau_x470sv94
3266  else
3267  tau550 = tau_x470sv94
3268  endif
3269  else
3270  if (tau_x470sv995 < 0.3) then
3271  tau550 = tau_x470sv995
3272  elseif (tau_x470sv995 < 0.6) then
3273  aod_frac = (tau_x470sv995-0.3)/(0.6-0.3)
3274  tau550 = (1.0-aod_frac)*tau_x470sv995+aod_frac*tau_x470sv96
3275  else
3276  tau550 = tau_x470sv96
3277  endif
3278  endif
3279  endif
3280 
3281  if (gzflg == 33) then! Palencia
3282  if (season == 1) then
3283  if (tau_x412sv96 < 0.3) then
3284  tau550 = tau_x412sv96
3285  elseif (tau_x412sv96 < 0.6) then
3286  aod_frac = (tau_x412sv96-0.3)/(0.6-0.3)
3287  tau550 = (1.0-aod_frac)*tau_x412sv96+aod_frac*tau_x412sv94
3288  else
3289  tau550 = tau_x412sv94
3290  endif
3291  else
3292  if (tau_x412sv98 < 0.3) then
3293  tau550 = tau_x412sv98
3294  elseif (tau_x412sv98 < 0.6) then
3295  aod_frac = (tau_x412sv98-0.3)/(0.6-0.3)
3296  tau550 = (1.0-aod_frac)*tau_x412sv98+aod_frac*tau_x412sv96
3297  else
3298  tau550 = tau_x412sv96
3299  endif
3300  endif
3301  endif
3302  endif
3303 
3304 ! 24 January 2018 JLee TEST
3305 ! print *, gzflg, px_elev
3306 ! if (gzflg == 15 .and. px_elev > 600.0) then! India
3307 ! tau550 = tau_x412sv94
3308 ! endif
3309 
3310 !end from MODIS code for swir-vis surface database
3311 
3312 ! jlee test 2.2 um AOD everywhere (global) 20170725
3313 ! tau550 = tau_x470sv
3314 ! end jlee test
3315 
3316 ! -- We only really want to extrapolate aot550 past 3.5. All other AOT values
3317 ! -- should be truncated at 3.5 to maintain current AE and SSA values.
3318  if (xtau(1) > 3.5) xtau(1) = 3.5
3319  if (xtau(2) > 3.5) xtau(2) = 3.5
3320  if (xtau(3) > 3.5) xtau(3) = 3.5
3321 
3322 ! -- To detect the smoke plumes for aerosol typing
3323  dda1 = 999.0
3324  dda2 = 0.0
3325  if (realbuf(22) > 0.0 .and. realbuf(22) < 50.0) then
3326  dda1 = (realbuf(22)*realbuf(6)) / (realbuf(23)*realbuf(23))
3327  endif
3328 
3329  ddx = 0.90
3330  if (xlong < -20.0 .or. (xlat <30.0 .and. xlat >-40.0 .and.
3331  1 xlong >-30.0 .and. xlong <60.0 .and. gzflg < 1) .or.
3332  1 gzflg == 26 .or. gzflg == 27 .or. (gzflg > 0 .and. gzflg < 6)) go to 868
3333  if (xphi >= 90.0 .and. xthet > 50.0 .and. realbuf(23) < 33.0) then
3334  ddx = 0.90 - (0.90-0.83)*(xthet-50.0)/15.0
3335  if (xthet > 65.0) ddx = 0.83
3336  endif
3337 868 continue
3338 
3339  if (tau550 > 0.4 .and. dda1 < ddx .and. bt11 > 286.0) dda2 = 1.0
3340 
3341  dda3 = 999.0
3342  if (realbuf(22) > 0.0 .and. realbuf(22) < 50.0) then
3343  dda3 = realbuf(23)-realbuf(22) - (realbuf(6)-realbuf(23))*(488.-412.)/(670.-488.)
3344  endif
3345  if (dda1 < 0.90 .and. dda1 > 0.86 .and. dda3 < 1.4) dda2 = 0.0
3346 
3347 ! -- check consistency when tau550 is undefined
3348  if (tau550 < -900.0) then
3349  xtau(1) = -999.0
3350  xtau(2) = -999.0
3351  xtau(3) = -999.0
3352  ssa(1) = -999.0
3353  ssa(2) = -999.0
3354  ssa(3) = -999.0
3355  tau550 = -999.0
3356  alpha = -999.0
3357  r412 = -999.0
3358  r470 = -999.0
3359  r650 = -999.0
3360  tau_x470_flag = -999
3361  tau_x650_flag = -999
3362  alg_flag = -999
3363  end if
3364 c-------------------------------------------------------------
3365 c set output buffer
3366 c-------------------------------------------------------------
3367  if (debug) print *, 'final spectra aot: ', xtau(1), xtau(2), xtau(3)
3368  if (debug) print *, 'final tau550, ae: ', tau550, alpha
3369  if (debug) print *, 'final SSA: ', ssa(1), ssa(2), ssa(3)
3370  do i=1,3
3371  outbuf(i) = xtau(i)
3372  outbuf(i+3) = ssa(i)
3373  enddo
3374 
3375  outbuf(7) = tau550
3376  outbuf(8) = alpha
3377  outbuf(9) = r412
3378  outbuf(10) = -999.0
3379  if (tau_x470_flag_veg == 1 .or. tau_x470_flag > 0) outbuf(10) = 1.0 ! out-of-bound flag
3380  !tau_x470_flag_veg : veg code outbound flag, tau_x470_flag : DB code outbound flag
3381 ! outbuf(10) = 1.0*tau_x650_flag
3382  outbuf(11) = r470
3383  outbuf(12) = r650
3384  outbuf(13) = xthet
3385  outbuf(14) = dda2
3386  outbuf(15) = sfc_typ
3387  outbuf(18) = float(alg_flag)
3388  outbuf(19) = outbufvg(21) ! veg bright surface flag
3389  outbuf(20) = outbufvg(7) ! veg AOD to be used for rcndvi_avg
3390  outbuf(21) = -999.0
3391  if (realbuf(22) > 0.0 .and. realbuf(22) < 50.0)
3392  1 outbuf(21) = realbuf(22)/realbuf(23) ! LER412/488nm ratio
3393 
3394 10 continue
3395  return
3396  end
3397 
3398 c --------------
3399 c
3400  subroutine search(dflag,xbar,x,n,i)
3401 c
3402 c purpose
3403 c locate position in table of point at which interpolation is
3404 c required
3405 c
3406 c usage
3407 c call search (xbar, x, n, i)
3408 c
3409 c description of parameters
3410 c xbar - point at which interpolation is required
3411 c x - array containing independent variable
3412 c n - number of points in x array
3413 c i - index specifying segment containing xbar
3414 c
3415  logical dflag
3416  dimension x(n)
3417  data b/.69314718/
3418  icnt = 0
3419  if (n.lt.2) go to 15
3420  if(x(1).gt.x(2)) go to 17
3421  m = int((log(float(n)))/b)
3422  i=2**m
3423  if (i .ge. n) i = n-1
3424  k=i
3425  10 k=k/2
3426  if (k .eq. 0) icnt = icnt + 1
3427  if (icnt .ge. 2) goto 27
3428  if (xbar.ge.x(i).and.xbar.lt.x(i+1))return
3429  if (xbar.gt.x(i)) go to 12
3430  i = i-k
3431  go to 10
3432  12 i = i+k
3433  if (i.lt.n) go to 10
3434  i=n-1
3435  go to 10
3436  15 print *, "Search n is less than 2."
3437  return
3438  17 print *, "Search table is not in increasing order."
3439  return
3440 
3441  27 continue
3442  do 22 i=1,n-1
3443  if (xbar.ge.x(i).and.xbar.le.x(i+1)) return
3444  22 continue
3445  write(6,*) "setting dflag = true"
3446  dflag = .true.
3447  return
3448 
3449  end
3450 
3451  subroutine search2(dflag,xbar,x,n,i,fx)
3452 c
3453 c call search (xbar, x, n, i)
3454 c
3455 c description of parameters
3456 c xbar - point at which interpolation is required
3457 c x - array containing independent variable
3458 c n - number of points in x array
3459 c i - index specifying segment containing xbar
3460 c fx - fraction from i(between i and i+1)
3461 c
3462  dimension x(n)
3463  logical dflag
3464 
3465  call search(dflag,xbar,x,n,i)
3466 
3467  fx = (xbar-x(i))/ (x(i+1)-x(i))
3468 
3469  end
3470 
3471 
3472 
3473 c------------------------------------------------------------
3474  subroutine aero_mod (tau_x412,tau_x470,tau_x650,
3475  1 tau_x412_91,aot_mod)
3476 
3477  real aot_mod(6)
3478 
3479  aot_mod(1) = tau_x650
3480  aot_mod(2) = tau_x470
3481  aot_mod(3) = tau_x412_91
3482  aot_mod(4) = tau_x412*1.9
3483  aot_mod(5) = tau_x470*2.
3484  aot_mod(6) = tau_x650*2.8
3485  return
3486  end
3487 
3488  integer function handle_lut_out_of_bounds(geo_zone, retr_flag, aot) result(status)
3489  implicit none
3490 
3491  integer, intent(in) :: geo_zone
3492  integer, intent(in) :: retr_flag
3493  real, intent(inout) :: aot
3494 
3495  status = 0
3496 
3497  if (retr_flag == -10) then ! -10 = return value from aero_412, etc indicating
3498  select case (geo_zone) ! refl < yy(1) scenario.
3499  case (5:11, 15, 16, 19, 20, 21, 23, 24, 26, 27)
3500  aot = 0.06
3501  case (1:4, 12:14, 17, 18, 22, 25, 28, 29, 30, 31, 32, 33, 34)
3502  aot = 0.02
3503  case (-999)
3504  aot = 0.02
3505  case default
3506  aot = 0.02
3507  print *, "ERROR: Invalid geographical zone in handle_lut_out_of_bounds: ", geo_zone
3508  print *, "Use default aot of 0.02 rather than return without output"
3509  !status = -1
3510  !return
3511  end select
3512  end if
3513 
3514  return
3515 
3516  end function handle_lut_out_of_bounds
3517 
3518 
3519  subroutine smoke_mod(view,scat_ang,dd412,dd470,dd650,
3520  1 tau_x412ss,tau_x412ss_91,tau_smoke)
3521  implicit none
3522 
3523  real, intent(in) :: view
3524  real, intent(in) :: scat_ang
3525  real, intent(in) :: dd412
3526  real, intent(in) :: dd470
3527  real, intent(in) :: dd650
3528  real, intent(in) :: tau_x412ss, tau_x412ss_91
3529  real :: tau_smoke
3530 
3531  if (dd412 <1.0 .and. dd470 <1.0) then
3532  tau_smoke = tau_x412ss * 2.85
3533  if (dd470 > dd412 .and. dd650 > dd412)
3534  1 tau_smoke = tau_x412ss_91 * 1.77
3535  if (abs(dd412-dd470) < 0.25) then
3536  if (dd412 <0.35) tau_smoke = tau_x412ss * 4.3
3537  if (scat_ang <158.0 .and. dd412 <0.8) then
3538  tau_smoke = tau_x412ss * 3.5
3539  if (view >45.0 .and. dd470 <0.2 .and. abs(dd412-dd470) < 0.24)
3540  1 tau_smoke = tau_x412ss * 6.0
3541  endif
3542  endif
3543  if (abs(dd412-dd650) < 0.1 .and. dd412 >0.7)
3544  1 tau_smoke = tau_x412ss_91
3545  if (dd650 > 1.0 .and. abs(dd412-dd470) < 0.25) then
3546  tau_smoke = tau_x412ss * 4.5
3547  if (scat_ang <158.0 .and. dd412 >0.8)
3548  1 tau_smoke = tau_x412ss * 3.5
3549  endif
3550  endif
3551 
3552  if (dd412 >1.0 .and. dd470 <1.0) then
3553  tau_smoke = tau_x412ss * 2.7
3554  if (abs(dd412-dd470) < 0.25 .and. dd650 <0.)
3555  1 tau_smoke = tau_x412ss * 2.4
3556  endif
3557 
3558  if (dd412 <1.0 .and. dd470 >1.0) then
3559  tau_smoke = tau_x412ss * 2.7
3560  if (abs(dd412-dd470) < 0.25 .and. dd650 <0.)
3561  1 tau_smoke = tau_x412ss * 2.4
3562  if (dd650 - dd470 > -0.2) then
3563  tau_smoke = tau_x412ss * 2.1
3564  if (view < 45.0) tau_smoke = tau_x412ss_91
3565  endif
3566  endif
3567 
3568  if (dd412 >1.0 .and. dd470 >1.0) then
3569  if (dd412 > dd470) tau_smoke = tau_x412ss_91 * 1.15
3570  if (dd470 > dd650 .and. dd650 > dd412 .and. dd650 >1.0) then
3571  if (view >= 45.0) tau_smoke = tau_x412ss * 2.1
3572  if (view < 45.0) tau_smoke = tau_x412ss_91
3573  endif
3574  if (dd412 <2.0 .and. dd470 >2.0)
3575  1 tau_smoke = tau_x412ss * 1.4
3576  if (dd650 > dd470 .and. dd470 > dd412) then
3577  if (view <= 25.0) tau_smoke = (tau_x412ss+tau_x412ss_91)/2.
3578  if (view > 25.0 .and. abs(dd412-dd470) < 0.2)
3579  1 tau_smoke = (tau_x412ss+tau_x412ss_91)/2.
3580  endif
3581  if (dd650 > 3.5) tau_smoke = tau_x412ss
3582  endif
3583 
3584  if (dd650 <-1.0 .and. abs(dd412-dd650)>1.9 .and. view >45.) then
3585  if (dd412 >= 0.88) tau_smoke = tau_x412ss * 4.0
3586  if (dd412 <0.88 .and. dd412 >= 0.6) tau_smoke = tau_x412ss *6.0
3587  if (dd412 <0.6 .and. dd412 >= 0.35) tau_smoke = tau_x412ss *8.0
3588  if (dd412 <0.35) tau_smoke = tau_x412ss_91 *8.0
3589  if (abs(dd412-dd470) < 0.14 .and. dd650 > -1.7 .and. dd412 <0.4)
3590  1 tau_smoke = tau_x412ss * 3.5
3591  if (dd412 <0.) tau_smoke = tau_x412ss_91 *9.0
3592  endif
3593 
3594  return
3595 
3596  end subroutine smoke_mod
3597 
dust
subroutine dust
Definition: 6sm1.f:3800
index
an array had not been initialized Several spelling and grammar corrections were which is read from the appropriate MCF the above metadata values were hard coded A problem calculating the average background DN for SWIR bands when the moon is in the space view port was corrected The new algorithm used to calculate the average background DN for all reflective bands when the moon is in the space view port is now the same as the algorithm employed by the thermal bands For non SWIR changes in the averages are typically less than Also for non SWIR the black body DNs remain a backup in case the SV DNs are not available For SWIR the changes in computed averages were larger because the old which used the black body suffered from contamination by the micron leak As a consequence of the if SV DNs are not available for the SWIR the EV pixels will not be the granule time is used to identify the appropriate tables within the set given for one LUT the first two or last two tables respectively will be used for the interpolation If there is only one LUT in the set of it will be treated as a constant LUT The manner in which Earth View data is checked for saturation was changed Previously the raw Earth View DNs and Space View DNs were checked against the lookup table values contained in the table dn_sat The change made is to check the raw Earth and Space View DNs to be sure they are less than the maximum saturation value and to check the Space View subtracted Earth View dns against a set of values contained in the new lookup table dn_sat_ev The metadata configuration and ASSOCIATEDINSTRUMENTSHORTNAME from the MOD02HKM product The same metatdata with extensions and were removed from the MOD021KM and MOD02OBC products ASSOCIATEDSENSORSHORTNAME was set to MODIS in all products These changes are reflected in new File Specification which users may consult for exact the pow functions were eliminated in Emissive_Cal and Emissive bands replaced by more efficient code Other calculations throughout the code were also made more efficient Aside from a few round off there was no difference to the product The CPU time decreased by about for a day case and for a night case A minor bug in calculating the uncertainty index for emissive bands was corrected The frame index(0-based) was previously being used the frame number(1-based) should have been used. There were only a few minor changes to the uncertainty index(maximum of 1 digit). 3. Some inefficient arrays(Sigma_RVS_norm_sq) were eliminated and some code lines in Preprocess_L1A_Data were moved into Process_OBCEng_Emiss. There were no changes to the product. Required RAM was reduced by 20 MB. Now
code
===========================================================================V4.1.3 12/18/2002============================================================================Changes which do not affect scientific output:1. The R *LUT was eliminated and the equivalent formulation for R *, i.e. 1/(m1 *e_sun_over_pi), was substituted for it in the only instance of its use, which is in the calculation of the RSB uncertainty index. This reduces the size of the Reflective LUT HDF file by approximately 1/4 to 1/3. The equivalent formulation of R *differed from the new by at most 0.056% in test granules and uncertainty differences of at most 1 count(out of a range of 0-15) were found in no more than 1 in 100, 000 pixels. 2. In Preprocess.c, a small error where the trailing dropped scan counter was incremented when the leading dropped scan counter should have been was fixed. This counter is internal only and is not yet used for any purpose. 3. NEW MYD02OBC Metadata Configuration Files. MCST wishes to have the OBC files archived even when the Orbit Number is recorded as "-1". Accordingly, ECS has delivered new MCF files for OBC output having all elements in the OrbitCalculatedSpatialDomain container set to "MANDATORY=FALSE". 4. pgs_in.version is now reset to "1" in Metadata.c before the call to look up the geolocation gringpoint information.============================================================================V4.1.1 CODE SPECIFIC TO MODIS/AQUA(FM1) 10/03/2002============================================================================Two changes were made to the code which do not affect scientific output:1. A bug which caused PGE02 to fail when scans were dropped between granules was fixed.(The length of the error message generated was shortened.) 2. Messages regarding an invalid MCST LUT Version or an invalid Write High Resolution Night Mode Output value in the PCF file were added.==============================================================================V4.1.0 CODE SPECIFIC TO MODIS/AQUA(FM1)(NEVER USED IN PRODUCTION) 07/30/2002==============================================================================Changes which impact scientific output of code:1. The LUT type of the RVS corrections was changed to piecewise linear. In addition the RVS LUTs were changed from listing the RVS corrections to listing the quadratic coefficients necessary to make the RVS corrections. The coefficients are now calculated by interpolating on the granule collection time and the RVS corrections are then generated using the interpolated coefficients. Previously used Emissive and Reflective RVS LUT tables were eliminated and new ones introduced. Several changes were made to the code which should not affect scientific output. They are:1. The ADC correction algorithm and related LUTs were stripped from the code.(The ADC correction has always been set to "0" so this has no scientific impact.) 2. Some small changes to the code, chiefly to casting of variables, were added to make it LINUX-compatible. Output of code run on LINUX machines displays differences of at most 1 scaled integer(SI) from output of code run on IRIX machines. The data type of the LUT "dn_sat_ev" was changed to float64 to avoid discrepancies seen between MOD_PR02 run on LINUX systems and IRIX systems where values were flagged under one operating system but not the other. 3. Checking for non-functioning detectors, sector rotation, incalculable values of the Emissive calibration factor "b1", and incalculable values of SV or BB averages was moved outside the loop over frames in Emissive_Cal.c since none of these quantities are frame-dependent. 4. The code was altered so that if up to five scans are dropped between the leading/middle or middle/trailing granules, the leading or trailing granule will still be used in emissive calibration to form a cross-granule average. QA bits 25 and 26 are set for a gap between the leading/middle and middle/trailing granules respectively. This may in rare instances lead to a change in emissive calibration coefficients for scans at the beginning or end of a granule. 5.(MODIS/AQUA ONLY) The name of the seed(error message) file was changed from "MODIS_36100.h" to "MODIS_36110.h". 6. Metadata.c was changed so that the source of the geolocation metadata is the input geolocation file rather than the L1A granule. 7. To reduce to overall size of the reflective LUT HDF files, fill values were eliminated from all LUTs previously dimensioned "BDSM"([NUM_REFLECTIVE_BANDS] *[MAX_DETECTORS_PER_BAND] *[MAX_SAMPLES_PER_BAND] *[NUM_MIRROR_SIDES]) in the LUT HDF files. Each table piece is stored in the HDF file with dimensions NUM_REFLECTIVE_INDICES, where NUM_REFLECTIVE_INDICES=[NUM_250M_BANDS *DETECTORS_PER_250M_BAND *SAMPLES_PER_250M_BAND *NUM_MIRROR_SIDES]+[NUM_500M_BANDS *DETECTORS_PER_500M_BAND *SAMPLES_PER_500M_BAND *NUM_MIRROR_SIDES]+[NUM_1000M_BANDS *DETECTORS_PER_1KM_BAND *SAMPLES_PER_1KM_BAND *NUM_MIRROR_SIDES] with SAMPLES_PER_250M_BAND=4, SAMPLES_PER_500M_BAND=2, and SAMPLES_PER_1KM_BAND=1. Values within each table piece appear in the order listed above. The overall dimensions of time dependent BDSM LUTs are now[NUM_TIMES] *[NUM_REFLECTIVE_INDICES], where NUM_TIMES is the number of time dependent table pieces. 8. Checking for non-functioning detectors, sector rotation, incalculable values of the Emissive calibration factor "b1", and incalculable values of SV or BB averages was moved outside the loop over frames in Emissive_Cal.c since none of these quantities are frame-dependent. 9. The code was altered so that if up to five scans are dropped between the leading/middle or middle/trailing granules, the leading or trailing granule will still be used in emissive calibration to form a cross-granule average. QA bits 25 and 26 are set for a gap between the leading/middle and middle/trailing granules respectively. This may in rare instances lead to a change in emissive calibration coefficients for scans at the beginning or end of a granule. 10. The array of b1s in Preprocess.c was being initialized to -1 outside the loop over bands, which meant that if b1 could not be computed, the value of b1 from the previous band for that scan/detector combination was used. The initialization was moved inside the band loop.============================================================================V3.1.0(Original Aqua-specific code version) 02/06/2002============================================================================AQUA-Specific changes made:1. A correction to a problem with blackbody warmup on bands 33, 35, and 36 was inserted. PC Bands 33, 35, and 36 on MODIS Aqua saturate on BB warmup before 310K, which means current code will not provide correct b1 calibration coefficients when the BB temperatures are above the saturation threshold. A LUT with default b1s and band-dependent temperature thresholds will be inserted in code. If the BB temperature is over the saturation threshold for the band, the default b1 from the table is used. 2. The number of possible wavelengths in the Emissive LUT RSR file was changed to 67 in order to accommodate the Aqua RSR tables. 3. Several changes to the upper and lower bound limits on LUT values were inserted. Changes to both Aqua and Terra Code:1. A check was put into Emissive_Cal.c to see whether the value of b1 being used to calibrate a pixel is negative. If so, the pixel is flagged with the newly created flag TEB_B1_NOT_CALCULATED, value 65526, and the number of pixels for which this occurs is counted in the QA_common table. 2. The array of b1s in Preprocess.c was being initialized to -1 outside the loop over bands, which meant that if b1 could not be computed, the value of b1 from the previous band for that scan/detector combination was used. The initialization was moved inside the band loop. 3. Minor code changes were made to eliminate compiler warnings when the code is compiled in 64-bit mode. 4. Temperature equations were upgraded to be MODIS/Aqua or MODIS/Terra specific and temperature conversion coefficients for Aqua were inserted.========================================================================================================================================================ALL CHANGES BELOW ARE TO COMMON TERRA/AQUA CODE USED BEFORE 02/06/2002========================================================================================================================================================v3.0.1 11/26/2001============================================================================Several small changes to the code were made, none of which changes the scientific output:1. The code was changed so that production of 250m and 500m resolution data when all scans of a granule are in night mode may be turned off/on through the PCF file. 2. A check on the times of the leading and trailing granules was inserted. If a leading or trailing granule does not immediately precede or follow(respectively) the middle granule, it is treated as a missing granule and a warning message is printed. 3. The code now reads the "MCST Version Number"(e.g. "3.0.1.0_Terra") from the PCF file and checks it against the MCST Version number contained in the LUT HDF files. This was done to allow the user to make sure the code is being run using the correct LUT files.(The designators "0_Terra", "1_Terra", etc.) refer to the LUT versions.) 4. A small bug in Preprocess.c was corrected code
Definition: HISTORY.txt:661
table
an array had not been initialized Several spelling and grammar corrections were which is read from the appropriate MCF the above metadata values were hard coded A problem calculating the average background DN for SWIR bands when the moon is in the space view port was corrected The new algorithm used to calculate the average background DN for all reflective bands when the moon is in the space view port is now the same as the algorithm employed by the thermal bands For non SWIR changes in the averages are typically less than Also for non SWIR the black body DNs remain a backup in case the SV DNs are not available For SWIR the changes in computed averages were larger because the old which used the black body suffered from contamination by the micron leak As a consequence of the if SV DNs are not available for the SWIR the EV pixels will not be the granule time is used to identify the appropriate tables within the set given for one LUT table(e.g. m1) to be used for interpolation. The table values will be linearly interpolated using the tables corresponding to the node times bracketing the granule time. If the granule time falls before the time of the first node or after the time of the last node
landcover::get_landcover
integer function, public get_landcover(lat, lon, status)
Definition: landcover_nc4.f95:100
find_v_viirs
subroutine find_v_viirs(realbuf, tmpvg, outbuf, qa_flag, px_elev, lm, windsp, wv)
Definition: find_v_viirs_nc4.f:10
rayleigh
void rayleigh(l1str *l1rec, int32_t ip)
Definition: rayleigh.c:130
real
#define real
Definition: DbAlgOcean.cpp:26
calendars::season_from_doy
integer function, public season_from_doy(yr, doy)
Definition: calendars.f95:539
find_v_veg
subroutine find_v_veg(month, season, realbuf, tmpvg, r412sv, r470sv, gzflg, outbufvg, tau_x470_flag)
Definition: find_v_vegset_nc4.f:3
calendars::gregorian_from_doy
type(gdatetime) function, public gregorian_from_doy(yr, doy)
Definition: calendars.f95:469
viirs_aerosol_luts::aero_412
subroutine, public aero_412(dflag, refl, x1, x2, x3, mm, nn, ll, ma, imod, r412, tau_x412, tau_x412_flag, trflg, model_frac, debug)
Definition: viirs_aerosol_luts_nc4.f95:1205
usage
void usage(char *progname)
Definition: main_biosmap.c:381
viirs_aerosol_luts
Definition: viirs_aerosol_luts_nc4.f95:1
search
subroutine search(dflag, xbar, x, n, i)
Definition: find_v_viirs_nc4.f:3401
search2
subroutine search2(dflag, xbar, x, n, i, fx)
Definition: find_v_viirs_nc4.f:3452
load
void load(float x1, float v[], float y[])
max
#define max(A, B)
Definition: main_biosmap.c:61
and
===========================================================================V5.0.48(Terra) 03/20/2015 Changes shown below are differences from MOD_PR02 V5.0.46(Terra)============================================================================Changes noted for V6.1.20(Terra) below were also instituted for this version.============================================================================V6.1.20(Terra) 03/12/2015 Changes shown below are differences from MOD_PR02 V6.1.18(Terra)============================================================================Changes from v6.1.18 which may affect scientific output:A situation can occur in which a scan which contains sector rotated data has a telemetry value indicating the completeness of the sector rotation. This issue is caused by the timing of the instrument command to perform the sector rotation and the recording of the telemetry point that reports the status of sector rotation. In this case a scan is considered valid by L1B and pass through the calibration - reporting extremely high radiances. Operationally the TEB calibration uses a 40 scan average coefficient, so the 20 scans(one mirror side) after the sector rotation are contaminated with anomalously high radiance values. A similar timing issue appeared before the sector rotation was fixed in V6.1.2. Our analysis indicates the ‘SET_FR_ENC_DELTA’ telemetry correlates well with the sector rotation encoder position. The use of this telemetry point to determine scans that are sector rotated should fix the anomaly occured before and after the sector rotation(usually due to the lunar roll maneuver). The fix related to the sector rotation in V6.1.2 is removed in this version.============================================================================V6.1.18(Terra) 10/01/2014 Changes shown below are differences from MOD_PR02 V6.1.16(Terra)============================================================================Added doi attributes to NRT(Near-Real-Time) product.============================================================================V6.1.16(Terra) 01/27/2014 Changes shown below are differences from MOD_PR02 V6.1.14(Terra)============================================================================Migrate to SDP Toolkit 5.2.17============================================================================V6.1.14(Terra) 06/26/2012 Changes shown below are differences from MOD_PR02 V6.1.12(Terra)============================================================================Added the doi metadata to L1B product============================================================================V6.1.12(Terra) 04/25/2011 Changes shown below are differences from MOD_PR02 V6.1.8(Terra)============================================================================1. The algorithm to calculate uncertainties for reflective solar bands(RSB) is updated. The current uncertainty in L1B code includes 9 terms from prelaunch analysis. The new algorithm regroups them with the new added contributions into 5 terms:u1:the common term(AOI and time independent) and
Definition: HISTORY.txt:126
newsfc412_arab
subroutine newsfc412_arab(dflag2, mod_sfc, xday, xlatp, xlonp, xthet, xphi, scat_ang, terrain_flag_new5, r412_135, r412new)
Definition: pack_412.f:226
smoke_mod
subroutine smoke_mod(view, scat_ang, dd412, dd470, dd650,
Definition: find_v_viirs_nc4.f:3520
modis_surface
Definition: modis_surface_nc4.f95:1
set
set(CMAKE_CXX_FLAGS "${CMAKE_CXX_FLAGS} -DWITH_MPI") target_link_libraries(afrt_nc4 $
Definition: CMakeLists.txt:16
min
#define min(A, B)
Definition: main_biosmap.c:62
locate
subroutine locate(xx, n, x, j)
Definition: atrem_app_refl_f90_cubeio.f:3114
eq
#define eq(a, b)
Definition: rice.h:167
modis_surface::get_brdfcorr_sr
integer function, public get_brdfcorr_sr(lat, lon, ra, sa, vza, amf, elev, month, ndvi, stdv, gzone, lc_type, bgaod, sr412, sr470, sr650, use_alternate_brdf, debug)
Definition: modis_surface_nc4.f95:749
aero_mod
subroutine aero_mod(tau_x412, tau_x470, tau_x650,
Definition: find_v_viirs_nc4.f:3475
copy
void copy(double **aout, double **ain, int n)
handle_lut_out_of_bounds
integer function handle_lut_out_of_bounds(geo_zone, retr_flag, aot)
Definition: find_v_viirs_nc4.f:3489
for
for(i=0;i< NROOTS;i++) s[i]
Definition: decode_rs.h:85
abs
#define abs(a)
Definition: misc.h:90
on
this program makes no use of any feature of the SDP Toolkit that could generate such a then geolocation is calculated at that and then aggregated up to Resolved feature request Bug by adding three new int8 SDSs for each high resolution offsets between the high resolution geolocation and a bi linear interpolation extrapolation of the positions This can be used to reconstruct the high resolution geolocation Resolved Bug by delaying cumulation of gflags until after validation of derived products Resolved Bug by setting Latitude and Longitude to the correct fill resolving to support Near Real Time because they may be unnecessary if use of entrained ephemeris and attitude data is turned on(as it will be in Near-Real-Time processing).
direct
subroutine direct(N, M, N1, M1, NN, MM, C)
Definition: tmd.lp.f:1299


ABOUT

Our Team

Contact Us

DATA

Get Started

Find Data

Use Data

TOOLS

SeaBASS

SeaDAS

Giovanni Data Visualizer

Worldview Mapping Interface

Level 3 Time Series Plotter

Request a SeaHawk Image

Mission Quality Monitor

HELP

Level 1 & 2 Browser Help

Level 3 & 4 Browser Help

File Search Help

Search & Download Help

Earthdata Forum

Join Mailing List

How to Cite Data

CONNECT WITH US
COMPLIANCE

Web Privacy Policy

Data & Information Policy

Communications Policy

Accessibility Policy

Freedom of Information Act

USA.gov

Responsible NASA Official: Sean Bailey
Curator: OceanColor Webmaster