module calc_fov_crosstrk !$$$ module documentation block ! . . . . ! module: calc_fov_crosstrk ! ! abstract: contains all routines necessary for fov calculations ! for cross track scanning instruments ! ! program history log: ! 2007-08-13 kleespies/gayno initial version ! 2009-03-05 Kleespies AMSU antenna pattern now channel ! specific. add channel specific ! MHS antenna patterns. add airs ! and iasi. ! 2011-09-13 gayno improve module's error handling ! 2015-03-24 ejones changed maxinstr to accomodate the addition of ! megha tropiques ! 2016-11-11 gayno add ATMS 3.3 ! ! subroutines included: ! sub fov_ellipse_crosstrk - calc lat/lon of fov polygon ! sub fovanglessizes - compute cross/along track angles/fov sizes ! sub get_sat_height - get satellite height ! sub fov_cleanup - free up memory ! sub fov_check - ensure fov number is valid ! sub instrument_init - initialize instrument fields ! sub inside_fov_crosstrk - determine if point inside fov ! ! Variable Definitions: ! def alongtrackangle - along track angle ! def alongtrackfovsize - along track fov size ! def amsucoeff - amsu-a antenna power coefficients ! for a specific satellite. ! def amsucoeff_XX - amsu-a antenna power coefficients ! for satellite number XX. ! def crosstrackangle - cross track angle ! def crosstrackfovsize - cross track fov size ! def eccen - fov eccentricity ! def maxinstr - maximum number of instruments ! def maxfov - maximum number of fovs for each instrument ! def mhscoeff - mhs antenna power coefficients ! for a specific satellite. ! def mhscoeff_XX - mhs antenna power coefficients ! for satellite number XX. ! def npoly - number of vertices in the fov polygon ! def psi - angle of each vertex of fov polygon ! def rmax - major axis of ellipse ! def rmin - minor axis of ellipse ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ end documentation block use kinds, only : i_kind, r_kind implicit none private integer(i_kind) , parameter, public :: npoly = 30 integer(i_kind) , parameter, private :: maxinstr = 20 integer(i_kind) , dimension(maxinstr), private:: maxfov = (/ 2048,2048,2048, & 56,56,56, & 56,56, 8, & 11,30,90, & 90,96,96, & 96,90,30, & 130,96 /) real(r_kind) , dimension(:), allocatable, private :: alongtrackangle real(r_kind) , dimension(:), allocatable, private :: crosstrackangle real(r_kind) , dimension(:), allocatable, private :: alongtrackfovsize real(r_kind) , dimension(:), allocatable, private :: crosstrackfovsize real(r_kind) , dimension(:), allocatable, private :: rmax real(r_kind) , dimension(:), allocatable, private :: rmin real(r_kind) , dimension(:), allocatable, private :: eccen real(r_kind) , private :: psi(npoly) public get_sat_height public instrument_init public fov_check public fov_cleanup public inside_fov_crosstrk public fov_ellipse_crosstrk real(r_kind), pointer, private, dimension(:,:,:) :: amsucoeff real(r_kind), pointer, private, dimension(:,:,:) :: atmscoeff ! Variables for AMSUA power determination. the number after the underscore is the ! satellite number. real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_15 = reshape( (/ & -1.00367358e-002_r_kind, -1.17344270e-002_r_kind, -9.72701728e-001_r_kind, -7.81457871e-003_r_kind, & -2.11871578e-003_r_kind, 6.00348110e-004_r_kind, -6.77819597e-004_r_kind, -2.68106451e-005_r_kind, & -2.42081527e-002_r_kind, 2.95991510e-001_r_kind, -9.96972978e-001_r_kind, -1.81507263e-002_r_kind, & 3.29453568e-003_r_kind, 3.66205233e-003_r_kind, -9.45091131e-004_r_kind, -1.26842526e-004_r_kind, & -6.65183738e-002_r_kind, -1.89245090e-001_r_kind, -8.38231504e-001_r_kind, -7.44134281e-003_r_kind, & -2.76436768e-002_r_kind, -1.69177575e-003_r_kind, 9.04653862e-004_r_kind, 6.48917630e-005_r_kind, & -9.20323208e-002_r_kind, -9.59552359e-003_r_kind, -8.66163492e-001_r_kind, -4.93938252e-002_r_kind, & -2.02207621e-002_r_kind, 2.83769006e-003_r_kind, 3.93795577e-004_r_kind, -6.02218970e-005_r_kind, & 4.64320648e-003_r_kind, 2.47414529e-001_r_kind, -8.36344659e-001_r_kind, -1.23036997e-002_r_kind, & -1.09757520e-002_r_kind, 1.78341044e-003_r_kind, 4.46469421e-005_r_kind, -5.94182457e-005_r_kind, & 4.34463955e-002_r_kind, 1.12749077e-001_r_kind, -9.76884484e-001_r_kind, -4.08062451e-002_r_kind, & 7.08068081e-004_r_kind, 5.47452038e-003_r_kind, -3.97553435e-004_r_kind, -1.02515020e-004_r_kind, & 3.91776152e-002_r_kind, -9.75272208e-002_r_kind, -9.23579991e-001_r_kind, -1.84527859e-002_r_kind, & -2.79494980e-003_r_kind, 4.74081701e-003_r_kind, -3.26868758e-005_r_kind, -6.42077284e-005_r_kind, & -1.61407292e-002_r_kind, -3.28314036e-001_r_kind, -8.62950087e-001_r_kind, 3.28221768e-002_r_kind, & -5.26583334e-003_r_kind, 7.51853397e-004_r_kind, -1.13260481e-004_r_kind, -2.68223033e-007_r_kind, & -1.07014296e-003_r_kind, 3.36564213e-001_r_kind, -8.70783627e-001_r_kind, -1.84429635e-003_r_kind, & -5.49669797e-003_r_kind, 3.81773541e-004_r_kind, -7.79173351e-005_r_kind, -8.48751370e-006_r_kind, & -3.09105590e-002_r_kind, 5.92924535e-001_r_kind, -9.27605867e-001_r_kind, 3.00982548e-003_r_kind, & -2.47880863e-003_r_kind, 1.88208383e-003_r_kind, -4.73521737e-004_r_kind, 4.02910246e-006_r_kind, & 1.13214113e-001_r_kind, 7.23310336e-002_r_kind, -9.99909580e-001_r_kind, 3.64811085e-002_r_kind, & 1.19802766e-002_r_kind, 4.69912251e-004_r_kind, -8.19198438e-004_r_kind, 6.82059472e-005_r_kind, & 3.13684642e-002_r_kind, -2.45412495e-002_r_kind, -8.36684287e-001_r_kind, 8.01033247e-003_r_kind, & -6.20706240e-003_r_kind, 1.85168872e-003_r_kind, -6.99376760e-005_r_kind, -1.83860229e-005_r_kind, & 3.33500355e-002_r_kind, 2.10814878e-001_r_kind, -8.89469445e-001_r_kind, -1.38472561e-002_r_kind, & -3.73062352e-003_r_kind, 2.76804063e-003_r_kind, -1.62501543e-004_r_kind, -5.73568032e-005_r_kind, & 4.07360271e-002_r_kind, -4.81393337e-002_r_kind, -9.35139358e-001_r_kind, -2.37985887e-002_r_kind, & -9.27939080e-004_r_kind, 3.88786383e-003_r_kind, -3.52895237e-004_r_kind, -8.52064040e-005_r_kind, & 3.51734906e-002_r_kind, -7.54902139e-002_r_kind, -9.31113124e-001_r_kind, -1.53986933e-002_r_kind, & -1.36302866e-003_r_kind, 4.01552487e-003_r_kind, -1.58052790e-004_r_kind, -7.67155652e-005_r_kind, & -1.59598049e-002_r_kind, -2.72733510e-001_r_kind, -9.45698082e-001_r_kind, 5.62173650e-002_r_kind, & -9.21276514e-004_r_kind, -6.84147817e-005_r_kind, -1.42305653e-004_r_kind, -8.02842897e-006_r_kind, & -3.72697525e-002_r_kind, 7.18167961e-001_r_kind, -9.05988336e-001_r_kind, -2.88108010e-002_r_kind, & -9.10314266e-003_r_kind, 6.75110379e-003_r_kind, -4.47749808e-005_r_kind, -1.43583864e-004_r_kind, & -1.57145746e-002_r_kind, 1.89278293e-002_r_kind, -8.31693709e-001_r_kind, 4.78307344e-002_r_kind, & -1.14242034e-002_r_kind, -4.96163173e-003_r_kind, -5.86726936e-004_r_kind, 1.64914774e-004_r_kind, & -3.72697525e-002_r_kind, 7.18167961e-001_r_kind, -9.05988336e-001_r_kind, -2.88108010e-002_r_kind, & -9.10314266e-003_r_kind, 6.75110379e-003_r_kind, -4.47749808e-005_r_kind, -1.43583864e-004_r_kind, & -1.57145746e-002_r_kind, 1.89278293e-002_r_kind, -8.31693709e-001_r_kind, 4.78307344e-002_r_kind, & -1.14242034e-002_r_kind, -4.96163173e-003_r_kind, -5.86726936e-004_r_kind, 1.64914774e-004_r_kind, & -3.72697525e-002_r_kind, 7.18167961e-001_r_kind, -9.05988336e-001_r_kind, -2.88108010e-002_r_kind, & -9.10314266e-003_r_kind, 6.75110379e-003_r_kind, -4.47749808e-005_r_kind, -1.43583864e-004_r_kind, & -1.57145746e-002_r_kind, 1.89278293e-002_r_kind, -8.31693709e-001_r_kind, 4.78307344e-002_r_kind, & -1.14242034e-002_r_kind, -4.96163173e-003_r_kind, -5.86726936e-004_r_kind, 1.64914774e-004_r_kind, & -3.72697525e-002_r_kind, 7.18167961e-001_r_kind, -9.05988336e-001_r_kind, -2.88108010e-002_r_kind, & -9.10314266e-003_r_kind, 6.75110379e-003_r_kind, -4.47749808e-005_r_kind, -1.43583864e-004_r_kind, & -1.57145746e-002_r_kind, 1.89278293e-002_r_kind, -8.31693709e-001_r_kind, 4.78307344e-002_r_kind, & -1.14242034e-002_r_kind, -4.96163173e-003_r_kind, -5.86726936e-004_r_kind, 1.64914774e-004_r_kind, & -3.72697525e-002_r_kind, 7.18167961e-001_r_kind, -9.05988336e-001_r_kind, -2.88108010e-002_r_kind, & -9.10314266e-003_r_kind, 6.75110379e-003_r_kind, -4.47749808e-005_r_kind, -1.43583864e-004_r_kind, & -1.57145746e-002_r_kind, 1.89278293e-002_r_kind, -8.31693709e-001_r_kind, 4.78307344e-002_r_kind, & -1.14242034e-002_r_kind, -4.96163173e-003_r_kind, -5.86726936e-004_r_kind, 1.64914774e-004_r_kind, & -3.72697525e-002_r_kind, 7.18167961e-001_r_kind, -9.05988336e-001_r_kind, -2.88108010e-002_r_kind, & -9.10314266e-003_r_kind, 6.75110379e-003_r_kind, -4.47749808e-005_r_kind, -1.43583864e-004_r_kind, & -1.57145746e-002_r_kind, 1.89278293e-002_r_kind, -8.31693709e-001_r_kind, 4.78307344e-002_r_kind, & -1.14242034e-002_r_kind, -4.96163173e-003_r_kind, -5.86726936e-004_r_kind, 1.64914774e-004_r_kind, & 3.87737527e-002_r_kind, -2.41649374e-001_r_kind, -9.37661827e-001_r_kind, 1.40356580e-002_r_kind, & 1.63714646e-003_r_kind, 1.03963783e-003_r_kind, -4.49913321e-004_r_kind, 2.07235353e-005_r_kind, & -1.22435763e-001_r_kind, -3.87322575e-001_r_kind, -4.49794233e-001_r_kind, 7.79016390e-002_r_kind, & -6.47236556e-002_r_kind, 4.48352861e-004_r_kind, 1.75243057e-003_r_kind, -7.81173076e-005_r_kind/), & (/8,2,15/) ) real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_16 = reshape( (/ & 2.36253720e-002_r_kind, -2.01956004e-001_r_kind, -9.42441463e-001_r_kind, 4.67489986e-003_r_kind, & -5.80562279e-003_r_kind, -6.12580625e-004_r_kind, -4.41984390e-004_r_kind, -3.52604752e-006_r_kind, & 2.19323412e-002_r_kind, -2.41464257e-001_r_kind, -8.76489699e-001_r_kind, 7.69139780e-003_r_kind, & -9.98940202e-004_r_kind, -6.81366888e-004_r_kind, -7.79129798e-004_r_kind, -9.83764585e-006_r_kind, & -1.03580669e-001_r_kind, -4.66340750e-001_r_kind, -8.48327339e-001_r_kind, 1.64177082e-002_r_kind, & -3.32605094e-002_r_kind, -2.10031378e-003_r_kind, 1.23471778e-003_r_kind, 7.77417736e-005_r_kind, & -1.74308643e-001_r_kind, -4.48501498e-001_r_kind, -7.84894824e-001_r_kind, -2.24881396e-002_r_kind, & -3.02584395e-002_r_kind, 3.28395003e-003_r_kind, 8.63748952e-004_r_kind, -8.32117294e-005_r_kind, & 1.41929025e-002_r_kind, -3.19731236e-001_r_kind, -8.94576848e-001_r_kind, -8.88386834e-003_r_kind, & -5.62882330e-003_r_kind, 1.30801124e-003_r_kind, -2.34643157e-004_r_kind, -3.08262061e-005_r_kind, & 3.66718285e-002_r_kind, 3.11904073e-001_r_kind, -9.76906300e-001_r_kind, -7.42634246e-003_r_kind, & -1.05722966e-002_r_kind, 2.43136333e-003_r_kind, 6.73694740e-006_r_kind, -5.72227436e-005_r_kind, & 2.43993644e-002_r_kind, 2.97852010e-001_r_kind, -8.69731128e-001_r_kind, 5.95848588e-003_r_kind, & -1.40882675e-002_r_kind, 1.54716417e-003_r_kind, 1.30630637e-004_r_kind, -1.47408900e-005_r_kind, & 1.95409320e-002_r_kind, 3.06262583e-001_r_kind, -9.00652170e-001_r_kind, -2.64048314e-004_r_kind, & -2.35236697e-002_r_kind, 4.14541783e-003_r_kind, 7.66475452e-004_r_kind, -1.78228700e-004_r_kind, & 3.27163115e-002_r_kind, 1.54714525e-001_r_kind, -9.09946263e-001_r_kind, -1.01834759e-002_r_kind, & -1.32468473e-002_r_kind, 2.72665126e-003_r_kind, 9.75327930e-005_r_kind, -8.41464716e-005_r_kind, & -9.23690666e-003_r_kind, 2.90396154e-001_r_kind, -9.56662178e-001_r_kind, -1.72896422e-002_r_kind, & -8.10661633e-003_r_kind, 3.93800158e-003_r_kind, -1.26169616e-004_r_kind, -1.14084214e-004_r_kind, & -3.62905040e-002_r_kind, 4.50896472e-001_r_kind, -9.21234190e-001_r_kind, 4.01789322e-004_r_kind, & -1.30359689e-002_r_kind, 3.19115492e-003_r_kind, 1.62691125e-004_r_kind, -5.46014235e-005_r_kind, & 7.36603187e-003_r_kind, 2.78184146e-001_r_kind, -8.94912243e-001_r_kind, 7.14887865e-003_r_kind, & -1.17667923e-002_r_kind, 3.00864619e-003_r_kind, 1.51504573e-004_r_kind, -8.34179082e-005_r_kind, & -6.21941835e-002_r_kind, 3.43763500e-001_r_kind, -9.09502566e-001_r_kind, -2.85539613e-003_r_kind, & -1.45036662e-002_r_kind, 2.73107341e-003_r_kind, 1.15043804e-004_r_kind, -1.07309730e-004_r_kind, & 4.28958330e-004_r_kind, 3.40785056e-001_r_kind, -1.04653192e+000_r_kind, -1.01089070e-003_r_kind, & -4.18298785e-003_r_kind, 1.99506595e-003_r_kind, -1.34124552e-004_r_kind, -5.35102045e-005_r_kind, & 3.57502624e-002_r_kind, 5.50667197e-002_r_kind, -9.01254117e-001_r_kind, -1.54003268e-002_r_kind, & -1.58270784e-002_r_kind, 3.67956469e-003_r_kind, 3.11302574e-004_r_kind, -7.75082663e-005_r_kind, & -9.44588482e-002_r_kind, 7.50300825e-001_r_kind, -1.01519608e+000_r_kind, -1.62286758e-002_r_kind, & -2.35847179e-002_r_kind, 5.65459207e-003_r_kind, 9.53018840e-004_r_kind, -2.01222618e-004_r_kind, & 3.39001790e-002_r_kind, 2.41042480e-001_r_kind, -9.38472331e-001_r_kind, -9.81334690e-003_r_kind, & -1.13598527e-002_r_kind, 2.97825597e-003_r_kind, -3.46505985e-005_r_kind, -8.45731338e-005_r_kind, & 1.22634219e-002_r_kind, 1.28929421e-001_r_kind, -1.26483619e+000_r_kind, 1.04203574e-001_r_kind, & 1.66568719e-002_r_kind, -9.81434155e-003_r_kind, -8.17351334e-004_r_kind, 2.54426064e-004_r_kind, & 3.39001790e-002_r_kind, 2.41042480e-001_r_kind, -9.38472331e-001_r_kind, -9.81334690e-003_r_kind, & -1.13598527e-002_r_kind, 2.97825597e-003_r_kind, -3.46505985e-005_r_kind, -8.45731338e-005_r_kind, & 1.22634219e-002_r_kind, 1.28929421e-001_r_kind, -1.26483619e+000_r_kind, 1.04203574e-001_r_kind, & 1.66568719e-002_r_kind, -9.81434155e-003_r_kind, -8.17351334e-004_r_kind, 2.54426064e-004_r_kind, & 3.39001790e-002_r_kind, 2.41042480e-001_r_kind, -9.38472331e-001_r_kind, -9.81334690e-003_r_kind, & -1.13598527e-002_r_kind, 2.97825597e-003_r_kind, -3.46505985e-005_r_kind, -8.45731338e-005_r_kind, & 1.22634219e-002_r_kind, 1.28929421e-001_r_kind, -1.26483619e+000_r_kind, 1.04203574e-001_r_kind, & 1.66568719e-002_r_kind, -9.81434155e-003_r_kind, -8.17351334e-004_r_kind, 2.54426064e-004_r_kind, & 3.39001790e-002_r_kind, 2.41042480e-001_r_kind, -9.38472331e-001_r_kind, -9.81334690e-003_r_kind, & -1.13598527e-002_r_kind, 2.97825597e-003_r_kind, -3.46505985e-005_r_kind, -8.45731338e-005_r_kind, & 1.22634219e-002_r_kind, 1.28929421e-001_r_kind, -1.26483619e+000_r_kind, 1.04203574e-001_r_kind, & 1.66568719e-002_r_kind, -9.81434155e-003_r_kind, -8.17351334e-004_r_kind, 2.54426064e-004_r_kind, & 3.39001790e-002_r_kind, 2.41042480e-001_r_kind, -9.38472331e-001_r_kind, -9.81334690e-003_r_kind, & -1.13598527e-002_r_kind, 2.97825597e-003_r_kind, -3.46505985e-005_r_kind, -8.45731338e-005_r_kind, & 1.22634219e-002_r_kind, 1.28929421e-001_r_kind, -1.26483619e+000_r_kind, 1.04203574e-001_r_kind, & 1.66568719e-002_r_kind, -9.81434155e-003_r_kind, -8.17351334e-004_r_kind, 2.54426064e-004_r_kind, & 3.39001790e-002_r_kind, 2.41042480e-001_r_kind, -9.38472331e-001_r_kind, -9.81334690e-003_r_kind, & -1.13598527e-002_r_kind, 2.97825597e-003_r_kind, -3.46505985e-005_r_kind, -8.45731338e-005_r_kind, & 1.22634219e-002_r_kind, 1.28929421e-001_r_kind, -1.26483619e+000_r_kind, 1.04203574e-001_r_kind, & 1.66568719e-002_r_kind, -9.81434155e-003_r_kind, -8.17351334e-004_r_kind, 2.54426064e-004_r_kind, & 2.49365754e-002_r_kind, 2.07554460e-001_r_kind, -9.28714156e-001_r_kind, -3.25016864e-003_r_kind, & -1.66032407e-002_r_kind, 3.08119995e-003_r_kind, 4.82223346e-004_r_kind, -8.70772346e-005_r_kind, & -8.07163268e-002_r_kind, 3.89334202e-001_r_kind, -1.21344244e+000_r_kind, 1.38754621e-001_r_kind, & -7.67704891e-003_r_kind, -9.02157184e-003_r_kind, 6.78973622e-004_r_kind, 1.44773308e-004_r_kind/),& (/ 8,2,15/) ) real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_17 = reshape( (/ & 1.04546631e-002_r_kind, -1.68463737e-002_r_kind, -9.45984364e-001_r_kind, 1.37877243e-003_r_kind, & -3.02611245e-003_r_kind, -4.58378316e-004_r_kind, -6.50028989e-004_r_kind, 1.92575317e-005_r_kind, & 8.65941588e-003_r_kind, -5.86456992e-002_r_kind, -1.05106378e+000_r_kind, -3.31778708e-003_r_kind, & 7.11675407e-003_r_kind, 3.46897374e-004_r_kind, -1.00040296e-003_r_kind, -2.32413208e-006_r_kind, & -1.14871182e-001_r_kind, -1.51089326e-001_r_kind, -7.85072088e-001_r_kind, 3.14582437e-002_r_kind, & -4.47488837e-002_r_kind, -3.65157402e-003_r_kind, 1.55843853e-003_r_kind, 9.86069936e-005_r_kind, & -1.04314178e-001_r_kind, -2.74469480e-002_r_kind, -8.81356478e-001_r_kind, -2.02818569e-002_r_kind, & -2.86703520e-002_r_kind, -5.17870300e-004_r_kind, 1.00237259e-003_r_kind, 4.57371680e-005_r_kind, & 2.61975098e-002_r_kind, 2.57663220e-001_r_kind, -9.20347750e-001_r_kind, 3.82701890e-003_r_kind, & -6.57390244e-003_r_kind, -4.66979400e-004_r_kind, -2.55515653e-004_r_kind, 2.64522496e-005_r_kind, & 2.15802751e-002_r_kind, -5.81635386e-002_r_kind, -9.84405994e-001_r_kind, 5.96869783e-003_r_kind, & -6.38028421e-003_r_kind, -7.33875029e-004_r_kind, -1.59733929e-004_r_kind, 2.81333669e-005_r_kind, & 3.78270298e-002_r_kind, 3.36511105e-001_r_kind, -8.78714681e-001_r_kind, 4.41203304e-002_r_kind, & -1.62802711e-002_r_kind, -3.45655321e-003_r_kind, 6.39992184e-004_r_kind, 5.08310477e-005_r_kind, & 1.01109396e-003_r_kind, 6.70918003e-002_r_kind, -9.33075070e-001_r_kind, -1.74348708e-002_r_kind, & -2.16634739e-002_r_kind, 5.47799282e-003_r_kind, 9.33322823e-004_r_kind, -1.77886817e-004_r_kind, & -2.73076892e-002_r_kind, 4.72141653e-001_r_kind, -9.53527510e-001_r_kind, -6.89434644e-004_r_kind, & -5.42717846e-003_r_kind, 5.59467240e-004_r_kind, -3.03966750e-004_r_kind, 1.11371974e-005_r_kind, & -3.98405595e-004_r_kind, 1.15465194e-001_r_kind, -9.67057765e-001_r_kind, 1.19083177e-003_r_kind, & -5.91184292e-003_r_kind, 6.99364347e-004_r_kind, -2.35179934e-004_r_kind, -1.14399882e-005_r_kind, & 3.54661047e-002_r_kind, 5.62102981e-002_r_kind, -9.27996159e-001_r_kind, 2.60837171e-002_r_kind, & -1.34692034e-002_r_kind, -3.63388471e-003_r_kind, 6.23947941e-004_r_kind, 8.62157394e-005_r_kind, & 1.91068146e-002_r_kind, 1.31960198e-001_r_kind, -8.89069080e-001_r_kind, -1.12417685e-002_r_kind, & -1.85574889e-002_r_kind, 2.10160110e-003_r_kind, 7.30395375e-004_r_kind, -6.73496033e-005_r_kind, & 1.24845589e-002_r_kind, 2.64344871e-001_r_kind, -9.50289309e-001_r_kind, -4.72221710e-003_r_kind, & -6.23426493e-003_r_kind, 1.04571111e-003_r_kind, -2.39053625e-004_r_kind, -2.11943225e-005_r_kind, & 4.29382026e-002_r_kind, 4.44974098e-003_r_kind, -1.08117700e+000_r_kind, 4.19260468e-003_r_kind, & 2.20504985e-003_r_kind, -3.38937796e-004_r_kind, -4.93370055e-004_r_kind, 7.60885268e-006_r_kind, & -1.83724221e-002_r_kind, 2.54897237e-001_r_kind, -8.93854260e-001_r_kind, 5.88120036e-002_r_kind, & -1.57443024e-002_r_kind, -5.48272207e-003_r_kind, 6.38363475e-004_r_kind, 1.04877683e-004_r_kind, & -7.62775689e-002_r_kind, 6.55387521e-001_r_kind, -1.09617865e+000_r_kind, 7.14722602e-003_r_kind, & -1.56472009e-002_r_kind, 8.87253773e-005_r_kind, 9.29544040e-004_r_kind, -4.59635412e-005_r_kind, & -1.86581500e-002_r_kind, 2.93479741e-001_r_kind, -9.95894194e-001_r_kind, 1.39501002e-002_r_kind, & -2.50316015e-003_r_kind, -1.03119342e-003_r_kind, -4.54510242e-004_r_kind, 6.13363445e-005_r_kind, & -2.11460143e-001_r_kind, 8.15824568e-001_r_kind, -8.98793340e-001_r_kind, 1.26764514e-002_r_kind, & -1.21458182e-002_r_kind, -9.57899145e-004_r_kind, -2.11658567e-004_r_kind, 9.96286144e-006_r_kind, & -1.86581500e-002_r_kind, 2.93479741e-001_r_kind, -9.95894194e-001_r_kind, 1.39501002e-002_r_kind, & -2.50316015e-003_r_kind, -1.03119342e-003_r_kind, -4.54510242e-004_r_kind, 6.13363445e-005_r_kind, & -2.11460143e-001_r_kind, 8.15824568e-001_r_kind, -8.98793340e-001_r_kind, 1.26764514e-002_r_kind, & -1.21458182e-002_r_kind, -9.57899145e-004_r_kind, -2.11658567e-004_r_kind, 9.96286144e-006_r_kind, & -1.86581500e-002_r_kind, 2.93479741e-001_r_kind, -9.95894194e-001_r_kind, 1.39501002e-002_r_kind, & -2.50316015e-003_r_kind, -1.03119342e-003_r_kind, -4.54510242e-004_r_kind, 6.13363445e-005_r_kind, & -2.11460143e-001_r_kind, 8.15824568e-001_r_kind, -8.98793340e-001_r_kind, 1.26764514e-002_r_kind, & -1.21458182e-002_r_kind, -9.57899145e-004_r_kind, -2.11658567e-004_r_kind, 9.96286144e-006_r_kind, & -1.86581500e-002_r_kind, 2.93479741e-001_r_kind, -9.95894194e-001_r_kind, 1.39501002e-002_r_kind, & -2.50316015e-003_r_kind, -1.03119342e-003_r_kind, -4.54510242e-004_r_kind, 6.13363445e-005_r_kind, & -2.11460143e-001_r_kind, 8.15824568e-001_r_kind, -8.98793340e-001_r_kind, 1.26764514e-002_r_kind, & -1.21458182e-002_r_kind, -9.57899145e-004_r_kind, -2.11658567e-004_r_kind, 9.96286144e-006_r_kind, & -1.86581500e-002_r_kind, 2.93479741e-001_r_kind, -9.95894194e-001_r_kind, 1.39501002e-002_r_kind, & -2.50316015e-003_r_kind, -1.03119342e-003_r_kind, -4.54510242e-004_r_kind, 6.13363445e-005_r_kind, & -2.11460143e-001_r_kind, 8.15824568e-001_r_kind, -8.98793340e-001_r_kind, 1.26764514e-002_r_kind, & -1.21458182e-002_r_kind, -9.57899145e-004_r_kind, -2.11658567e-004_r_kind, 9.96286144e-006_r_kind, & -1.86581500e-002_r_kind, 2.93479741e-001_r_kind, -9.95894194e-001_r_kind, 1.39501002e-002_r_kind, & -2.50316015e-003_r_kind, -1.03119342e-003_r_kind, -4.54510242e-004_r_kind, 6.13363445e-005_r_kind, & -2.11460143e-001_r_kind, 8.15824568e-001_r_kind, -8.98793340e-001_r_kind, 1.26764514e-002_r_kind, & -1.21458182e-002_r_kind, -9.57899145e-004_r_kind, -2.11658567e-004_r_kind, 9.96286144e-006_r_kind, & -4.41915691e-002_r_kind, 1.59778103e-001_r_kind, -9.23613429e-001_r_kind, 6.46106228e-002_r_kind, & -2.46498063e-002_r_kind, -7.18605192e-003_r_kind, 9.84023209e-004_r_kind, 1.64684461e-004_r_kind, & -8.41936469e-002_r_kind, 2.42628723e-001_r_kind, -9.72940922e-001_r_kind, -2.08250862e-002_r_kind, & -2.04343498e-002_r_kind, -3.75839561e-004_r_kind, 8.31458718e-004_r_kind, 2.66584175e-005_r_kind/),& (/ 8,2,15/) ) real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_18 = reshape( (/ & 9.22924746e-003_r_kind, 1.49358943e-001_r_kind, -1.00322044e+000_r_kind, -3.36418790e-003_r_kind, & 7.25480262e-003_r_kind, 5.94239740e-004_r_kind, -1.06518064e-003_r_kind, 8.64808044e-006_r_kind, & -1.35190440e-002_r_kind, -1.06490999e-001_r_kind, -1.07750916e+000_r_kind, -5.29923243e-003_r_kind, & 8.65704659e-003_r_kind, 8.76962658e-005_r_kind, -1.11807382e-003_r_kind, 2.69893872e-006_r_kind, & -1.26408741e-001_r_kind, 2.66547471e-001_r_kind, -8.51005554e-001_r_kind, -1.33056249e-002_r_kind, & -3.28817517e-002_r_kind, 4.60545672e-003_r_kind, 1.13587396e-003_r_kind, -1.36787363e-004_r_kind, & -9.95343253e-002_r_kind, 7.50881582e-002_r_kind, -8.77490580e-001_r_kind, 1.01857305e-001_r_kind, & -3.39203440e-002_r_kind, -2.77569343e-004_r_kind, 1.21702568e-003_r_kind, -1.26924133e-004_r_kind, & -1.08985245e-001_r_kind, -1.21106640e-001_r_kind, -9.80138302e-001_r_kind, -1.55268027e-003_r_kind, & -3.36335576e-003_r_kind, -4.89384867e-004_r_kind, -7.03390921e-004_r_kind, 1.62836750e-005_r_kind, & 2.32525673e-002_r_kind, -9.72827375e-002_r_kind, -9.76187348e-001_r_kind, -1.21772941e-002_r_kind, & -8.53854977e-003_r_kind, 6.21675106e-004_r_kind, -2.20879287e-004_r_kind, -1.55190664e-005_r_kind, & -3.20610404e-002_r_kind, 1.59049734e-001_r_kind, -8.63536596e-001_r_kind, -3.94894071e-002_r_kind, & -2.71239337e-002_r_kind, 6.07018499e-003_r_kind, 7.67658174e-004_r_kind, -1.19957927e-004_r_kind, & -6.88230246e-002_r_kind, 9.77551863e-002_r_kind, -9.07400787e-001_r_kind, 1.13464310e-003_r_kind, & -2.58067437e-002_r_kind, 4.45144856e-003_r_kind, 9.31801915e-004_r_kind, -1.49155152e-004_r_kind, & -7.73048252e-002_r_kind, -6.37678578e-002_r_kind, -9.94577944e-001_r_kind, 1.28959669e-002_r_kind, & -7.75872730e-003_r_kind, -1.52835681e-003_r_kind, -3.36256810e-004_r_kind, 4.88392907e-005_r_kind, & -5.15876301e-002_r_kind, -1.67232484e-001_r_kind, -1.02346241e+000_r_kind, 7.78744638e-004_r_kind, & -3.73176741e-003_r_kind, -4.07326850e-004_r_kind, -4.87616286e-004_r_kind, 1.14948148e-006_r_kind, & -1.01931587e-001_r_kind, -7.53218979e-002_r_kind, -9.00053680e-001_r_kind, -1.89454183e-002_r_kind, & -2.32969522e-002_r_kind, 4.95998468e-003_r_kind, 7.01531244e-004_r_kind, -1.08523644e-004_r_kind, & -9.89292935e-002_r_kind, -1.67843536e-001_r_kind, -1.01929891e+000_r_kind, 3.55926808e-004_r_kind, & -7.59821432e-003_r_kind, 4.09116549e-003_r_kind, 3.45942281e-005_r_kind, -5.99469204e-005_r_kind, & -9.26474035e-002_r_kind, 1.80861250e-001_r_kind, -1.02432120e+000_r_kind, -4.59155673e-003_r_kind, & -3.57208913e-003_r_kind, 8.09858844e-004_r_kind, -6.03338587e-004_r_kind, -1.04107357e-005_r_kind, & -1.41061982e-002_r_kind, 1.04682565e-001_r_kind, -1.02307773e+000_r_kind, 1.11318217e-003_r_kind, & -4.14737221e-003_r_kind, 2.20937582e-004_r_kind, -4.10551816e-004_r_kind, -2.26539196e-005_r_kind, & -1.22036435e-001_r_kind, -4.96090576e-002_r_kind, -9.67436492e-001_r_kind, -1.91587433e-002_r_kind, & -1.80738792e-002_r_kind, 4.41139704e-003_r_kind, 4.79028065e-004_r_kind, -6.78068609e-005_r_kind, & -1.14758760e-001_r_kind, -9.71220583e-002_r_kind, -9.63069022e-001_r_kind, -4.80926456e-003_r_kind, & -1.62106343e-002_r_kind, 5.26481727e-003_r_kind, 5.52151294e-004_r_kind, -1.55562826e-004_r_kind, & -2.97650043e-002_r_kind, -6.76572695e-002_r_kind, -9.97930944e-001_r_kind, 6.65108487e-003_r_kind, & -4.90541570e-003_r_kind, -1.19191140e-003_r_kind, -4.43050201e-004_r_kind, 4.20195611e-005_r_kind, & -1.01213023e-001_r_kind, 2.81079948e-001_r_kind, -9.86194432e-001_r_kind, 5.90509661e-002_r_kind, & -1.20200068e-002_r_kind, -1.81926333e-003_r_kind, -4.47495928e-004_r_kind, 7.30048778e-005_r_kind, & -2.97650043e-002_r_kind, -6.76572695e-002_r_kind, -9.97930944e-001_r_kind, 6.65108487e-003_r_kind, & -4.90541570e-003_r_kind, -1.19191140e-003_r_kind, -4.43050201e-004_r_kind, 4.20195611e-005_r_kind, & -1.01213023e-001_r_kind, 2.81079948e-001_r_kind, -9.86194432e-001_r_kind, 5.90509661e-002_r_kind, & -1.20200068e-002_r_kind, -1.81926333e-003_r_kind, -4.47495928e-004_r_kind, 7.30048778e-005_r_kind, & -2.97650043e-002_r_kind, -6.76572695e-002_r_kind, -9.97930944e-001_r_kind, 6.65108487e-003_r_kind, & -4.90541570e-003_r_kind, -1.19191140e-003_r_kind, -4.43050201e-004_r_kind, 4.20195611e-005_r_kind, & -1.01213023e-001_r_kind, 2.81079948e-001_r_kind, -9.86194432e-001_r_kind, 5.90509661e-002_r_kind, & -1.20200068e-002_r_kind, -1.81926333e-003_r_kind, -4.47495928e-004_r_kind, 7.30048778e-005_r_kind, & -2.97650043e-002_r_kind, -6.76572695e-002_r_kind, -9.97930944e-001_r_kind, 6.65108487e-003_r_kind, & -4.90541570e-003_r_kind, -1.19191140e-003_r_kind, -4.43050201e-004_r_kind, 4.20195611e-005_r_kind, & -1.01213023e-001_r_kind, 2.81079948e-001_r_kind, -9.86194432e-001_r_kind, 5.90509661e-002_r_kind, & -1.20200068e-002_r_kind, -1.81926333e-003_r_kind, -4.47495928e-004_r_kind, 7.30048778e-005_r_kind, & -2.97650043e-002_r_kind, -6.76572695e-002_r_kind, -9.97930944e-001_r_kind, 6.65108487e-003_r_kind, & -4.90541570e-003_r_kind, -1.19191140e-003_r_kind, -4.43050201e-004_r_kind, 4.20195611e-005_r_kind, & -1.01213023e-001_r_kind, 2.81079948e-001_r_kind, -9.86194432e-001_r_kind, 5.90509661e-002_r_kind, & -1.20200068e-002_r_kind, -1.81926333e-003_r_kind, -4.47495928e-004_r_kind, 7.30048778e-005_r_kind, & -2.97650043e-002_r_kind, -6.76572695e-002_r_kind, -9.97930944e-001_r_kind, 6.65108487e-003_r_kind, & -4.90541570e-003_r_kind, -1.19191140e-003_r_kind, -4.43050201e-004_r_kind, 4.20195611e-005_r_kind, & -1.01213023e-001_r_kind, 2.81079948e-001_r_kind, -9.86194432e-001_r_kind, 5.90509661e-002_r_kind, & -1.20200068e-002_r_kind, -1.81926333e-003_r_kind, -4.47495928e-004_r_kind, 7.30048778e-005_r_kind, & -1.05498888e-001_r_kind, 1.38700381e-001_r_kind, -9.18839693e-001_r_kind, 5.23423683e-003_r_kind, & -2.72325296e-002_r_kind, 4.52917721e-003_r_kind, 1.03306631e-003_r_kind, -1.61810633e-004_r_kind, & -1.89013511e-001_r_kind, -5.92692830e-002_r_kind, -9.70508218e-001_r_kind, 2.08470389e-001_r_kind, & -3.25677581e-002_r_kind, -7.15242233e-003_r_kind, 1.50414580e-003_r_kind, -4.58021022e-005_r_kind/),& (/ 8,2,15/) ) real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_19 = reshape( (/ & 2.88968161e-002_r_kind, 5.72211891e-002_r_kind, -8.97529483e-001_r_kind, 1.36518311e-002_r_kind, & -6.98799780e-003_r_kind, -1.08093256e-003_r_kind, -4.31537890e-004_r_kind, 2.89580021e-005_r_kind, & -8.22370946e-002_r_kind, 1.84141174e-001_r_kind, -9.10946906e-001_r_kind, -1.95769705e-002_r_kind, & -1.63720059e-003_r_kind, 2.61152891e-004_r_kind, -8.77573970e-004_r_kind, 6.08989940e-005_r_kind, & -1.22222744e-001_r_kind, 7.16604441e-002_r_kind, -7.71629572e-001_r_kind, -4.01977152e-002_r_kind, & -3.85947339e-002_r_kind, 6.56226650e-003_r_kind, 1.26806693e-003_r_kind, -1.73192995e-004_r_kind, & -1.74645901e-001_r_kind, -6.10361695e-002_r_kind, -6.59068048e-001_r_kind, 1.98529903e-002_r_kind, & -5.53864725e-002_r_kind, 6.79773558e-003_r_kind, 1.79760670e-003_r_kind, -3.09538184e-004_r_kind, & -1.04661090e-002_r_kind, 9.94574130e-002_r_kind, -9.82419491e-001_r_kind, -4.86512057e-004_r_kind, & -3.80066503e-003_r_kind, -1.88199090e-004_r_kind, -6.61909115e-004_r_kind, 1.81017949e-005_r_kind, & -8.66462104e-003_r_kind, 6.19573519e-002_r_kind, -1.00814879e+000_r_kind, 2.42177467e-003_r_kind, & -5.26584499e-003_r_kind, 2.81416695e-004_r_kind, -3.27982794e-004_r_kind, -1.35570344e-005_r_kind, & -1.29350603e-001_r_kind, -1.55049682e-001_r_kind, -8.12684000e-001_r_kind, -2.43408754e-002_r_kind, & -3.69030349e-002_r_kind, 3.85890855e-003_r_kind, 1.24171830e-003_r_kind, -9.62883350e-005_r_kind, & -1.13065429e-001_r_kind, 1.36988282e-001_r_kind, -8.88776898e-001_r_kind, 1.10446138e-003_r_kind, & -3.07487398e-002_r_kind, 5.69969090e-003_r_kind, 1.13500340e-003_r_kind, -2.08897967e-004_r_kind, & -8.47922172e-003_r_kind, 7.31471414e-003_r_kind, -1.00471723e+000_r_kind, 8.77908152e-003_r_kind, & -1.87136326e-003_r_kind, -1.34208635e-003_r_kind, -5.31051133e-004_r_kind, 6.27772606e-005_r_kind, & -3.90902273e-002_r_kind, 3.12932044e-001_r_kind, -1.03079987e+000_r_kind, -6.53349608e-003_r_kind, & -1.71896431e-003_r_kind, 1.65031443e-003_r_kind, -5.45169576e-004_r_kind, -2.35154166e-005_r_kind, & -1.06491432e-001_r_kind, 3.12845199e-003_r_kind, -8.48457575e-001_r_kind, -1.00595159e-002_r_kind, & -3.40486281e-002_r_kind, 4.14518919e-003_r_kind, 1.17543642e-003_r_kind, -1.27937790e-004_r_kind, & -6.93388954e-002_r_kind, 2.57634282e-001_r_kind, -9.94598567e-001_r_kind, 6.70428481e-003_r_kind, & -1.74006820e-002_r_kind, 4.02611215e-003_r_kind, 4.36087867e-004_r_kind, -8.74491816e-005_r_kind, & -2.71350276e-002_r_kind, -2.07411796e-001_r_kind, -1.00102842e+000_r_kind, -3.59892892e-003_r_kind, & -3.32985423e-003_r_kind, -1.80928022e-004_r_kind, -5.37585292e-004_r_kind, 1.40997736e-005_r_kind, & -6.48028627e-002_r_kind, -9.58651211e-003_r_kind, -1.05065870e+000_r_kind, -6.63892739e-003_r_kind, & -4.75495541e-003_r_kind, 3.86347063e-004_r_kind, -3.41132225e-004_r_kind, -1.49872949e-005_r_kind, & -5.01554534e-002_r_kind, -1.95738703e-001_r_kind, -8.77366960e-001_r_kind, -4.23899442e-002_r_kind, & -2.77511068e-002_r_kind, 6.38218131e-003_r_kind, 9.16048361e-004_r_kind, -1.49990592e-004_r_kind, & -1.54721767e-001_r_kind, 1.00917861e-001_r_kind, -1.01476681e+000_r_kind, -1.35878744e-002_r_kind, & -1.70912687e-002_r_kind, 6.88493345e-003_r_kind, 6.32103591e-004_r_kind, -2.11467093e-004_r_kind, & -1.13729835e-001_r_kind, -6.83060661e-002_r_kind, -9.83069777e-001_r_kind, 6.39191689e-003_r_kind, & -5.39532816e-003_r_kind, -1.23508472e-003_r_kind, -3.51379276e-004_r_kind, 4.07678272e-005_r_kind, & -3.75476144e-002_r_kind, 1.21206477e-001_r_kind, -1.05813158e+000_r_kind, 2.07087342e-002_r_kind, & -4.67697438e-003_r_kind, -1.32091378e-003_r_kind, -5.81131026e-004_r_kind, 3.38749996e-005_r_kind, & -1.13729835e-001_r_kind, -6.83060661e-002_r_kind, -9.83069777e-001_r_kind, 6.39191689e-003_r_kind, & -5.39532816e-003_r_kind, -1.23508472e-003_r_kind, -3.51379276e-004_r_kind, 4.07678272e-005_r_kind, & -3.75476144e-002_r_kind, 1.21206477e-001_r_kind, -1.05813158e+000_r_kind, 2.07087342e-002_r_kind, & -4.67697438e-003_r_kind, -1.32091378e-003_r_kind, -5.81131026e-004_r_kind, 3.38749996e-005_r_kind, & -1.13729835e-001_r_kind, -6.83060661e-002_r_kind, -9.83069777e-001_r_kind, 6.39191689e-003_r_kind, & -5.39532816e-003_r_kind, -1.23508472e-003_r_kind, -3.51379276e-004_r_kind, 4.07678272e-005_r_kind, & -3.75476144e-002_r_kind, 1.21206477e-001_r_kind, -1.05813158e+000_r_kind, 2.07087342e-002_r_kind, & -4.67697438e-003_r_kind, -1.32091378e-003_r_kind, -5.81131026e-004_r_kind, 3.38749996e-005_r_kind, & -1.13729835e-001_r_kind, -6.83060661e-002_r_kind, -9.83069777e-001_r_kind, 6.39191689e-003_r_kind, & -5.39532816e-003_r_kind, -1.23508472e-003_r_kind, -3.51379276e-004_r_kind, 4.07678272e-005_r_kind, & -3.75476144e-002_r_kind, 1.21206477e-001_r_kind, -1.05813158e+000_r_kind, 2.07087342e-002_r_kind, & -4.67697438e-003_r_kind, -1.32091378e-003_r_kind, -5.81131026e-004_r_kind, 3.38749996e-005_r_kind, & -1.13729835e-001_r_kind, -6.83060661e-002_r_kind, -9.83069777e-001_r_kind, 6.39191689e-003_r_kind, & -5.39532816e-003_r_kind, -1.23508472e-003_r_kind, -3.51379276e-004_r_kind, 4.07678272e-005_r_kind, & -3.75476144e-002_r_kind, 1.21206477e-001_r_kind, -1.05813158e+000_r_kind, 2.07087342e-002_r_kind, & -4.67697438e-003_r_kind, -1.32091378e-003_r_kind, -5.81131026e-004_r_kind, 3.38749996e-005_r_kind, & -1.13729835e-001_r_kind, -6.83060661e-002_r_kind, -9.83069777e-001_r_kind, 6.39191689e-003_r_kind, & -5.39532816e-003_r_kind, -1.23508472e-003_r_kind, -3.51379276e-004_r_kind, 4.07678272e-005_r_kind, & -3.75476144e-002_r_kind, 1.21206477e-001_r_kind, -1.05813158e+000_r_kind, 2.07087342e-002_r_kind, & -4.67697438e-003_r_kind, -1.32091378e-003_r_kind, -5.81131026e-004_r_kind, 3.38749996e-005_r_kind, & -1.09833397e-001_r_kind, 8.32149014e-002_r_kind, -9.15308714e-001_r_kind, -9.75017436e-003_r_kind, & -2.72793751e-002_r_kind, 5.97133115e-003_r_kind, 9.92798945e-004_r_kind, -1.96212786e-004_r_kind, & -1.12866320e-001_r_kind, 7.99755678e-002_r_kind, -1.04918265e+000_r_kind, 2.13572606e-001_r_kind, & -2.87145432e-002_r_kind, -6.76067220e-003_r_kind, 1.32346945e-003_r_kind, -6.15322861e-005_r_kind/),& (/ 8,2,15/) ) real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_20 = reshape( (/ & -5.53324930e-002_r_kind, 8.96679163e-002_r_kind, -9.26657975e-001_r_kind, 9.66761820e-003_r_kind, & 2.83907837e-004_r_kind, -6.78065931e-004_r_kind, -6.63810060e-004_r_kind, 2.54925835e-005_r_kind, & -5.02091786e-003_r_kind, -1.77881233e-002_r_kind, -1.09147155e+000_r_kind, -4.96294815e-004_r_kind, & 1.18014142e-002_r_kind, -1.61008938e-005_r_kind, -1.16086623e-003_r_kind, -4.97733663e-006_r_kind, & -1.16351970e-001_r_kind, 1.33443147e-001_r_kind, -7.46442199e-001_r_kind, 1.89830922e-002_r_kind, & -3.65450568e-002_r_kind, 5.66841278e-004_r_kind, 1.10465661e-003_r_kind, -3.34551769e-005_r_kind, & -4.48111631e-002_r_kind, -3.16455811e-002_r_kind, -9.74579811e-001_r_kind, 9.11674425e-002_r_kind, & -1.72947571e-002_r_kind, -1.92617078e-003_r_kind, 3.38681246e-004_r_kind, 1.86966990e-005_r_kind, & -2.66249888e-002_r_kind, 4.51199748e-002_r_kind, -1.00546598e+000_r_kind, 9.39478632e-004_r_kind, & -3.53710214e-003_r_kind, -1.79096489e-004_r_kind, -6.82093261e-004_r_kind, -4.31082117e-006_r_kind, & -5.99183887e-002_r_kind, 1.30984396e-001_r_kind, -9.32524204e-001_r_kind, -1.27927994e-003_r_kind, & -1.63929146e-002_r_kind, 4.12638718e-003_r_kind, 4.37202369e-004_r_kind, -9.56427102e-005_r_kind, & -1.55220898e-002_r_kind, -6.09744387e-003_r_kind, -9.12356019e-001_r_kind, -1.93568710e-002_r_kind, & -1.92841608e-002_r_kind, 3.54133872e-003_r_kind, 3.17921542e-004_r_kind, -9.52724877e-005_r_kind, & -6.82753250e-002_r_kind, 1.66306257e-001_r_kind, -9.80789363e-001_r_kind, 8.84867087e-003_r_kind, & -5.49255032e-003_r_kind, 1.29161309e-003_r_kind, -1.96004970e-004_r_kind, -2.83482495e-005_r_kind, & -1.36695012e-001_r_kind, -1.94087148e-001_r_kind, -9.40238893e-001_r_kind, -1.52252773e-002_r_kind, & -2.30028909e-002_r_kind, 2.24579475e-003_r_kind, 5.63844922e-004_r_kind, -2.53392736e-005_r_kind, & -7.49450773e-002_r_kind, 1.62285313e-001_r_kind, -9.37521935e-001_r_kind, -6.94234809e-003_r_kind, & -1.34419724e-002_r_kind, 3.78534500e-003_r_kind, 1.85084413e-004_r_kind, -1.01173340e-004_r_kind, & -1.64406467e-002_r_kind, 1.01123825e-001_r_kind, -9.76054192e-001_r_kind, 2.07524165e-003_r_kind, & -8.59644078e-003_r_kind, 9.05853347e-004_r_kind, -2.01289775e-004_r_kind, -1.85748886e-005_r_kind, & -1.02546252e-002_r_kind, -1.13465235e-001_r_kind, -9.96470332e-001_r_kind, 1.32228865e-003_r_kind, & -6.50392054e-003_r_kind, 5.66910137e-004_r_kind, -1.64041601e-004_r_kind, -1.60055206e-005_r_kind, & -1.17206797e-001_r_kind, -1.87622290e-002_r_kind, -9.34394658e-001_r_kind, -2.15133242e-002_r_kind, & -2.12136190e-002_r_kind, 4.24697716e-003_r_kind, 5.30900143e-004_r_kind, -1.18717391e-004_r_kind, & -9.86387860e-003_r_kind, 2.19956830e-001_r_kind, -1.00827980e+000_r_kind, 2.09442037e-003_r_kind, & -2.11799294e-002_r_kind, 4.06272896e-003_r_kind, 5.31111960e-004_r_kind, -1.27632244e-004_r_kind, & -8.87944084e-003_r_kind, 2.99282800e-002_r_kind, -1.01551545e+000_r_kind, 4.39770659e-003_r_kind, & -9.93058924e-003_r_kind, -7.66939920e-005_r_kind, -1.68053346e-004_r_kind, 4.35948896e-005_r_kind, & -2.19885521e-002_r_kind, 1.00838229e-001_r_kind, -1.00064743e+000_r_kind, -1.63910221e-002_r_kind, & -1.11649465e-002_r_kind, 4.00572736e-003_r_kind, 1.80401798e-006_r_kind, -1.32877278e-004_r_kind, & -6.20569475e-002_r_kind, 1.00097619e-002_r_kind, -9.41858709e-001_r_kind, 1.25842076e-003_r_kind, & -1.62702035e-002_r_kind, 1.49221520e-003_r_kind, 2.20242451e-004_r_kind, -1.20322366e-005_r_kind, & -5.10414988e-002_r_kind, -1.43528238e-001_r_kind, -1.00225520e+000_r_kind, -6.67478470e-003_r_kind, & -6.59647910e-003_r_kind, 3.82588268e-003_r_kind, 1.02092381e-005_r_kind, -1.06037653e-004_r_kind, & -6.20569475e-002_r_kind, 1.00097619e-002_r_kind, -9.41858709e-001_r_kind, 1.25842076e-003_r_kind, & -1.62702035e-002_r_kind, 1.49221520e-003_r_kind, 2.20242451e-004_r_kind, -1.20322366e-005_r_kind, & -5.10414988e-002_r_kind, -1.43528238e-001_r_kind, -1.00225520e+000_r_kind, -6.67478470e-003_r_kind, & -6.59647910e-003_r_kind, 3.82588268e-003_r_kind, 1.02092381e-005_r_kind, -1.06037653e-004_r_kind, & -6.20569475e-002_r_kind, 1.00097619e-002_r_kind, -9.41858709e-001_r_kind, 1.25842076e-003_r_kind, & -1.62702035e-002_r_kind, 1.49221520e-003_r_kind, 2.20242451e-004_r_kind, -1.20322366e-005_r_kind, & -5.10414988e-002_r_kind, -1.43528238e-001_r_kind, -1.00225520e+000_r_kind, -6.67478470e-003_r_kind, & -6.59647910e-003_r_kind, 3.82588268e-003_r_kind, 1.02092381e-005_r_kind, -1.06037653e-004_r_kind, & -6.20569475e-002_r_kind, 1.00097619e-002_r_kind, -9.41858709e-001_r_kind, 1.25842076e-003_r_kind, & -1.62702035e-002_r_kind, 1.49221520e-003_r_kind, 2.20242451e-004_r_kind, -1.20322366e-005_r_kind, & -5.10414988e-002_r_kind, -1.43528238e-001_r_kind, -1.00225520e+000_r_kind, -6.67478470e-003_r_kind, & -6.59647910e-003_r_kind, 3.82588268e-003_r_kind, 1.02092381e-005_r_kind, -1.06037653e-004_r_kind, & -6.20569475e-002_r_kind, 1.00097619e-002_r_kind, -9.41858709e-001_r_kind, 1.25842076e-003_r_kind, & -1.62702035e-002_r_kind, 1.49221520e-003_r_kind, 2.20242451e-004_r_kind, -1.20322366e-005_r_kind, & -5.10414988e-002_r_kind, -1.43528238e-001_r_kind, -1.00225520e+000_r_kind, -6.67478470e-003_r_kind, & -6.59647910e-003_r_kind, 3.82588268e-003_r_kind, 1.02092381e-005_r_kind, -1.06037653e-004_r_kind, & -6.20569475e-002_r_kind, 1.00097619e-002_r_kind, -9.41858709e-001_r_kind, 1.25842076e-003_r_kind, & -1.62702035e-002_r_kind, 1.49221520e-003_r_kind, 2.20242451e-004_r_kind, -1.20322366e-005_r_kind, & -5.10414988e-002_r_kind, -1.43528238e-001_r_kind, -1.00225520e+000_r_kind, -6.67478470e-003_r_kind, & -6.59647910e-003_r_kind, 3.82588268e-003_r_kind, 1.02092381e-005_r_kind, -1.06037653e-004_r_kind, & -9.63953212e-002_r_kind, 5.01724146e-002_r_kind, -9.92312431e-001_r_kind, -9.74688586e-003_r_kind, & -7.74239656e-003_r_kind, 1.17365213e-003_r_kind, -2.19702270e-004_r_kind, -4.11848814e-005_r_kind, & -5.69399334e-002_r_kind, -2.14651451e-001_r_kind, -8.92655075e-001_r_kind, 2.09093615e-002_r_kind, & -1.54947042e-002_r_kind, -1.60737406e-003_r_kind, -9.41211692e-005_r_kind, 8.22987568e-005_r_kind/),& (/ 8,2,15/) ) real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_21 = reshape( (/ & -5.93122318e-002_r_kind, 1.51160493e-001_r_kind, -9.29837704e-001_r_kind, 1.11907059e-002_r_kind, & -2.92930426e-003_r_kind, -5.37499378e-004_r_kind, -4.92306892e-004_r_kind, 2.17214656e-005_r_kind, & 9.18469334e-004_r_kind, -9.17405635e-002_r_kind, -1.02836871e+000_r_kind, -1.24654919e-002_r_kind, & 1.02571091e-002_r_kind, 1.73965876e-003_r_kind, -1.22978329e-003_r_kind, -9.14341945e-005_r_kind, & -2.52755154e-002_r_kind, 1.59815118e-001_r_kind, -7.83926487e-001_r_kind, -3.39760706e-002_r_kind, & -2.61132009e-002_r_kind, 5.50437253e-003_r_kind, 6.93825656e-004_r_kind, -1.08948530e-004_r_kind, & -1.53507411e-001_r_kind, -8.26360360e-002_r_kind, -7.57060170e-001_r_kind, 8.65705162e-002_r_kind, & -4.13534306e-002_r_kind, 3.42966465e-004_r_kind, 1.35908509e-003_r_kind, -1.21754485e-004_r_kind, & -3.82020660e-002_r_kind, -1.15718156e-001_r_kind, -1.00587749e+000_r_kind, -5.29566593e-003_r_kind, & -4.67831036e-003_r_kind, 1.97040470e-004_r_kind, -6.34025782e-004_r_kind, -2.47025910e-005_r_kind, & 2.69212294e-003_r_kind, 5.70879392e-002_r_kind, -1.01243186e+000_r_kind, -3.38277477e-003_r_kind, & -4.87565668e-003_r_kind, 4.40029398e-004_r_kind, -3.20241059e-004_r_kind, -1.64663506e-005_r_kind, & -7.73688331e-002_r_kind, 1.66169301e-001_r_kind, -8.98579121e-001_r_kind, -2.15497073e-002_r_kind, & -2.49314401e-002_r_kind, 4.49157134e-003_r_kind, 7.54829322e-004_r_kind, -9.55032301e-005_r_kind, & -1.75007015e-001_r_kind, 3.45267087e-001_r_kind, -8.83400679e-001_r_kind, 8.71584751e-003_r_kind, & -3.31111401e-002_r_kind, 5.04212501e-003_r_kind, 1.33993302e-003_r_kind, -2.30460922e-004_r_kind, & -9.48398039e-002_r_kind, 1.06142253e-001_r_kind, -1.01354051e+000_r_kind, 1.22573189e-002_r_kind, & -6.30498072e-003_r_kind, -1.53515034e-003_r_kind, -4.57574526e-004_r_kind, 5.10238351e-005_r_kind, & -6.70555159e-002_r_kind, 2.83843815e-001_r_kind, -1.04732060e+000_r_kind, 1.32900607e-002_r_kind, & -4.08430025e-003_r_kind, -1.00952887e-003_r_kind, -5.14037441e-004_r_kind, 5.12171864e-005_r_kind, & -1.05284780e-001_r_kind, 2.67886724e-002_r_kind, -8.93791139e-001_r_kind, -1.46267274e-002_r_kind, & -2.56743617e-002_r_kind, 4.89959354e-003_r_kind, 8.27621494e-004_r_kind, -1.19767203e-004_r_kind, & -1.38652995e-001_r_kind, -5.68602562e-001_r_kind, -9.37645018e-001_r_kind, -1.81222521e-002_r_kind, & -1.73796881e-002_r_kind, 6.10630121e-003_r_kind, 5.63356152e-004_r_kind, -1.54547743e-004_r_kind, & -2.72973537e-001_r_kind, 1.08176279e+000_r_kind, -1.04552341e+000_r_kind, -1.02911890e-002_r_kind, & -7.50666298e-003_r_kind, 4.90204804e-003_r_kind, -1.11065072e-003_r_kind, -1.94243585e-005_r_kind, & -3.24335620e-002_r_kind, -1.26401782e-001_r_kind, -1.02597928e+000_r_kind, -4.94325487e-003_r_kind, & -7.71400006e-003_r_kind, 3.96533520e-004_r_kind, -1.37512645e-004_r_kind, -1.33189415e-005_r_kind, & -7.26768821e-002_r_kind, 2.59305775e-001_r_kind, -9.59965587e-001_r_kind, -1.57393645e-002_r_kind, & -1.87098477e-002_r_kind, 4.90527926e-003_r_kind, 4.25422186e-004_r_kind, -7.98054389e-005_r_kind, & -8.71768296e-002_r_kind, -8.85518044e-002_r_kind, -9.95310962e-001_r_kind, -1.64069664e-002_r_kind, & -1.84447765e-002_r_kind, 6.60658721e-003_r_kind, 7.05660321e-004_r_kind, -2.16119399e-004_r_kind, & -3.70209366e-002_r_kind, 8.34630653e-002_r_kind, -9.68965590e-001_r_kind, -6.49541570e-003_r_kind, & -4.04640893e-003_r_kind, 1.03803957e-003_r_kind, -4.51976492e-004_r_kind, -3.95033567e-005_r_kind, & -1.12812385e-001_r_kind, -1.95549324e-001_r_kind, -9.88134205e-001_r_kind, -4.85489611e-003_r_kind, & -2.32097469e-002_r_kind, 2.60839798e-003_r_kind, 4.25688369e-004_r_kind, -1.21078338e-004_r_kind, & -3.70209366e-002_r_kind, 8.34630653e-002_r_kind, -9.68965590e-001_r_kind, -6.49541570e-003_r_kind, & -4.04640893e-003_r_kind, 1.03803957e-003_r_kind, -4.51976492e-004_r_kind, -3.95033567e-005_r_kind, & -1.12812385e-001_r_kind, -1.95549324e-001_r_kind, -9.88134205e-001_r_kind, -4.85489611e-003_r_kind, & -2.32097469e-002_r_kind, 2.60839798e-003_r_kind, 4.25688369e-004_r_kind, -1.21078338e-004_r_kind, & -3.70209366e-002_r_kind, 8.34630653e-002_r_kind, -9.68965590e-001_r_kind, -6.49541570e-003_r_kind, & -4.04640893e-003_r_kind, 1.03803957e-003_r_kind, -4.51976492e-004_r_kind, -3.95033567e-005_r_kind, & -1.12812385e-001_r_kind, -1.95549324e-001_r_kind, -9.88134205e-001_r_kind, -4.85489611e-003_r_kind, & -2.32097469e-002_r_kind, 2.60839798e-003_r_kind, 4.25688369e-004_r_kind, -1.21078338e-004_r_kind, & -3.70209366e-002_r_kind, 8.34630653e-002_r_kind, -9.68965590e-001_r_kind, -6.49541570e-003_r_kind, & -4.04640893e-003_r_kind, 1.03803957e-003_r_kind, -4.51976492e-004_r_kind, -3.95033567e-005_r_kind, & -1.12812385e-001_r_kind, -1.95549324e-001_r_kind, -9.88134205e-001_r_kind, -4.85489611e-003_r_kind, & -2.32097469e-002_r_kind, 2.60839798e-003_r_kind, 4.25688369e-004_r_kind, -1.21078338e-004_r_kind, & -3.70209366e-002_r_kind, 8.34630653e-002_r_kind, -9.68965590e-001_r_kind, -6.49541570e-003_r_kind, & -4.04640893e-003_r_kind, 1.03803957e-003_r_kind, -4.51976492e-004_r_kind, -3.95033567e-005_r_kind, & -1.12812385e-001_r_kind, -1.95549324e-001_r_kind, -9.88134205e-001_r_kind, -4.85489611e-003_r_kind, & -2.32097469e-002_r_kind, 2.60839798e-003_r_kind, 4.25688369e-004_r_kind, -1.21078338e-004_r_kind, & -3.70209366e-002_r_kind, 8.34630653e-002_r_kind, -9.68965590e-001_r_kind, -6.49541570e-003_r_kind, & -4.04640893e-003_r_kind, 1.03803957e-003_r_kind, -4.51976492e-004_r_kind, -3.95033567e-005_r_kind, & -1.12812385e-001_r_kind, -1.95549324e-001_r_kind, -9.88134205e-001_r_kind, -4.85489611e-003_r_kind, & -2.32097469e-002_r_kind, 2.60839798e-003_r_kind, 4.25688369e-004_r_kind, -1.21078338e-004_r_kind, & -1.46404251e-001_r_kind, 1.02264553e-001_r_kind, -8.96511912e-001_r_kind, -8.85267928e-003_r_kind, & -2.83047054e-002_r_kind, 5.31384861e-003_r_kind, 1.06912781e-003_r_kind, -1.96049135e-004_r_kind, & -1.04192808e-001_r_kind, -3.60217482e-001_r_kind, -9.75351274e-001_r_kind, 1.85505807e-001_r_kind, & -2.61770226e-002_r_kind, -6.75599044e-003_r_kind, 1.03888079e-003_r_kind, -2.92561799e-005_r_kind/),& (/ 8,2,15/) ) real(r_kind), target, private, dimension(0:7,2,15) :: amsucoeff_22 = reshape( (/ & 9.75104049e-003_r_kind, -2.67954972e-002_r_kind, -1.00441885e+000_r_kind, -4.13552392e-003_r_kind, & 4.83619282e-003_r_kind, 2.23311290e-004_r_kind, -9.30484326e-004_r_kind, -1.82476233e-005_r_kind, & 8.98046792e-003_r_kind, -1.36037081e-001_r_kind, -1.01462388e+000_r_kind, 2.42774701e-003_r_kind, & 1.15878545e-002_r_kind, -9.88824759e-004_r_kind, -1.09468901e-003_r_kind, 3.26434965e-005_r_kind, & -1.26818016e-001_r_kind, 5.55154160e-002_r_kind, -8.69268537e-001_r_kind, -1.17066959e-002_r_kind, & -2.80337539e-002_r_kind, 3.85268894e-003_r_kind, 9.08550981e-004_r_kind, -9.60082616e-005_r_kind, & -5.62816784e-002_r_kind, -2.38852397e-001_r_kind, -8.00485492e-001_r_kind, 8.28931779e-002_r_kind, & -2.70346627e-002_r_kind, -1.18868623e-003_r_kind, 5.74906473e-004_r_kind, 7.22288223e-006_r_kind, & -8.62975866e-002_r_kind, -6.18472397e-001_r_kind, -9.35284019e-001_r_kind, -1.79203991e-002_r_kind, & -6.34222571e-003_r_kind, -5.96674741e-004_r_kind, -6.55032927e-004_r_kind, -3.02471090e-005_r_kind, & -7.36862645e-002_r_kind, -5.46863616e-001_r_kind, -9.86678123e-001_r_kind, -5.25613595e-003_r_kind, & -3.89240379e-003_r_kind, 1.29066524e-003_r_kind, -3.13781522e-004_r_kind, -3.81810169e-005_r_kind, & -1.33809060e-001_r_kind, 3.30931365e-001_r_kind, -7.22728193e-001_r_kind, 1.24727814e-002_r_kind, & -4.63315435e-002_r_kind, 3.82323237e-003_r_kind, 1.52931688e-003_r_kind, -1.58564944e-004_r_kind, & -5.35287187e-002_r_kind, -1.78305641e-001_r_kind, -9.91991997e-001_r_kind, -2.25094184e-002_r_kind, & -1.11914258e-002_r_kind, 3.07255681e-003_r_kind, -2.32107377e-005_r_kind, -1.36529721e-004_r_kind, & -1.16859958e-001_r_kind, -6.38576865e-001_r_kind, -9.91750658e-001_r_kind, -1.42509583e-002_r_kind, & -5.06967632e-003_r_kind, -1.73841341e-004_r_kind, -4.88598947e-004_r_kind, -3.51224116e-005_r_kind, & -6.07711682e-003_r_kind, 1.11489547e-002_r_kind, -1.00144279e+000_r_kind, 1.48511874e-002_r_kind, & -1.19725578e-002_r_kind, 1.70405302e-003_r_kind, -3.94240451e-005_r_kind, -3.69073427e-006_r_kind, & -1.43366650e-001_r_kind, 5.62716901e-001_r_kind, -8.58107805e-001_r_kind, 3.04583721e-002_r_kind, & -3.56974378e-002_r_kind, 3.12131504e-003_r_kind, 1.21328817e-003_r_kind, -1.34941030e-004_r_kind, & -3.02330963e-003_r_kind, 9.48470738e-003_r_kind, -1.04689682e+000_r_kind, 7.54171535e-002_r_kind, & -1.16670355e-002_r_kind, -3.58890556e-003_r_kind, 3.09759111e-004_r_kind, -3.82375001e-005_r_kind, & -1.23668313e-002_r_kind, -1.36103466e-001_r_kind, -9.72598374e-001_r_kind, -1.09806759e-002_r_kind, & -3.20895668e-003_r_kind, 9.01909953e-004_r_kind, -6.85705047e-004_r_kind, -6.03781882e-005_r_kind, & 7.78510468e-003_r_kind, -1.76615104e-001_r_kind, -1.00559688e+000_r_kind, -7.88976159e-003_r_kind, & -8.16779677e-003_r_kind, 4.48832056e-003_r_kind, 3.61518105e-005_r_kind, -1.19117700e-004_r_kind, & -1.47040308e-001_r_kind, 5.81465840e-001_r_kind, -8.27971518e-001_r_kind, 4.39170562e-002_r_kind, & -4.99625541e-002_r_kind, 3.92394979e-003_r_kind, 1.96118420e-003_r_kind, -2.30956604e-004_r_kind, & -8.52618814e-002_r_kind, 2.70994157e-001_r_kind, -1.01976037e+000_r_kind, 5.55564500e-002_r_kind, & -1.30984532e-002_r_kind, -7.25675840e-004_r_kind, 1.76902213e-005_r_kind, 2.08494976e-006_r_kind, & -2.37837387e-003_r_kind, -2.99516082e-001_r_kind, -9.47261035e-001_r_kind, 1.36433728e-002_r_kind, & -7.83403311e-003_r_kind, -1.57008832e-003_r_kind, -2.59426859e-004_r_kind, 4.62863754e-005_r_kind, & -1.11635305e-001_r_kind, -4.89074111e-001_r_kind, -9.05803621e-001_r_kind, 4.85747866e-002_r_kind, & -2.35798005e-002_r_kind, -2.64092139e-003_r_kind, 4.40482574e-004_r_kind, 2.19099729e-005_r_kind, & -2.37837387e-003_r_kind, -2.99516082e-001_r_kind, -9.47261035e-001_r_kind, 1.36433728e-002_r_kind, & -7.83403311e-003_r_kind, -1.57008832e-003_r_kind, -2.59426859e-004_r_kind, 4.62863754e-005_r_kind, & -1.11635305e-001_r_kind, -4.89074111e-001_r_kind, -9.05803621e-001_r_kind, 4.85747866e-002_r_kind, & -2.35798005e-002_r_kind, -2.64092139e-003_r_kind, 4.40482574e-004_r_kind, 2.19099729e-005_r_kind, & -2.37837387e-003_r_kind, -2.99516082e-001_r_kind, -9.47261035e-001_r_kind, 1.36433728e-002_r_kind, & -7.83403311e-003_r_kind, -1.57008832e-003_r_kind, -2.59426859e-004_r_kind, 4.62863754e-005_r_kind, & -1.11635305e-001_r_kind, -4.89074111e-001_r_kind, -9.05803621e-001_r_kind, 4.85747866e-002_r_kind, & -2.35798005e-002_r_kind, -2.64092139e-003_r_kind, 4.40482574e-004_r_kind, 2.19099729e-005_r_kind, & -2.37837387e-003_r_kind, -2.99516082e-001_r_kind, -9.47261035e-001_r_kind, 1.36433728e-002_r_kind, & -7.83403311e-003_r_kind, -1.57008832e-003_r_kind, -2.59426859e-004_r_kind, 4.62863754e-005_r_kind, & -1.11635305e-001_r_kind, -4.89074111e-001_r_kind, -9.05803621e-001_r_kind, 4.85747866e-002_r_kind, & -2.35798005e-002_r_kind, -2.64092139e-003_r_kind, 4.40482574e-004_r_kind, 2.19099729e-005_r_kind, & -2.37837387e-003_r_kind, -2.99516082e-001_r_kind, -9.47261035e-001_r_kind, 1.36433728e-002_r_kind, & -7.83403311e-003_r_kind, -1.57008832e-003_r_kind, -2.59426859e-004_r_kind, 4.62863754e-005_r_kind, & -1.11635305e-001_r_kind, -4.89074111e-001_r_kind, -9.05803621e-001_r_kind, 4.85747866e-002_r_kind, & -2.35798005e-002_r_kind, -2.64092139e-003_r_kind, 4.40482574e-004_r_kind, 2.19099729e-005_r_kind, & -2.37837387e-003_r_kind, -2.99516082e-001_r_kind, -9.47261035e-001_r_kind, 1.36433728e-002_r_kind, & -7.83403311e-003_r_kind, -1.57008832e-003_r_kind, -2.59426859e-004_r_kind, 4.62863754e-005_r_kind, & -1.11635305e-001_r_kind, -4.89074111e-001_r_kind, -9.05803621e-001_r_kind, 4.85747866e-002_r_kind, & -2.35798005e-002_r_kind, -2.64092139e-003_r_kind, 4.40482574e-004_r_kind, 2.19099729e-005_r_kind, & -8.75972584e-002_r_kind, -1.67156681e-001_r_kind, -9.24742103e-001_r_kind, -1.26165589e-002_r_kind, & -1.95397697e-002_r_kind, 5.43559855e-003_r_kind, 6.36994955e-004_r_kind, -1.51825137e-004_r_kind, & -9.06753615e-002_r_kind, 6.05665110e-002_r_kind, -1.01096153e+000_r_kind, 1.62782580e-001_r_kind, & -2.76400503e-002_r_kind, -4.19872068e-003_r_kind, 1.38002378e-003_r_kind, -9.29419766e-005_r_kind/) , & (/ 8,2,15/) ) real(r_kind), pointer, private, dimension(:,:,:) :: mhscoeff ! coefficients for mhs. each satellite in a separate data statement. real(r_kind), target, private, dimension(0:7,2,5) :: mhscoeff_18 = Reshape( (/ & 3.4965403e-002_r_kind, 2.1628516e+000_r_kind, -9.5464048e+000_r_kind, -1.6210480e+000_r_kind, & -1.4857870e+000_r_kind, 1.3150197e+000_r_kind, 5.7198799e-001_r_kind, -2.8890452e-001_r_kind, & -5.9196603e-002_r_kind, -6.7636013e-001_r_kind, -8.6714487e+000_r_kind, -1.6742078e+000_r_kind, & -2.1325920e+000_r_kind, 1.3199428e+000_r_kind, 4.0995410e-001_r_kind, -3.3879450e-001_r_kind, & -6.1441816e-002_r_kind, 1.4736807e+000_r_kind, -9.5058241e+000_r_kind, 1.3065364e+000_r_kind, & -2.4255955e+000_r_kind, -1.2079014e+000_r_kind, 1.7266749e-001_r_kind, 5.6346798e-001_r_kind, & -2.3348415e-002_r_kind, -1.1687171e+000_r_kind, -1.0514970e+001_r_kind, -7.4772894e-001_r_kind, & -9.8419696e-002_r_kind, 8.1826812e-001_r_kind, -1.2856618e+000_r_kind, -4.6397689e-001_r_kind, & -1.5548328e-001_r_kind, 3.4516063e+000_r_kind, -9.3134050e+000_r_kind, -3.9432240e+000_r_kind, & -1.0626595e+000_r_kind, 2.2548215e+000_r_kind, 4.0096134e-001_r_kind, -3.6071670e-001_r_kind, & 2.1672262e-002_r_kind, -1.2760131e+000_r_kind, -1.0493879e+001_r_kind, 5.9861401e-003_r_kind, & 2.6603895e-001_r_kind, -2.6706794e-001_r_kind, 1.1921029e-001_r_kind, 1.0875686e-001_r_kind, & -1.5548328e-001_r_kind, 3.4516063e+000_r_kind, -9.3134050e+000_r_kind, -3.9432240e+000_r_kind, & -1.0626595e+000_r_kind, 2.2548215e+000_r_kind, 4.0096134e-001_r_kind, -3.6071670e-001_r_kind, & 2.1672262e-002_r_kind, -1.2760131e+000_r_kind, -1.0493879e+001_r_kind, 5.9861401e-003_r_kind, & 2.6603895e-001_r_kind, -2.6706794e-001_r_kind, 1.1921029e-001_r_kind, 1.0875686e-001_r_kind, & -8.0707878e-002_r_kind, 2.7671516e+000_r_kind, -1.0256569e+001_r_kind, -1.5148641e+000_r_kind, & 2.3168449e-001_r_kind, 6.6968900e-001_r_kind, 1.5130611e-001_r_kind, -1.0532290e-001_r_kind, & -1.4505735e-001_r_kind, -1.4259889e+000_r_kind, -9.2915115e+000_r_kind, 9.5977956e-001_r_kind, & 2.6349644e-003_r_kind, -6.5305108e-001_r_kind, 1.1358295e-001_r_kind, 1.2522188e-001_r_kind/), & (/ 8,2,5 /) ) real(r_kind), target, private, dimension(0:7,2,5) :: mhscoeff_19 = Reshape( (/ & 1.6703354e-001_r_kind, -8.1146163e-001_r_kind, -1.0738719e+001_r_kind, 2.2400708e+000_r_kind, & -5.1594537e-001_r_kind, -1.6325672e+000_r_kind, 4.3368569e-001_r_kind, 3.2888728e-001_r_kind, & -6.6353470e-002_r_kind, -1.1122364e+000_r_kind, -9.9416256e+000_r_kind, 1.1250391e+000_r_kind, & -1.1115935e-001_r_kind, -1.6518873e+000_r_kind, -6.3417786e-001_r_kind, 4.7320408e-001_r_kind, & -7.8627609e-002_r_kind, 6.9862741e-001_r_kind, -6.9884124e+000_r_kind, -1.5533140e+000_r_kind, & -4.4324641e+000_r_kind, 1.1733495e+000_r_kind, 1.0958253e+000_r_kind, -2.0334625e-001_r_kind, & -5.4213259e-002_r_kind, 3.9298350e-001_r_kind, -7.3326321e+000_r_kind, 6.7302763e-001_r_kind, & -2.8431201e+000_r_kind, -6.1369413e-001_r_kind, -3.0428585e-001_r_kind, 4.5030949e-001_r_kind, & 1.4621608e-001_r_kind, 1.6264635e+000_r_kind, -7.9710269e+000_r_kind, -3.4492483e+000_r_kind, & -1.9521054e+000_r_kind, 9.8081595e-001_r_kind, 5.2375311e-001_r_kind, 1.9963808e-002_r_kind, & 1.7905952e-002_r_kind, -1.4191066e+000_r_kind, -7.2373867e+000_r_kind, -8.4922844e-001_r_kind, & -1.7871433e+000_r_kind, 6.7800480e-001_r_kind, 4.4923934e-001_r_kind, -1.0867891e-001_r_kind, & 1.4621608e-001_r_kind, 1.6264635e+000_r_kind, -7.9710269e+000_r_kind, -3.4492483e+000_r_kind, & -1.9521054e+000_r_kind, 9.8081595e-001_r_kind, 5.2375311e-001_r_kind, 1.9963808e-002_r_kind, & 1.7905952e-002_r_kind, -1.4191066e+000_r_kind, -7.2373867e+000_r_kind, -8.4922844e-001_r_kind, & -1.7871433e+000_r_kind, 6.7800480e-001_r_kind, 4.4923934e-001_r_kind, -1.0867891e-001_r_kind, & 3.6919810e-002_r_kind, 1.3084960e+000_r_kind, -7.9118624e+000_r_kind, -2.3499284e+000_r_kind, & -1.2343282e+000_r_kind, 1.0643171e+000_r_kind, 3.3098811e-001_r_kind, -1.2246436e-001_r_kind, & 2.0455542e-001_r_kind, -4.9476732e-002_r_kind, -1.0600542e+001_r_kind, 8.5216874e-001_r_kind, & 2.0095403e+000_r_kind, -7.0319116e-001_r_kind, -4.6245363e-001_r_kind, 1.1573103e-001_r_kind/), & (/ 8,2,5 /) ) real(r_kind), target, private, dimension(0:7,2,5) :: mhscoeff_20 = Reshape( (/ & -2.1541499e-002_r_kind, 3.0639741e-001_r_kind, -8.4143925e+000_r_kind, -3.7867212e+000_r_kind, & -5.9394364e+000_r_kind, 4.2701573e+000_r_kind, 7.0584254e+000_r_kind, 1.9181076e+000_r_kind, & -1.2439569e-001_r_kind, 2.6713184e-001_r_kind, -8.2538710e+000_r_kind, -6.7376202e-001_r_kind, & -3.0175061e+000_r_kind, 7.5293374e-001_r_kind, 6.8172354e-001_r_kind, -1.9056873e-001_r_kind, & -1.7424159e-001_r_kind, 1.0344873e-001_r_kind, -5.8175201e+000_r_kind, -1.0307833e+000_r_kind, & -5.4402261e+000_r_kind, 7.4550921e-001_r_kind, 1.3684468e+000_r_kind, -1.2857871e-001_r_kind, & 1.1077634e-001_r_kind, 1.3520555e-001_r_kind, -9.1045370e+000_r_kind, 6.7574865e-001_r_kind, & 2.4825224e-001_r_kind, -8.1358206e-001_r_kind, -1.3783551e+000_r_kind, 4.2219239e-001_r_kind, & -1.4964459e-002_r_kind, 6.2674385e-001_r_kind, -1.0447960e+001_r_kind, -4.3160647e-001_r_kind, & 1.4803721e+000_r_kind, -1.6797391e+000_r_kind, -3.4519408e+000_r_kind, -1.1398560e+000_r_kind, & 1.0782687e-001_r_kind, 1.9570646e-001_r_kind, -1.2181684e+001_r_kind, 1.7590965e-001_r_kind, & 2.5933206e+000_r_kind, -1.2687303e-001_r_kind, -4.8907807e-001_r_kind, 8.3254566e-003_r_kind, & -1.4964459e-002_r_kind, 6.2674385e-001_r_kind, -1.0447960e+001_r_kind, -4.3160647e-001_r_kind, & 1.4803721e+000_r_kind, -1.6797391e+000_r_kind, -3.4519408e+000_r_kind, -1.1398560e+000_r_kind, & 1.0782687e-001_r_kind, 1.9570646e-001_r_kind, -1.2181684e+001_r_kind, 1.7590965e-001_r_kind, & 2.5933206e+000_r_kind, -1.2687303e-001_r_kind, -4.8907807e-001_r_kind, 8.3254566e-003_r_kind, & 3.6925436e-003_r_kind, 8.2820499e-001_r_kind, -1.0636771e+001_r_kind, -3.3515611e+000_r_kind, & 2.0096780e-001_r_kind, 2.4985476e+000_r_kind, 1.1403252e+000_r_kind, 5.8821380e-002_r_kind, & -6.9040768e-002_r_kind, 6.4781201e-001_r_kind, -1.1727573e+001_r_kind, -1.3895519e-001_r_kind, & 3.8286030e+000_r_kind, -5.2255921e-002_r_kind, -1.6591578e+000_r_kind, -4.7322434e-001_r_kind/), & (/ 8,2,5 /) ) real(r_kind), target, private, dimension(0:7,2,5) :: mhscoeff_21 = Reshape( (/ & 7.6232433e-002_r_kind, 1.3942748e+000_r_kind, -9.8446512e+000_r_kind, -2.0756340e+000_r_kind, & -1.3671181e+000_r_kind, 1.6274825e+000_r_kind, 5.9784073e-001_r_kind, -3.4564933e-001_r_kind, & -3.1608678e-002_r_kind, 2.4189080e-001_r_kind, -9.2069426e+000_r_kind, -2.5967264e-001_r_kind, & -7.4456704e-001_r_kind, 3.9401662e-001_r_kind, -1.9786409e-001_r_kind, -3.8355362e-002_r_kind, & -3.6970297e-001_r_kind, 7.6490498e-001_r_kind, -5.0033422e+000_r_kind, -5.3878754e-001_r_kind, & -6.3954039e+000_r_kind, 7.5371617e-001_r_kind, 1.6196293e+000_r_kind, -1.9679393e-001_r_kind, & 8.1099616e-003_r_kind, 3.2490087e-001_r_kind, -8.8245935e+000_r_kind, -1.2075572e-001_r_kind, & -3.0251434e-001_r_kind, -3.4918312e-002_r_kind, -1.1037029e+000_r_kind, 7.6611511e-002_r_kind, & -2.0474860e-002_r_kind, 2.8463020e+000_r_kind, -9.2072830e+000_r_kind, -4.7145214e+000_r_kind, & -1.7613492e+000_r_kind, 2.9370275e+000_r_kind, 5.9442198e-001_r_kind, -4.9637201e-001_r_kind, & 1.9742247e-001_r_kind, 1.0064976e+000_r_kind, -1.1721347e+001_r_kind, -2.4134441e-001_r_kind, & 1.6217790e+000_r_kind, -1.9909312e-001_r_kind, -2.2648722e-001_r_kind, 6.2695511e-002_r_kind, & -2.0474860e-002_r_kind, 2.8463020e+000_r_kind, -9.2072830e+000_r_kind, -4.7145214e+000_r_kind, & -1.7613492e+000_r_kind, 2.9370275e+000_r_kind, 5.9442198e-001_r_kind, -4.9637201e-001_r_kind, & 1.9742247e-001_r_kind, 1.0064976e+000_r_kind, -1.1721347e+001_r_kind, -2.4134441e-001_r_kind, & 1.6217790e+000_r_kind, -1.9909312e-001_r_kind, -2.2648722e-001_r_kind, 6.2695511e-002_r_kind, & 1.2025704e-002_r_kind, 1.2674695e+000_r_kind, -1.0516845e+001_r_kind, -1.0019062e+000_r_kind, & 1.4141989e-001_r_kind, 6.2272865e-001_r_kind, 1.7704187e-001_r_kind, -1.0959747e-001_r_kind, & -4.0295798e-002_r_kind, 2.1012750e+000_r_kind, -9.8839283e+000_r_kind, -9.8200977e-001_r_kind, & 8.2772732e-001_r_kind, 5.4999664e-003_r_kind, -3.2910269e-002_r_kind, 3.2401454e-002_r_kind/) , & (/ 8,2,5 /) ) real(r_kind), target, private, dimension(0:7,2,5) :: mhscoeff_22 = Reshape( (/ & 2.5050598e-001_r_kind, -6.2730297e-002_r_kind, -8.8656340e+000_r_kind, -6.3468504e-001_r_kind, & -1.9718552e+000_r_kind, 4.9594674e-001_r_kind, 6.9558799e-001_r_kind, -8.0025569e-002_r_kind, & -3.3202425e-002_r_kind, -2.3962040e-001_r_kind, -8.4419193e+000_r_kind, 1.2665473e+000_r_kind, & -1.9571604e+000_r_kind, -1.1674975e+000_r_kind, 6.2022537e-001_r_kind, 2.3576872e-001_r_kind, & -1.7486788e-001_r_kind, 6.0032076e-001_r_kind, -5.4257035e+000_r_kind, -1.2123948e+000_r_kind, & -6.4185920e+000_r_kind, 1.2191445e+000_r_kind, 1.7030684e+000_r_kind, -2.8053701e-001_r_kind, & -3.5663303e-002_r_kind, -4.3715101e-001_r_kind, -7.2040515e+000_r_kind, -7.2689578e-002_r_kind, & -2.9207318e+000_r_kind, 3.4304473e-001_r_kind, 6.0343611e-001_r_kind, -1.3739409e-001_r_kind, & -1.8886766e-001_r_kind, 2.4629128e+000_r_kind, -6.8254299e+000_r_kind, -5.2823763e+000_r_kind, & -2.7499585e+000_r_kind, 3.0000246e+000_r_kind, 6.5303785e-001_r_kind, -4.5435145e-001_r_kind, & 5.5550888e-002_r_kind, -1.1213347e+000_r_kind, -1.0163505e+001_r_kind, -2.3168841e-001_r_kind, & 5.5254114e-001_r_kind, 1.6528653e-001_r_kind, 6.1143074e-002_r_kind, -3.9106724e-003_r_kind, & -1.8886766e-001_r_kind, 2.4629128e+000_r_kind, -6.8254299e+000_r_kind, -5.2823763e+000_r_kind, & -2.7499585e+000_r_kind, 3.0000246e+000_r_kind, 6.5303785e-001_r_kind, -4.5435145e-001_r_kind, & 5.5550888e-002_r_kind, -1.1213347e+000_r_kind, -1.0163505e+001_r_kind, -2.3168841e-001_r_kind, & 5.5254114e-001_r_kind, 1.6528653e-001_r_kind, 6.1143074e-002_r_kind, -3.9106724e-003_r_kind, & -6.4852767e-002_r_kind, 9.2157054e-001_r_kind, -9.4032297e+000_r_kind, -2.8345444e+000_r_kind, & -7.5975519e-001_r_kind, 2.1118293e+000_r_kind, 2.9422322e-001_r_kind, -3.9856154e-001_r_kind, & -6.8119340e-002_r_kind, -1.2434633e+000_r_kind, -9.6530972e+000_r_kind, 9.4927631e-002_r_kind, & 8.7302977e-001_r_kind, -1.0281616e-001_r_kind, -2.4406007e-002_r_kind, 3.3993088e-002_r_kind/) , & (/ 8,2,5 /) ) contains subroutine instrument_init(instr, satid, expansion, valid) !$$$ subprogram documentation block ! . . . . ! subprogram: instrument_init initalize instrument fields ! ! prgmmr: kleespies org: nesdis date: 2008-08-09 ! ! abstract: initialize variables required by this module. ! ! program history log: ! 2008-08-09 kleespies ! 2008-11-06 gayno - modified for gsi software standards ! 2009-03-05 Kleespies - Add airs and iasi ! 2011-09-13 gayno - pass back bad status (variable valid) ! if inputs are incorrect. use amsua ! noaa19 coefficients as default for aqua. ! 2015-03-24 ejones - Add SAPHIR ! 2016-11-11 gayno - add ATMS 3.3 ! ! input argument list: ! instr - Instrument number ! 1 = AVHRR-2 LAC/HRPT ! 2 = AVHRR-3 LAC/HRPT ! 3 = AVHRR-3 LAC/HRPT on NOAA-16 ! 4 = HIRS-2 ! 5 = HIRS-2I ! 6 = HIRS-3 NOAA-K ! 7 = HIRS-3 NOAA-L,M ! 8 = HIRS-4 ! 9 = SSU ! 10 = MSU ! 11 = AMSU-A ! 12 = AMSU-B, HSB ! 13 = MHS ! 14 = ATMS 5.2 DEG ! 15 = ATMS 2.2 DEG ! 16 = ATMS 1.1 DEG ! 17 = AIRS ! 18 = IASI ! 19 = SAPHIR ! 20 = ATMS 3.3 DEG ! satid - satellite id ! expansion - expansion factor. Must be 1.0 for accurate renderine, ! > 1.0 makes bigger ellipses, < 1.0 makes smaller ellipses. ! do not make bigger than 3.0. use 1.0 for ir and 3.0 for ! microwave. ! ! output argument list: ! valid - error status; set to false when inputs are incorrect. ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ use constants, only : pi, half, one, two implicit none ! Declare passed variables. character(len=*), intent(in ) :: satid integer(i_kind) , intent(in ) :: instr logical , intent( out) :: valid real(r_kind) , intent(in ) :: expansion ! Declare local variables. integer(i_kind) :: i, ifov real(r_kind) :: ata, cta, atf, ctf, ratio, height valid=.true. if (instr < 1 .or. instr > maxinstr) then write(6,*) "INSTRUMENT_INIT: INSTRUMENT NUMBER OF: ", instr, " IS OUT OF RANGE." valid=.false. return end if if (expansion < 0.95_r_kind .or. expansion > 3.0_r_kind) then write(6,*) "INSTRUMENT_INIT: EXPANSION FACTOR IS OUT OF RANGE." valid=.false. return endif call get_sat_height(satid, height, valid) if (.not.valid) return call fov_cleanup allocate (alongtrackangle(1:maxfov(instr))) allocate (crosstrackangle(1:maxfov(instr))) allocate (alongtrackfovsize(1:maxfov(instr))) allocate (crosstrackfovsize(1:maxfov(instr))) allocate (rmax(1:maxfov(instr))) allocate (rmin(1:maxfov(instr))) allocate (eccen(1:maxfov(instr))) do i = 1, npoly psi(i) = two*pi*float(i-1)/float(npoly-1) ! Will connect Npoly points enddo ! Precompute angles and sizes for speed. For accurate representation of fov, ! this computation should go with the height from the 1B. do ifov = 1 , maxfov(instr) call fovanglessizes(instr,height,ifov,ata,cta,atf,ctf) alongtrackangle(ifov) = ata crosstrackangle(ifov) = cta alongtrackfovsize(ifov) = atf crosstrackfovsize(ifov) = ctf enddo ! ifov rmax = half * crosstrackangle * expansion ! remember, these are semiaxes rmin = half * alongtrackangle * expansion do i = 1 , maxfov(instr) ratio = rmin(i)**2/rmax(i)**2 if(ratio > one ) ratio = one ! this takes care of some precision issues eccen(i) = sqrt(one - ratio) enddo ! set antenna power coefficients for amsu-a sensor. each satellite ! has different coefficients. i had to create separate data ! statements becuase a single data statement had too many elements ! according to the aix compiler. nullify(amsucoeff) if (instr == 11) then ! amsu-a sensor select case (trim(satid)) case('n15') amsucoeff=>amsucoeff_15 case('n16') amsucoeff=>amsucoeff_16 case('n17') amsucoeff=>amsucoeff_17 case('n18') amsucoeff=>amsucoeff_18 case('n19') amsucoeff=>amsucoeff_19 case('metop-a') amsucoeff=>amsucoeff_20 case('metop-b') amsucoeff=>amsucoeff_21 case('metop-c') amsucoeff=>amsucoeff_22 case('aqua') amsucoeff=>amsucoeff_19 case default write(6,*) "INSTRUMENT_INIT: NO AMSUA COEFFICIENTS FOR SATELLITE ",trim(satid) valid=.false. return end select endif nullify(atmscoeff) if (instr == 20) then ! atms 3.3. we don't have coefficients for atms. ! but it is similar to amsua. atmscoeff=>amsucoeff_20 endif ! set antenna power coefficients for mhs sensor. each satellite ! has different coefficients. i had to create separate data ! statements becuase a single data statement had too many elements ! according to the aix compiler. nullify(mhscoeff) if (instr == 13) then ! mhs sensor select case (trim(satid)) case('n18') mhscoeff=>mhscoeff_18 case('n19') mhscoeff=>mhscoeff_19 case('metop-a') mhscoeff=>mhscoeff_20 case('metop-b') mhscoeff=>mhscoeff_21 case('metop-c') mhscoeff=>mhscoeff_22 case default write(6,*) "INSTRUMENT_INIT: NO MHS COEFFICIENTS FOR SATELLITE ",trim(satid) valid=.false. return end select endif end subroutine instrument_init subroutine fov_cleanup !$$$ subprogram documentation block ! . . . . ! subprogram: fov cleanup deallocate module arrays ! ! prgmmr: gayno org: emc date: 2008-08-09 ! ! abstract: free up memory by deallocating module arrays ! ! program history log: ! 2008-08-09 gayno ! ! input argument list: ! n/a ! ! output argument list: ! n/a ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ implicit none if(allocated(alongtrackangle)) deallocate (alongtrackangle) if(allocated(crosstrackangle)) deallocate (crosstrackangle) if(allocated(alongtrackfovsize)) deallocate (alongtrackfovsize) if(allocated(crosstrackfovsize)) deallocate (crosstrackfovsize) if(allocated(rmax)) deallocate (rmax) if(allocated(rmin)) deallocate (rmin) if(allocated(eccen)) deallocate (eccen) end subroutine fov_cleanup subroutine fov_ellipse_crosstrk (ifov, instr, satellite_azimuth, & lat, lon, elats, elons ) !$$$ subprogram documentation block ! . . . . ! subprogram: fov_ellipse_crosstrk computes fov ellipses ! ! prgmmr: kleespies org: nesdis date: 2007-08-09 ! ! abstract: computes FOV ellipses in latitude/longitude coordinates ! for a scan line of a cross track scanning instrument. ! ! program history log: ! 2007-08-09 kleespies ! 2007-08-13 gayno - modified for gsi software standards ! 2009-03-05 kleespies - Add AMSU-MHS Antenna Patterns ! 2009-09-21 kleespies - Add AIRS-IASI ! 2016-11-11 gayno - add ATMS 3.3 ! ! input argument list: ! ifov - fov number min = 1, max as below ! 2048 = AVHRR-2 LAC/HRPT ! 2048 = AVHRR-3 LAC/HRPT ! 2048 = AVHRR-3 LAC/HRPT on NOAA-16 ! 56 = HIRS-2 ! 56 = HIRS-2I ! 56 = HIRS-3 ! 56 = HIRS-4 ! 8 = SSU ! 11 = MSU ! 13 = MHS ! 30 = AMSU-A ! 90 = AMSU-B ! 90 = AIRS ! 96 = ATMS ! 120 = IASI ! instr - Instrument number ! 1 = AVHRR-2 LAC/HRPT ! 2 = AVHRR-3 LAC/HRPT ! 3 = AVHRR-3 LAC/HRPT on NOAA-16 ! 4 = HIRS-2 ! 5 = HIRS-2I ! 6 = HIRS-3 NOAA-K ! 7 = HIRS-3 NOAA-L,M ! 8 = HIRS-4 ! 9 = SSU ! 10 = MSU ! 11 = AMSU-A ! 12 = AMSU-B, HSB ! 13 = MHS ! 14 = ATMS 5.2 DEG ! 15 = ATMS 2.2 DEG ! 16 = ATMS 1.1 DEG ! 17 = AIRS ! 18 = IASI ! 20 = ATMS 3.3 DEG ! satellite_azimuth - satellite azimuth angle ! lat - latitude of center of fov ! lon - longitude of center of fov ! ! output argument list: ! elats - ellipse latitudes centered about lat,lon ! elons - ellipse longitudes centered about lat,lon ! ! remarks: ! ! Set npoly for the order of the polygon that represents the ellipse. ! ! elats and elon must be dimensioned to npoly in the calling routine. ! ! There are several engineering checks to handle things like ! arcsin( x>1.0 or x<-1.0). These things can happen because POES navigation ! uses an oblate spheroid earth, while here we are using a spherical ! earth, so there is a small inconsistency in computing arc angles. ! ! No provisions are made for spacecraft roll-pitch-yaw errors, ! which are presumed to be small. (SC attitude is usually held to within 0.1 deg) ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ use constants, only : rad2deg, one implicit none ! Declare passed variables integer(i_kind), intent(in ) :: ifov integer(i_kind), intent(in ) :: instr real(r_kind) , intent(in ) :: satellite_azimuth real(r_kind) , intent(in ) :: lat real(r_kind) , intent(in ) :: lon real(r_kind) , intent( out) :: elats(npoly) real(r_kind) , intent( out) :: elons(npoly) ! Declare local variables integer(i_kind):: i ! loop counters integer(i_kind):: fov real(r_kind):: pos_ang ! rotation angle of the ellipse real(r_kind):: psip(npoly), r(npoly), cosc(npoly), c(npoly), sinb(npoly), b(npoly) fov = ifov if(instr == 18) then ! iasi fov = (ifov-1)/4 + 1 endif pos_ang = satellite_azimuth if(pos_ang > 180._r_kind) pos_ang = pos_ang-360._r_kind if(pos_ang < -180._r_kind) pos_ang = 360._r_kind + pos_ang psip = psi + pos_ang/rad2deg r = rmax(fov) * sqrt( (one - eccen(fov)**2)/(one - eccen(fov)**2 *cos(psi)**2) ) cosc = cos((90._r_kind-lat)/rad2deg)*cos(r/rad2deg) + & sin((90._r_kind-lat)/rad2deg)*sin(r/rad2deg)*cos(psip) c = acos(cosc)*rad2deg elats(1:npoly) = 90._r_kind - c sinb = sin(r/rad2deg)*sin(psip)/sin(c/rad2deg) do i = 1 , npoly ! handle numeric imprecision if(sinb(i) > one) sinb(i) = one if(sinb(i) < -(one)) sinb(i) = -(one) enddo b = asin(sinb)*rad2deg elons(1:npoly) = lon + b end subroutine fov_ellipse_crosstrk subroutine fovanglessizes(instr,height,fov,alongtrackangle,crosstrackangle,& alongtrackfovsize,crosstrackfovsize) !$$$ subprogram documentation block ! . . . . ! subprogram: fovanglessizes return cross and along track angles/sizes. ! ! prgmmr: kleespies org: nesdis date: 2007-08-09 ! ! abstract: computes the cross track and along track angles of an ! instrument FOV as viewed from the center of the earth, and ! cross track and along track FOV size in km. presumes a spherical earth. ! ! program history log: ! 2007-08-09 kleespies ! 2007-08-13 gayno - modified for gsi software standards ! 2009-09-21 kleespies - add AIRS and IASI ! 2015-03-24 ejones - add SAPHIR ! 2016-11-11 gayno - add ATMS 3.3 ! ! input argument list: ! instr - Instrument number ! 1 = AVHRR-2 LAC/HRPT ! 2 = AVHRR-3 LAC/HRPT ! 3 = AVHRR-3 LAC/HRPT on NOAA-16 ! 4 = HIRS-2 ! 5 = HIRS-2I ! 6 = HIRS-3 NOAA-K ! 7 = HIRS-3 NOAA-L,M ! 8 = HIRS-4 ! 9 = SSU ! 10 = MSU ! 11 = AMSU-A ! 12 = AMSU-B, HSB ! 13 = MHS ! 14 = ATMS 5.2 DEG ! 15 = ATMS 2.2 DEG ! 16 = ATMS 1.1 DEG ! 17 = AIRS ! 18 = IASI ! 19 = SAPHIR ! 20 = ATMS 3.3 DEG ! height - height of satellite in km ! fov - fov number min = 1, max as below ! 2048 = AVHRR-2 LAC/HRPT ! 2048 = AVHRR-3 LAC/HRPT ! 2048 = AVHRR-3 LAC/HRPT on NOAA-16 ! 56 = HIRS-2 ! 56 = HIRS-2I ! 56 = HIRS-3 ! 56 = HIRS-4 ! 8 = SSU ! 11 = MSU ! 90 = MHS ! 30 = AMSU-A ! 90 = AMSU-B ! 96 = ATMS ! 90 = AIRS ! 120 = IASI ! 130 = SAPHIR ! ! output argument list: ! AlongTrackAngle - along track angle of a fov ! CrossTrackAngle - cross track angle of a fov ! AlongTrackFOVSize - along track fov size in km ! CrossTrackFOVSize - cross track fov size in km ! ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ use constants, only : pi, rad2deg, zero, half, one, two, three implicit none ! Declare passed variables integer(i_kind), intent(in ) :: fov, instr real(r_kind) , intent(in ) :: height real(r_kind) , intent( out) :: alongtrackangle, crosstrackangle real(r_kind) , intent( out) :: alongtrackfovsize, crosstrackfovsize ! Declare local parameters real(r_kind), parameter:: radius = 6371.22_r_kind ! assymetry degrees (+-) that the scanner is assymetric about nadir. ! This thusfar is important only for NOAA-16&17 HIRS. ! degscan FOV step angle ! fovangle angular dimension of the FOV ! halfscan degrees from nadir to the center of outer earth FOV ! i think AIRS is same as AMSU-B... anyway, have put MHS into old AIRS slot ! noaa-n parameters presently set to klm. may need different numbers for ! AQUA AMSU-A and HSB amsub=1.1, mhs=10/9. real(r_kind) , dimension(maxinstr) :: degscan = (/ 0.0541e0_r_kind, 0.0541e0_r_kind, 0.053981436e0_r_kind, & 1.8e0_r_kind, 1.8e0_r_kind, 1.8e0_r_kind, & 1.8e0_r_kind, 1.8e0_r_kind, 10.0e0_r_kind, & 9.4737e0_r_kind, 3.e0_r_kind+one/three, & 1.1e0_r_kind , 10.e0_r_kind/9.e0_r_kind , 1.1e0_r_kind, & 1.1e0_r_kind, 1.1e0_r_kind, & 1.1e0_r_kind, 3.e0_r_kind+one/three, 0.6660465_r_kind, & 1.1e0_r_kind /) real(r_kind) , dimension(maxinstr) :: fovangle = (/ 0.0745_r_kind, 0.0745_r_kind, 0.0745_r_kind, & 1.22_r_kind, 1.40_r_kind, 1.40_r_kind, & 1.40_r_kind, 0.70_r_kind, 10.0_r_kind, & 7.5_r_kind, 3.3_r_kind, 1.1_r_kind, & 1.1_r_kind , 5.2_r_kind, 2.2_r_kind, 1.1_r_kind, & 1.1_r_kind, 0.839383_r_kind, 0.6660465_r_kind, & 3.3_r_kind /) real(r_kind) , dimension(maxinstr) :: halfscan = (/ 55.37_r_kind, 55.37_r_kind, 55.25_r_kind, 49.5_r_kind, & 49.5_r_kind, 49.5_r_kind, 49.5_r_kind, 49.5_r_kind, & 35.0_r_kind, 47.3685_r_kind, & 48._r_kind+one/three, 48.95_r_kind, & 48.95_r_kind, 52.73_r_kind, 52.73_r_kind, 52.73_r_kind, & 44.5_r_kind*10.0_r_kind/9.0_r_kind, & 48._r_kind+one/three, 42.96_r_kind, 52.73_r_kind /) real(r_kind) , dimension(maxinstr) :: assymetry = (/ zero, zero, zero, zero, & zero, zero, -1.8_r_kind, zero, & zero, zero, zero, zero, & zero, zero, zero, zero, & zero, zero, zero, zero /) ! declare local variables real(r_kind) nadirangle, nadirangle_m, nadirangle_p real(r_kind) compzacenter, compza_m, compza_p, distancetofov ! initialize to bad value alongtrackangle = -999._r_kind crosstrackangle = -999._r_kind !Nadir angles of center and crosstrack extremes of fov nadirangle = halfscan(instr) - (fov-1)*degscan(instr) + assymetry(instr) nadirangle_m = nadirangle - fovangle(instr)*half nadirangle_p = nadirangle + fovangle(Instr)*half !Complement of zenith angle for center and crosstrack extremes of fov compzacenter = 180._r_kind - asin((radius+height)/radius * sin(nadirangle /rad2deg))*rad2deg compza_m = 180._r_kind - asin((radius+height)/radius * sin(nadirangle_m/rad2deg))*rad2deg compza_p = 180._r_kind - asin((radius+height)/radius * sin(nadirangle_p/rad2deg))*rad2deg !cross track angle of the fov as viewed from center of the earth crosstrackangle = abs(nadirangle_p + compza_p - nadirangle_m - compza_m) !cross track fov size in km crosstrackfovsize = abs(crosstrackangle*two*pi*radius/360._r_kind) !distance from satellite to the center of the fov in km distancetofov = (radius+height) * & sin( (180._r_kind-nadirangle-compzacenter)/rad2deg)/sin((compzacenter)/rad2deg) !along track fov size in km ! the following is an approximation, but it is close. It underestimates the FOV by a smidge alongtrackfovsize = two*distancetofov*tan(fovangle(instr)*half/rad2deg) !along track angle of the fov as viewed from center of the earth alongtrackangle = 360._r_kind * alongtrackfovsize/(two*pi*radius) end subroutine fovanglessizes subroutine get_sat_height(satid, height, valid) !$$$ subprogram documentation block ! . . . . ! subprogram: get_sat_height return height of satellite in km ! prgmmr: gayno org: np23 date: 2007-08-09 ! ! abstract: return height of desired satellite in km. ! information from: ! 1) www.itc.nl/research/products/sensordb ! 2) a nasa/noaa booklet entitled, "Advanced TIROS-N ! (ATN): NOAA J" TL798.M4A38 (1993) ! ! program history log: ! 2007-08-09 gayno ! 2009-12-20 gayno - add metop-b/c ! 2011-09-13 gayno - pass back bad status for undefined ! satellites ! 2015-03-24 ejones - add megha tropiques ! ! input argument list: ! satid - satellite id ! ! output argument list: ! height - height of satellite in km ! valid - status flag, set to false for ! undefined satellites ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ implicit none ! Declare passed variables character(len=*), intent(in ) :: satid logical, intent( out) :: valid real(r_kind) , intent( out) :: height valid=.true. select case (trim(satid)) case('tirosn') height=870._r_kind case('dmsp', 'f13', 'f14', 'f15', 'f16', 'f17') height=830._r_kind case('trmm') height=350._r_kind case('n05') height=1511._r_kind case('n07', 'n09', 'n11', 'n14', 'n16', 'n18', 'n19') height=870._r_kind case('n06', 'n08', 'n10', 'n12', 'n15', 'n17') height=833._r_kind case('aura', 'aqua') height=705._r_kind case('metop-a', 'metop-b', 'metop-c') height=817._r_kind case('meghat') height=866._r_kind case('npp') height=840._r_kind case('n20') height=840._r_kind case default write(6,*) 'GET_SAT_HEIGHT: ERROR, unrecognized satellite id: ', trim(satid) valid=.false. return end select end subroutine get_sat_height subroutine inside_fov_crosstrk(instr, ifov, satellite_azimuth, & lat, lon, testlat, testlon, & expansion, ichan, inside) !$$$ subprogram documentation block ! . . . . ! subprogram: inside_fov_crosstrk determine antenna power ! ! prgmmr: kleespies org: nesdis date: 2008-08-09 ! ! abstract: determine if a point is inside a fov. if it is, return ! the relative antenna power. ! ! program history log: ! 2008-08-09 kleespies ! 2008-11-06 gayno - modified for gsi software standards ! 2009-03-05 kleespies - modify to include all amsu-a and mhs instruments and all channels ! 2009-09-21 kleespies - add AIRS and IASI ! 2016-11-11 gayno - add ATMS 3.3 ! ! input argument list: ! instr - instrument number ! 1 = AVHRR-2 LAC/HRPT ! 2 = AVHRR-3 LAC/HRPT ! 3 = AVHRR-3 LAC/HRPT on NOAA-16 ! 4 = HIRS-2 ! 5 = HIRS-2I ! 6 = HIRS-3 NOAA-K ! 7 = HIRS-3 NOAA-L,M ! 8 = HIRS-4 ! 9 = SSU ! 10 = MSU ! 11 = AMSU-A ! 12 = AMSU-B, HSB ! 13 = MHS ! 14 = ATMS 5.2 DEG ! 15 = ATMS 2.2 DEG ! 16 = ATMS 1.1 DEG ! 17 = AIRS ! 18 - IASI ! 19 = SAPHIR ! 20 = ATMS 3.3 DEG ! ifov - fov number ! satellite_azimuth - satellite azimuth angle (degrees) ! lat - latitude of fov center ! lon - longitude of fov center ! testlat - latitude of point to be tested ! testlon - longitude of point to be tested ! expansion - expansion factor. ! 1.0 for 50% power ! 2.0 for 95% power ! 3.0 for 99% power ! > 1.0 makes bigger ellipses, < 1.0 makes smaller ellipses. ! do not exceed 3.0. ! in general try to use 3.0 for microwave if you can afford it. ! use 1.0 for infrared ! ichan - channel number for microwave instruments, 1-15 for amsua, ! 1-5 for mhs, na for all else ! ! output argument list: ! inside - 0.0-1.0 relative antenna power if [testLat,testLon] is ! inside the FOV for amsu and mhs, 0 if outside. For other ! instruments, is 1.0 if inside or 0.0 if outside. ! ! remarks: ! there are several engineering checks to handle things like arcsin( x>1.0 or x<-1.0) ! these things can happen because POES navigation uses an oblate spheroid earth, while here we are ! using a spherical earth, so there is a small inconsistency in computing arc angles. ! ! no provisions are made for spacecraft roll-pitch-yaw errors, which are presumed to be small ! (SC attitude is usually held to within 0.1 deg) ! ! attributes: ! language: f90 ! machine: ibm RS/6000 SP ! !$$$ use constants, only : one, zero, half, two, pi, deg2rad, rad2deg, one_tenth implicit none ! Declare passed variables. integer(i_kind),intent(in ) :: instr integer(i_kind),intent(in ) :: ifov integer(i_kind),intent(in ) :: ichan real(r_kind) ,intent(in ) :: lat real(r_kind) ,intent(in ) :: lon real(r_kind) ,intent(in ) :: testlat real(r_kind) ,intent(in ) :: testlon real(r_kind) ,intent(in ) :: expansion real(r_kind) ,intent(in ) :: satellite_azimuth real(r_kind) ,intent( out) :: inside ! Declare local parameters. real(r_kind), parameter :: r1 = one ! equatorial radius. work in angular distance, not km (otherwise r1=6371) real(r_kind), parameter :: r2 = r1 ! assume spherical earth (otherwise r2 = polar radius) real(r_kind), dimension(maxinstr) :: fovangle = (/ 0.0745_r_kind, 0.0745_r_kind, 0.0745_r_kind, & 1.22_r_kind, 1.40_r_kind, 1.40_r_kind, & 1.40_r_kind, 0.70_r_kind, 10.0_r_kind, & 7.5_r_kind, 3.3_r_kind, 1.1_r_kind, & 1.1_r_kind, 5.2_r_kind, 2.2_r_kind, 1.1_r_kind, & 1.1_r_kind, 0.839383_r_kind, 0.6660465_r_kind, & 3.3_r_kind /) ! Declare local variables. integer(i_kind) :: fov real(r_kind) :: dellon ! longitude difference from fov center to test location real(r_kind) :: dellat ! latitude difference from fov center to test location real(r_kind) :: fovanglesize ! these two are what gets compared to determine if test location is within fov. real(r_kind) :: d ! angular distance from fov center to test location real(r_kind) :: r ! angular distance from fov center to ellipse along d real(r_kind) :: psi ! angle from 3:00 on fov to line from fov center to test location real(r_kind) :: psip ! angle from latitude parallel through fov center to ! line from fov center to test location ! these are the flat earth variables real(r_kind) :: distance_north ! north angular distance from fov center to test location real(r_kind) :: distance_east ! east angular distance from fov center to test location real(r_kind) :: bearing_to_test ! same as psip real(r_kind) :: bearing_to_test_deg ! in degrees real(r_kind) :: x,y,px,py,p,rat real(r_kind) :: sataz ! satellite azimuth, used in computing psi fov = ifov if(instr == 18) then fov = (ifov-1)/4 + 1 endif ! Get satellite az where we want it. ! 1st, convert +- to 0-360 sataz = satellite_azimuth if(sataz < zero) sataz = 360.0_r_kind + sataz ! 2nd, Shift rotation direction sataz = mod((450.0_r_kind-sataz),360.0_r_kind) ! delta lat and lon dellat = (testlat - lat)*deg2rad dellon = testlon - lon if (dellon > 180._r_kind) dellon = dellon - 360._r_kind dellon = dellon*deg2rad ! distance north and east in degrees distance_north = r1*dellat distance_east = r2*cos(lat*deg2rad)*dellon ! angle to the test point bearing_to_test = mod(atan2(distance_north,distance_east),2*pi ) bearing_to_test_deg = bearing_to_test*rad2deg ! convert to degrees ! this is the arc distance to the test point d=two*asin(sqrt((sin(dellat/two))**2 + & cos(testlat*deg2rad)*cos(lat*deg2rad)*(sin(dellon/two))**2)) d = d*rad2deg ! convert to degrees psip = bearing_to_test_deg psi = (psip - sataz) psi = psi*deg2rad ! convert to radians ! r is the angular distance from the fov center to the edge of the ellipse in degrees r = rmax(fov) * sqrt( (one - eccen(fov)**2)/(one - eccen(fov)**2 *cos(psi)**2) ) inside = zero if (d one) p = one inside = p endif ! amsuA if(instr == 20) then ! ATMS 3.3 DEG. Use amsua channel 15 coeffs, which is sfc sensitive. fovanglesize = fovangle(instr) rat = d / r * expansion * fovanglesize * half x = rat * cos(psi) y = rat * sin(psi) px = atmscoeff(0,1,15) + atmscoeff(1,1,15)*x + atmscoeff(2,1,15)*x**2 & + atmscoeff(3,1,15)*x**3 + atmscoeff(4,1,15)*x**4 & + atmscoeff(5,1,15)*x**5 + atmscoeff(6,1,15)*x**6 & + atmscoeff(7,1,15)*x**7 py = atmscoeff(0,2,15) + atmscoeff(1,2,15)*y + atmscoeff(2,2,15)*y**2 & + atmscoeff(3,2,15)*y**3 + atmscoeff(4,2,15)*y**4 & + atmscoeff(5,2,15)*y**5 + atmscoeff(6,2,15)*y**6 & + atmscoeff(7,2,15)*y**7 p = -(px+py) ! power in dB (positive) ! convert to fraction of max power p = 10._r_kind**(-p*one_tenth) ! convert to fraction of max power if(p > one) p = one inside = p endif ! atms if(instr == 13) then ! mhs fovanglesize = fovangle(instr) rat = d / r * expansion * fovanglesize * half x = rat * cos(psi) y = rat * sin(psi) px = mhscoeff(0,1,ichan) + mhscoeff(1,1,ichan)*x + mhscoeff(2,1,ichan)*x**2 & + mhscoeff(3,1,ichan)*x**3 + mhscoeff(4,1,ichan)*x**4 & + mhscoeff(5,1,ichan)*x**5 + mhscoeff(6,1,ichan)*x**6 & + mhscoeff(7,1,ichan)*x**7 py = mhscoeff(0,2,ichan) + mhscoeff(1,2,ichan)*y + mhscoeff(2,2,ichan)*y**2 & + mhscoeff(3,2,ichan)*y**3 + mhscoeff(4,2,ichan)*y**4 & + mhscoeff(5,2,ichan)*y**5 + mhscoeff(6,2,ichan)*y**6 & + mhscoeff(7,2,ichan)*y**7 p = -(px+py) ! power in dB (positive) p = 10._r_kind**(-p*one_tenth) ! convert to fraction of max power if(p > one) then p = one endif inside = p endif ! mhs endif ! is (d maxfov(instr) ) then write(6,*) "FOV_CHECK: ERROR, FOV GREATER THAN ",maxfov(instr), "FOR INSTRUMENT ", instr valid=.false. return endif if (instr == 11) then if (ichan < 1 .or. ichan > 15) then write(6,*) "FOV_CHECK: ERROR, invalid amsua channel number" valid=.false. return endif else if (instr == 13) then if (ichan < 1 .or. ichan > 5) then write(6,*) "FOV_CHECK: ERROR, invalid mhs channel number" valid=.false. return endif endif return end subroutine fov_check end module calc_fov_crosstrk