subroutine DBMP1 (X, AMPL, THETA)
c Copyright (c) 1996 California Institute of Technology, Pasadena, CA.
c ALL RIGHTS RESERVED.
c Based on Government Sponsored Research NAS7-03001.
c>> 2002-03-29 DBMP1 Krogh  Very minor clean up of code.
c>> 1996-04-27 DBMP1 Krogh  Changes to use .C. and C%%.
c>> 1996-03-30 DBMP1 Krogh  Added external statement.
C>> 1995-11-28 DBMP1 Krogh  Changes to simplify conversion to C.
C>> 1995-11-03 DBMP1 Krogh  Removed blanks in numbers for C conversion.
C>> 1994-11-11 DBMP1 Krogh   Declared all vars.
C>> 1994-10-20 DBMP1 Krogh  Changes to use M77CON
C>> 1991-01-14 DBMP1 CLL Changed to use generic name SQRT
C>> 1990-11-29 CLL Changed subroutine name to DBMP1
C>> 1985-08-02 D9B1MP Lawson  Initial code.
C JULY 1977 EDITION.  W. FULLERTON, C3, LOS ALAMOS SCIENTIFIC LAB.
C C.L.LAWSON & S.CHAN, JPL, 1984 FEB ADAPTED TO JPL MATH77 LIBRARY.
C
C EVALUATE THE MODULUS AND PHASE FOR THE BESSEL J1 AND Y1 FUNCTIONS.
C
C     ------------------------------------------------------------------
c--D replaces "?": ?BMP1, ?INITS, ?CSEVL, ?ERM1, ?ERV1
C     ------------------------------------------------------------------
      external D1MACH, DCSEVL
      integer NBM1, NBT12, NBM12, NBTH1
      double precision X, AMPL, THETA, BM1CS(37), BT12CS(39), ETA
      double precision BM12CS(40), BTH1CS(44), XMAX, PI4
      double precision Z, D1MACH, DCSEVL
C
      save ETA, NBM1, NBT12, NBM12, NBTH1, XMAX
C
C SERIES FOR BM1        ON THE INTERVAL  1.56250D-02 TO  6.25000D-02
C                                        WITH WEIGHTED ERROR   4.91D-32
C                                         LOG WEIGHTED ERROR  31.31
C                               SIGNIFICANT FIGURES REQUIRED  30.04
C                                    DECIMAL PLACES REQUIRED  32.09
C
c++ Save data by elements if ~.C.
      data BM1CS(1) /   +.1069845452618063014969985308538D+0 /
      data BM1CS(2) /   +.3274915039715964900729055143445D-2 /
      data BM1CS(3) /   -.2987783266831698592030445777938D-4 /
      data BM1CS(4) /   +.8331237177991974531393222669023D-6 /
      data BM1CS(5) /   -.4112665690302007304896381725498D-7 /
      data BM1CS(6) /   +.2855344228789215220719757663161D-8 /
      data BM1CS(7) /   -.2485408305415623878060026596055D-9 /
      data BM1CS(8) /   +.2543393338072582442742484397174D-10 /
      data BM1CS(9) /   -.2941045772822967523489750827909D-11 /
      data BM1CS(10) /  +.3743392025493903309265056153626D-12 /
      data BM1CS(11) /  -.5149118293821167218720548243527D-13 /
      data BM1CS(12) /  +.7552535949865143908034040764199D-14 /
      data BM1CS(13) /  -.1169409706828846444166290622464D-14 /
      data BM1CS(14) /  +.1896562449434791571721824605060D-15 /
      data BM1CS(15) /  -.3201955368693286420664775316394D-16 /
      data BM1CS(16) /  +.5599548399316204114484169905493D-17 /
      data BM1CS(17) /  -.1010215894730432443119390444544D-17 /
      data BM1CS(18) /  +.1873844985727562983302042719573D-18 /
      data BM1CS(19) /  -.3563537470328580219274301439999D-19 /
      data BM1CS(20) /  +.6931283819971238330422763519999D-20 /
      data BM1CS(21) /  -.1376059453406500152251408930133D-20 /
      data BM1CS(22) /  +.2783430784107080220599779327999D-21 /
      data BM1CS(23) /  -.5727595364320561689348669439999D-22 /
      data BM1CS(24) /  +.1197361445918892672535756799999D-22 /
      data BM1CS(25) /  -.2539928509891871976641440426666D-23 /
      data BM1CS(26) /  +.5461378289657295973069619199999D-24 /
      data BM1CS(27) /  -.1189211341773320288986289493333D-24 /
      data BM1CS(28) /  +.2620150977340081594957824000000D-25 /
      data BM1CS(29) /  -.5836810774255685901920938666666D-26 /
      data BM1CS(30) /  +.1313743500080595773423615999999D-26 /
      data BM1CS(31) /  -.2985814622510380355332778666666D-27 /
      data BM1CS(32) /  +.6848390471334604937625599999999D-28 /
      data BM1CS(33) /  -.1584401568222476721192960000000D-28 /
      data BM1CS(34) /  +.3695641006570938054301013333333D-29 /
      data BM1CS(35) /  -.8687115921144668243012266666666D-30 /
      data BM1CS(36) /  +.2057080846158763462929066666666D-30 /
      data BM1CS(37) /  -.4905225761116225518523733333333D-31 /
C
C SERIES FOR BT12       ON THE INTERVAL  1.56250D-02 TO  6.25000D-02
C                                        WITH WEIGHTED ERROR   3.33D-32
C                                         LOG WEIGHTED ERROR  31.48
C                               SIGNIFICANT FIGURES REQUIRED  31.05
C                                    DECIMAL PLACES REQUIRED  32.27
C
c++ Save data by elements if ~.C.
      data BT12CS(1) /   +.73823860128742974662620839792764D+0 /
      data BT12CS(2) /   -.33361113174483906384470147681189D-2 /
      data BT12CS(3) /   +.61463454888046964698514899420186D-4 /
      data BT12CS(4) /   -.24024585161602374264977635469568D-5 /
      data BT12CS(5) /   +.14663555577509746153210591997204D-6 /
      data BT12CS(6) /   -.11841917305589180567005147504983D-7 /
      data BT12CS(7) /   +.11574198963919197052125466303055D-8 /
      data BT12CS(8) /   -.13001161129439187449366007794571D-9 /
      data BT12CS(9) /   +.16245391141361731937742166273667D-10 /
      data BT12CS(10) /  -.22089636821403188752155441770128D-11 /
      data BT12CS(11) /  +.32180304258553177090474358653778D-12 /
      data BT12CS(12) /  -.49653147932768480785552021135381D-13 /
      data BT12CS(13) /  +.80438900432847825985558882639317D-14 /
      data BT12CS(14) /  -.13589121310161291384694712682282D-14 /
      data BT12CS(15) /  +.23810504397147214869676529605973D-15 /
      data BT12CS(16) /  -.43081466363849106724471241420799D-16 /
      data BT12CS(17) /  +.80202544032771002434993512550400D-17 /
      data BT12CS(18) /  -.15316310642462311864230027468799D-17 /
      data BT12CS(19) /  +.29928606352715568924073040554666D-18 /
      data BT12CS(20) /  -.59709964658085443393815636650666D-19 /
      data BT12CS(21) /  +.12140289669415185024160852650666D-19 /
      data BT12CS(22) /  -.25115114696612948901006977706666D-20 /
      data BT12CS(23) /  +.52790567170328744850738380799999D-21 /
      data BT12CS(24) /  -.11260509227550498324361161386666D-21 /
      data BT12CS(25) /  +.24348277359576326659663462400000D-22 /
      data BT12CS(26) /  -.53317261236931800130038442666666D-23 /
      data BT12CS(27) /  +.11813615059707121039205990399999D-23 /
      data BT12CS(28) /  -.26465368283353523514856789333333D-24 /
      data BT12CS(29) /  +.59903394041361503945577813333333D-25 /
      data BT12CS(30) /  -.13690854630829503109136383999999D-25 /
      data BT12CS(31) /  +.31576790154380228326413653333333D-26 /
      data BT12CS(32) /  -.73457915082084356491400533333333D-27 /
      data BT12CS(33) /  +.17228081480722747930705920000000D-27 /
      data BT12CS(34) /  -.40716907961286507941068800000000D-28 /
      data BT12CS(35) /  +.96934745136779622700373333333333D-29 /
      data BT12CS(36) /  -.23237636337765716765354666666666D-29 /
      data BT12CS(37) /  +.56074510673522029406890666666666D-30 /
      data BT12CS(38) /  -.13616465391539005860522666666666D-30 /
      data BT12CS(39) /  +.33263109233894654388906666666666D-31 /
C
C SERIES FOR BM12       ON THE INTERVAL  0.          TO  1.56250D-02
C                                        WITH WEIGHTED ERROR   5.01D-32
C                                         LOG WEIGHTED ERROR  31.30
C                               SIGNIFICANT FIGURES REQUIRED  29.99
C                                    DECIMAL PLACES REQUIRED  32.10
C
c++ Save data by elements if ~.C.
      data BM12CS(1) /   +.9807979156233050027272093546937D-1 /
      data BM12CS(2) /   +.1150961189504685306175483484602D-2 /
      data BM12CS(3) /   -.4312482164338205409889358097732D-5 /
      data BM12CS(4) /   +.5951839610088816307813029801832D-7 /
      data BM12CS(5) /   -.1704844019826909857400701586478D-8 /
      data BM12CS(6) /   +.7798265413611109508658173827401D-10 /
      data BM12CS(7) /   -.4958986126766415809491754951865D-11 /
      data BM12CS(8) /   +.4038432416421141516838202265144D-12 /
      data BM12CS(9) /   -.3993046163725175445765483846645D-13 /
      data BM12CS(10) /  +.4619886183118966494313342432775D-14 /
      data BM12CS(11) /  -.6089208019095383301345472619333D-15 /
      data BM12CS(12) /  +.8960930916433876482157048041249D-16 /
      data BM12CS(13) /  -.1449629423942023122916518918925D-16 /
      data BM12CS(14) /  +.2546463158537776056165149648068D-17 /
      data BM12CS(15) /  -.4809472874647836444259263718620D-18 /
      data BM12CS(16) /  +.9687684668292599049087275839124D-19 /
      data BM12CS(17) /  -.2067213372277966023245038117551D-19 /
      data BM12CS(18) /  +.4646651559150384731802767809590D-20 /
      data BM12CS(19) /  -.1094966128848334138241351328339D-20 /
      data BM12CS(20) /  +.2693892797288682860905707612785D-21 /
      data BM12CS(21) /  -.6894992910930374477818970026857D-22 /
      data BM12CS(22) /  +.1830268262752062909890668554740D-22 /
      data BM12CS(23) /  -.5025064246351916428156113553224D-23 /
      data BM12CS(24) /  +.1423545194454806039631693634194D-23 /
      data BM12CS(25) /  -.4152191203616450388068886769801D-24 /
      data BM12CS(26) /  +.1244609201503979325882330076547D-24 /
      data BM12CS(27) /  -.3827336370569304299431918661286D-25 /
      data BM12CS(28) /  +.1205591357815617535374723981835D-25 /
      data BM12CS(29) /  -.3884536246376488076431859361124D-26 /
      data BM12CS(30) /  +.1278689528720409721904895283461D-26 /
      data BM12CS(31) /  -.4295146689447946272061936915912D-27 /
      data BM12CS(32) /  +.1470689117829070886456802707983D-27 /
      data BM12CS(33) /  -.5128315665106073128180374017796D-28 /
      data BM12CS(34) /  +.1819509585471169385481437373286D-28 /
      data BM12CS(35) /  -.6563031314841980867618635050373D-29 /
      data BM12CS(36) /  +.2404898976919960653198914875834D-29 /
      data BM12CS(37) /  -.8945966744690612473234958242979D-30 /
      data BM12CS(38) /  +.3376085160657231026637148978240D-30 /
      data BM12CS(39) /  -.1291791454620656360913099916966D-30 /
      data BM12CS(40) /  +.5008634462958810520684951501254D-31 /
C
C SERIES FOR BTH1       ON THE INTERVAL  0.          TO  1.56250D-02
C                                        WITH WEIGHTED ERROR   2.82D-32
C                                         LOG WEIGHTED ERROR  31.55
C                               SIGNIFICANT FIGURES REQUIRED  31.12
C                                    DECIMAL PLACES REQUIRED  32.37
C
c++ Save data by elements if ~.C.
      data BTH1CS(1) /   +.74749957203587276055443483969695D+0 /
      data BTH1CS(2) /   -.12400777144651711252545777541384D-2 /
      data BTH1CS(3) /   +.99252442404424527376641497689592D-5 /
      data BTH1CS(4) /   -.20303690737159711052419375375608D-6 /
      data BTH1CS(5) /   +.75359617705690885712184017583629D-8 /
      data BTH1CS(6) /   -.41661612715343550107630023856228D-9 /
      data BTH1CS(7) /   +.30701618070834890481245102091216D-10 /
      data BTH1CS(8) /   -.28178499637605213992324008883924D-11 /
      data BTH1CS(9) /   +.30790696739040295476028146821647D-12 /
      data BTH1CS(10) /  -.38803300262803434112787347554781D-13 /
      data BTH1CS(11) /  +.55096039608630904934561726208562D-14 /
      data BTH1CS(12) /  -.86590060768383779940103398953994D-15 /
      data BTH1CS(13) /  +.14856049141536749003423689060683D-15 /
      data BTH1CS(14) /  -.27519529815904085805371212125009D-16 /
      data BTH1CS(15) /  +.54550796090481089625036223640923D-17 /
      data BTH1CS(16) /  -.11486534501983642749543631027177D-17 /
      data BTH1CS(17) /  +.25535213377973900223199052533522D-18 /
      data BTH1CS(18) /  -.59621490197413450395768287907849D-19 /
      data BTH1CS(19) /  +.14556622902372718620288302005833D-19 /
      data BTH1CS(20) /  -.37022185422450538201579776019593D-20 /
      data BTH1CS(21) /  +.97763074125345357664168434517924D-21 /
      data BTH1CS(22) /  -.26726821639668488468723775393052D-21 /
      data BTH1CS(23) /  +.75453300384983271794038190655764D-22 /
      data BTH1CS(24) /  -.21947899919802744897892383371647D-22 /
      data BTH1CS(25) /  +.65648394623955262178906999817493D-23 /
      data BTH1CS(26) /  -.20155604298370207570784076869519D-23 /
      data BTH1CS(27) /  +.63417768556776143492144667185670D-24 /
      data BTH1CS(28) /  -.20419277885337895634813769955591D-24 /
      data BTH1CS(29) /  +.67191464220720567486658980018551D-25 /
      data BTH1CS(30) /  -.22569079110207573595709003687336D-25 /
      data BTH1CS(31) /  +.77297719892989706370926959871929D-26 /
      data BTH1CS(32) /  -.26967444512294640913211424080920D-26 /
      data BTH1CS(33) /  +.95749344518502698072295521933627D-27 /
      data BTH1CS(34) /  -.34569168448890113000175680827627D-27 /
      data BTH1CS(35) /  +.12681234817398436504211986238374D-27 /
      data BTH1CS(36) /  -.47232536630722639860464993713445D-28 /
      data BTH1CS(37) /  +.17850008478186376177858619796417D-28 /
      data BTH1CS(38) /  -.68404361004510395406215223566746D-29 /
      data BTH1CS(39) /  +.26566028671720419358293422672212D-29 /
      data BTH1CS(40) /  -.10450402527914452917714161484670D-29 /
      data BTH1CS(41) /  +.41618290825377144306861917197064D-30 /
      data BTH1CS(42) /  -.16771639203643714856501347882887D-30 /
      data BTH1CS(43) /  +.68361997776664389173535928028528D-31 /
      data BTH1CS(44) /  -.28172247861233641166739574622810D-31 /
C
      data PI4 / 0.785398163397448309615660845819876D0 /
      data NBM1, NBT12, NBM12, NBTH1, XMAX / 4*0, 0.D0 /
C     ------------------------------------------------------------------
      if (NBM1 .eq. 0) then
         ETA = 0.1D0 * D1MACH(3)
         call DINITS (BM1CS, 37, ETA,  NBM1 )
         call DINITS (BT12CS, 39, ETA, NBT12)
         call DINITS (BM12CS, 40, ETA, NBM12)
         call DINITS (BTH1CS, 44, ETA, NBTH1)
C
         XMAX = 0.04D0/D1MACH(4)
      end if
C
      if (X .lt. 4.D0) then
        AMPL = 0.D0
        THETA = 0.D0
        call DERM1 ('DBMP1',1,0,'X MUST BE .GE. 4','X',X,'.')
      end if
C
      if (X .le. 8.D0) then
         Z = (128.D0/(X*X) - 5.D0)/3.D0
         AMPL = (0.75D0 + DCSEVL (Z, BM1CS, NBM1))/SQRT(X)
         THETA = X - 3*PI4 + DCSEVL (Z, BT12CS, NBT12)/X
         RETURN
      end if
C
      if (X .gt. XMAX) then
        AMPL = 0.D0
        THETA = 0.D0
        call DERM1 ('DBMP1',2,0,'NO PRECISION BECAUSE X .GT. XMAX',
     *              'X',X,',')
        call DERV1 ('XMAX',XMAX,'.')
      end if
C
      Z = 128.D0/(X*X) - 1.D0
      AMPL = (0.75D0 + DCSEVL (Z, BM12CS, NBM12))/SQRT(X)
      THETA = X - 3*PI4 + DCSEVL (Z, BTH1CS, NBTH1)/X
      return
C
      end

.