! ************************************************************************* ! * * ! * Function dphinv * ! * * ! ************************************************************************* ! Double Precision Version 1.1 ! Written by Gordon A. Fenton, TUNS, 1990 ! Latest Update: May 18, 1997 ! ! PURPOSE to return the inverse of the standard normal cumulative ! distribution function ! ! For a given probability q lying between 0 and 1, this function returns ! the corresponding standard normal variate value, ie z such that ! ! q = P[ Z <= z ] ! ! The value z is obtained through inverse interpolation of the tabulated ! probability values (see U.S. Dept. of Commerce, Handbook of Mathematical ! Functions, pg 954). ! Notes: ! - the following is based on NPROBS = 768, so if this is changed, this ! routine will require some major modifications. ! - due to catastrophic cancellation, this routine is much more accurate ! for q values close to 0 than for q values close to 1. ! ! This module sets up a table of cumulative probabilities for ! the standard normal distribution, that is the following probabilities ! ! p = P[ Z <= z ] ! ! are stored in p(.) for z = -0.01, -0.02, ..., -7.68. ! ! REVISION HISTORY: ! 1.1 brought contents of dgblck.f directly into this routine; some compilers ! don't like externally defined common blocks. (May 18/97) !--------------------------------------------------------------------------- real*8 function dphinv(q) implicit real*8 (a-h,o-z) parameter (NPROBS = 768) dimension p(NPROBS) ! powers of two integer itwo(8) data itwo/128,64,32,16,8,4,2,1/ ! (empirical) fudge factors data f1, f2, f3/ & 0.245575499584d0, 0.297337076337d0, -0.003148948457d0/ ! some constants data pt01, sixth, half, one, two/ & 0.01d0, 0.1666666666666666667d0, 0.5d0, 1.d0, 2.d0/ data rt2pi/2.5066282746310005024d0/, big/30.d0/, zero/0.d0/ data a1/2.30753d0/, a2/0.27061d0/, b1/0.99229d0/, b2/0.04481d0/ ! probability tables data (p(i),i=1,15)/ & 4.960106436853684d-01,4.920216862830980d-01,4.880335265858873d-01, & 4.840465631471693d-01,4.800611941616275d-01,4.760778173458932d-01, & 4.720968298194789d-01,4.681186279860126d-01,4.641436074148279d-01, & 4.601721627229710d-01,4.562046874576832d-01,4.522415739794162d-01, & 4.482832133454389d-01,4.443299951940936d-01,4.403823076297575d-01/ data (p(i),i=16,30)/ & 4.364405371085672d-01,4.325050683249616d-01,4.285762840990993d-01, & 4.246545652652046d-01,4.207402905608970d-01,4.168338365175577d-01, & 4.129355773517854d-01,4.090458848579941d-01,4.051651283022042d-01, & 4.012936743170763d-01,3.974318867982395d-01,3.935801268019605d-01, & 3.897387524442028d-01,3.859081188011227d-01,3.820885778110474d-01/ data (p(i),i=31,45)/ & 3.782804781779807d-01,3.744841652766800d-01,3.706999810593465d-01, & 3.669282639639719d-01,3.631693488243809d-01,3.594235667820088d-01, & 3.556912451994533d-01,3.519727075758372d-01,3.482682734640177d-01, & 3.445782583896759d-01,3.409029737723226d-01,3.372427268482495d-01, & 3.335978205954577d-01,3.299685536605936d-01,3.263552202879201d-01/ data (p(i),i=46,60)/ & 3.227581102503477d-01,3.191775087825558d-01,3.156136965162226d-01, & 3.120669494173905d-01,3.085375387259869d-01,3.050257308975195d-01, & 3.015317875469662d-01,2.980559653948764d-01,2.945985162156980d-01, & 2.911596867883464d-01,2.877397188490270d-01,2.843388490463242d-01, & 2.809573088985644d-01,2.775953247534649d-01,2.742531177500736d-01/ data (p(i),i=61,75)/ & 2.709309037830058d-01,2.676288934689831d-01,2.643472921156776d-01, & 2.610862996928617d-01,2.578461108058648d-01,2.546269146713361d-01, & 2.514288950953101d-01,2.482522304535706d-01,2.450970936743094d-01, & 2.419636522230730d-01,2.388520680899867d-01,2.357624977792512d-01, & 2.326950923008975d-01,2.296499971647906d-01,2.266273523768682d-01/ data (p(i),i=76,90)/ & 2.236272924375995d-01,2.206499463426497d-01,2.176954375857332d-01, & 2.147638841636371d-01,2.118553985833966d-01,2.089700878716017d-01, & 2.061080535858131d-01,2.032693918280685d-01,2.004541932604497d-01, & 1.976625431226924d-01,1.948945212518084d-01,1.921502021036962d-01, & 1.894296547767121d-01,1.867329430371727d-01,1.840601253467595d-01/ data (p(i),i=91,105)/ & 1.814112548917973d-01,1.787863796143717d-01,1.761855422452580d-01, & 1.736087803386246d-01,1.710561263084818d-01,1.685276074668378d-01, & 1.660232460635296d-01,1.635430593276923d-01,1.610870595108309d-01, & 1.586552539314571d-01,1.562476450212546d-01,1.538642303727349d-01, & 1.515050027883437d-01,1.491699503309815d-01,1.468590563758959d-01/ data (p(i),i=106,120)/ & 1.445722996639096d-01,1.423096543559392d-01,1.400710900887691d-01, & 1.378565720320355d-01,1.356660609463827d-01,1.334995132427472d-01, & 1.313568810427307d-01,1.292381122400178d-01,1.271431505627982d-01, & 1.250719356371502d-01,1.230244030513434d-01,1.210004844210182d-01, & 1.190001074552007d-01,1.170231960231087d-01,1.150696702217083d-01/ data (p(i),i=121,135)/ & 1.131394464439773d-01,1.112324374478346d-01,1.093485524256919d-01, & 1.074876970745869d-01,1.056497736668553d-01,1.038346811213005d-01, & 1.020423150748192d-01,1.002725679544421d-01,9.852532904974787d-02, & 9.680048458561030d-02,9.509791779523907d-02,9.341750899347179d-02, & 9.175913565028082d-02,9.012267246445249d-02,8.850799143740207d-02/ data (p(i),i=136,150)/ & 8.691496194708503d-02,8.534345082196698d-02,8.379332241501425d-02, & 8.226443867766897d-02,8.075665923377101d-02,7.926984145339244d-02, & 7.780384052654643d-02,7.635850953673912d-02,7.493369953432705d-02, & 7.352925960964835d-02,7.214503696589381d-02,7.078087699168556d-02, & 6.943662333333178d-02,6.811211796672545d-02,6.680720126885809d-02/ data (p(i),i=151,165)/ & 6.552171208891650d-02,6.425548781893586d-02,6.300836446397839d-02, & 6.178017671181191d-02,6.057075800205902d-02,5.937994059479307d-02, & 5.820755563855307d-02,5.705343323775425d-02,5.591740251946942d-02, & 5.479929169955799d-02,5.369892814811972d-02,5.261613845425206d-02, & 5.155074849008934d-02,5.050258347410369d-02,4.947146803364810d-02/ data (p(i),i=166,180)/ & 4.845722626672283d-02,4.745968180294735d-02,4.647865786372007d-02, & 4.551397732154983d-02,4.456546275854306d-02,4.363293652403194d-02, & 4.271622079132897d-02,4.181513761359495d-02,4.092950897880737d-02, & 4.005915686381706d-02,3.920390328748263d-02,3.836357036287125d-02, & 3.753798034851680d-02,3.672695569872630d-02,3.593031911292582d-02/ data (p(i),i=181,195)/ & 3.514789358403880d-02,3.437950244589000d-02,3.362496941962834d-02, & 3.288411865916385d-02,3.215677479561369d-02,3.144276298075271d-02, & 3.074190892946599d-02,3.005403896119979d-02,2.937898004040940d-02, & 2.871655981600180d-02,2.806660665977251d-02,2.742894970383680d-02, & 2.680341887705495d-02,2.618984494045268d-02,2.558805952163862d-02/ data (p(i),i=196,210)/ & 2.499789514822043d-02,2.441918528022258d-02,2.385176434150854d-02, & 2.329546775021185d-02,2.275013194817921d-02,2.221559442943144d-02, & 2.169169376764679d-02,2.117826964267228d-02,2.067516286607007d-02, & 2.018221540570442d-02,1.969927040937691d-02,1.922617222751732d-02, & 1.876276643493774d-02,1.830889985165896d-02,1.786442056281656d-02/ data (p(i),i=211,225)/ & 1.742917793765708d-02,1.700302264763282d-02,1.658580668360504d-02, & 1.617738337216612d-02,1.577760739109052d-02,1.538633478392548d-02, & 1.500342297373219d-02,1.462873077598925d-02,1.426211841066888d-02, & 1.390344751349859d-02,1.355258114642000d-02,1.320938380725628d-02, & 1.287372143860205d-02,1.254546143594659d-02,1.222447265504473d-02/ data (p(i),i=226,240)/ & 1.191062541854704d-02,1.160379152190355d-02,1.130384423855280d-02, & 1.101065832441139d-02,1.072411002167578d-02,1.044407706195111d-02, & 1.017043866871969d-02,9.903075559164254d-03,9.641869945358317d-03, & 9.386705534838558d-03,9.137467530572652d-03,8.894042630336774d-03, & 8.656319025516557d-03,8.424186399345668d-03,8.197535924596155d-03/ data (p(i),i=241,255)/ & 7.976260260733725d-03,7.760253550553653d-03,7.549411416309215d-03, & 7.343630955348346d-03,7.142810735271399d-03,6.946850788624337d-03, & 6.755652607140672d-03,6.569119135546753d-03,6.387154764943170d-03, & 6.209665325776159d-03,6.036558080412646d-03,5.867741715332553d-03, & 5.703126332950670d-03,5.542623443082595d-03,5.386145954066668d-03/ data (p(i),i=256,270)/ & 5.233608163555781d-03,5.084925748991054d-03,4.940015757770644d-03, & 4.798796597126176d-03,4.661188023718732d-03,4.527111132967332d-03, & 4.396488348121286d-03,4.269243409089352d-03,4.145301361036025d-03, & 4.024588542758334d-03,3.907032574852809d-03,3.792562347685491d-03, & 3.681108009174983d-03,3.572600952399752d-03,3.466973803040674d-03/ data (p(i),i=271,285)/ & 3.364160406669203d-03,3.264095815891321d-03,3.166716277357817d-03, & 3.071959218650500d-03,2.979763235054556d-03,2.890068076226160d-03, & 2.802814632765049d-03,2.717944922701276d-03,2.635402077904969d-03, & 2.555130330427924d-03,2.477074998785855d-03,2.401182474189245d-03, & 2.327400206731556d-03,2.255676691542308d-03,2.185961454913232d-03/ data (p(i),i=286,300)/ & 2.118205040404608d-03,2.052358994939774d-03,1.988375854894309d-03, & 1.926209132187884d-03,1.865813300384045d-03,1.807143780806431d-03, & 1.750156928676083d-03,1.694810019277238d-03,1.641061234157026d-03, & 1.588869647364877d-03,1.538195211738036d-03,1.488998745237446d-03, & 1.441241917340019d-03,1.394887235492248d-03,1.349898031630096d-03/ data (p(i),i=301,315)/ & 1.306238448769467d-03,1.263873427672299d-03,1.222768693592260d-03, & 1.182890743104407d-03,1.144206831022698d-03,1.106684957409247d-03, & 1.070293854678923d-03,1.035002974802841d-03,1.000782476614011d-03, & 9.676032132183563d-04,9.354367195140999d-04,9.042551998223409d-04, & 8.740315156315670d-04,8.447391734586284d-04,8.163523128285638d-04/ data (p(i),i=316,330)/ & 7.888456943755737d-04,7.621946880672361d-04,7.363752615539311d-04, & 7.113639686453650d-04,6.871379379158485d-04,6.636748614399681d-04, & 6.409529836600560d-04,6.189510903868352d-04,5.976484979344154d-04, & 5.770250423907672d-04,5.570610690246212d-04,5.377374218296949d-04, & 5.190354332069722d-04,5.009369137857219d-04,4.834241423837776d-04/ data (p(i),i=331,345)/ & 4.664798561075492d-04,4.500872405921174d-04,4.342299203816563d-04, & 4.188919494503698d-04,4.040578018640217d-04,3.897123625820324d-04, & 3.758409184000837d-04,3.624291490330445d-04,3.494631183379715d-04, & 3.369292656768813d-04,3.248143974188780d-04,3.131056785812003d-04, & 3.017906246086373d-04,2.908570932907434d-04,2.802932768161773d-04/ data (p(i),i=346,360)/ & 2.700876939634748d-04,2.602291824274666d-04,2.507068912805378d-04, & 2.415102735678360d-04,2.326290790355250d-04,2.240533469910931d-04, & 2.157733992947175d-04,2.077798334806213d-04,2.000635160073205d-04, & 1.926155756356333d-04,1.854273969332782d-04,1.784906139048473d-04, & 1.717971037459309d-04,1.653389807201100d-04,1.591085901575340d-04/ data (p(i),i=361,375)/ & 1.530985025737555d-04,1.473015079074726d-04,1.417106098758194d-04, & 1.363190204458020d-04,1.311201544204847d-04,1.261076241384867d-04, & 1.212752342853580d-04,1.166169768153681d-04,1.121270259822471d-04, & 1.077997334773883d-04,1.036296236740311d-04,9.961138897591672d-05, & 9.573988526891469d-05,9.201012747410561d-05,8.841728520080404d-05/ data (p(i),i=376,390)/ & 8.495667849799789d-05,8.162377370268624d-05,7.841417938358505d-05, & 7.532364237868341d-05,7.234804392511995d-05,6.948339587986524d-05, & 6.672583702968470d-05,6.407162948887459d-05,6.151715518325534d-05, & 5.905891241892255d-05,5.669351253425669d-05,5.441767663369976d-05, & 5.222823240182019d-05,5.012211099618837d-05,4.809634401760274d-05/ data (p(i),i=391,405)/ & 4.614806055620888d-05,4.427448431207072d-05,4.247293078876124d-05, & 4.074080455855082d-05,3.907559659778755d-05,3.747488169107352d-05, & 3.593631590285383d-05,3.445763411505314d-05,3.303664762940245d-05, & 3.167124183311996d-05,3.035937392661827d-05,2.909907071193095d-05, & 2.788842644056393d-05,2.672560071949210d-05,2.560881647404153d-05/ data (p(i),i=406,420)/ & 2.453635796640967d-05,2.350656886859557d-05,2.251785038852544d-05, & 2.156865944818060d-05,2.065750691254679d-05,1.978295586822407d-05, & 1.894361995055329d-05,1.813816171813091d-05,1.736529107360408d-05, & 1.662376372965224d-05,1.591237971908220d-05,1.522998194797792d-05, & 1.457545479086707d-05,1.394772272688124d-05,1.334574901590634d-05/ data (p(i),i=421,435)/ & 1.276853441373497d-05,1.221511592525306d-05,1.168456559470741d-05, & 1.117598933212056d-05,1.068852577493443d-05,1.022134518398408d-05, & 9.773648372917573d-06,9.344665670196364d-06,8.933655912827005d-06, & 8.539905470991814d-06,8.162727302763068d-06,7.801460038101355d-06, & 7.455467091355145d-06,7.124135801495341d-06,6.806876599334045d-06/ data (p(i),i=436,450)/ & 6.503122200992801d-06,6.212326826901514d-06,5.933965445624679d-06, & 5.667533041826751d-06,5.412543907703856d-06,5.168530957224142d-06, & 4.935045062533278d-06,4.711654411897247d-06,4.497943888567909d-06, & 4.293514469971870d-06,4.097982646636362d-06,3.910979860280711d-06, & 3.732151960514484d-06,3.561158679597556d-06,3.397673124730062d-06/ data (p(i),i=451,465)/ & 3.241381287353394d-06,3.091981568956177d-06,2.949184322891521d-06, & 2.812711411724217d-06,2.682295779638856d-06,2.557681039451524d-06, & 2.438621073779427d-06,2.324879649934414d-06,2.216230048117548d-06, & 2.112454702502847d-06,2.013344854809340d-06,1.918700219970900d-06, & 1.828328663524165d-06,1.742045890344663d-06,1.659675144371463d-06/ data (p(i),i=466,480)/ & 1.581046918970512d-06,1.505998677596157d-06,1.434374584420136d-06, & 1.366025244606141d-06,1.300807453917281d-06,1.238583957352471d-06, & 1.179223216516399d-06,1.122599185436174d-06,1.068591094545936d-06, & 1.017083242568706d-06,9.679647960327358d-07,9.211295961671412d-07, & 8.764759729292055d-07,8.339065659229126d-07,7.933281519755972d-07/ data (p(i),i=481,495)/ & 7.546514791463692d-07,7.177911069469003d-07,6.826652525616647d-07, & 6.491956428613364d-07,6.173073720091949d-07,5.869287644666382d-07, & 5.579912432097829d-07,5.304292029750950d-07,5.041798883575367d-07, & 4.791832765903206d-07,4.553819648407320d-07,4.327210618617021d-07, & 4.111480838439311d-07,3.906128543183266d-07,3.710674079633336d-07/ data (p(i),i=496,510)/ & 3.524658981764252d-07,3.347645082736184d-07,3.179213661852820d-07, & 3.018964625208491d-07,2.866515718791945d-07,2.721501772855827d-07, & 2.583573976399725d-07,2.452399180653704d-07,2.327659230486003d-07, & 2.209050322695440d-07,2.096282390183694d-07,1.989078511037129d-07, & 1.887174341580604d-07,1.790317572498343d-07,1.698267407147598d-07/ data (p(i),i=511,525)/ & 1.610794061221380d-07,1.527678282945667d-07,1.448710893025085d-07, & 1.373692343578420d-07,1.302432295332016d-07,1.234749212365168d-07, & 1.170469973726320d-07,1.109429501263468d-07,1.051470403035407d-07, & 9.964426316933494d-08,9.442031572442990d-08,8.946156536290779d-08, & 8.475501985682844d-08,8.028829861495895d-08,7.604960516488729d-08/ data (p(i),i=526,540)/ & 7.202770080965977d-08,6.821187941186212d-08,6.459194325982505d-08, & 6.115817997230602d-08,5.790134039964602d-08,5.481261748095645d-08, & 5.188362601842435d-08,4.910638333128551d-08,4.647329075344129d-08, & 4.397711594005889d-08,4.161097594981976d-08,3.936832107075916d-08, & 3.724291935887129d-08,3.522884185984314d-08,3.332044848542857d-08/ data (p(i),i=541,555)/ & 3.151237451708229d-08,2.979951771053636d-08,2.817702597603999d-08, & 2.664028560996721d-08,2.518491005446115d-08,2.380672916270041d-08, & 2.250177894826861d-08,2.126629179795917d-08,2.009668712817647d-08, & 1.898956246588774d-08,1.794168493584716d-08,1.694998313655083d-08, & 1.601153938809098d-08,1.512358233576103d-08,1.428347989392278d-08/ data (p(i),i=556,570)/ & 1.348873251527842d-08,1.273696677129993d-08,1.202592923015495d-08, & 1.135348061903221d-08,1.071759025831089d-08,1.011633075554139d-08, & 9.547872947704290d-09,9.010481080699081d-09,8.502508215475082d-09, & 8.022391850663509d-09,7.568649751997725d-09,7.139875979218420d-09, & 6.734737101557546d-09,6.351968593271951d-09,5.990371401063532d-09/ data (p(i),i=571,585)/ & 5.648808675570940d-09,5.326202659455512d-09,5.021531724924521d-09, & 4.733827553845581d-09,4.462172453901613d-09,4.205696804522029d-09, & 3.963576626597628d-09,3.735031270249742d-09,3.519321215174624d-09, & 3.315745978326164d-09,3.123642123930022d-09,2.942381371044379d-09, & 2.771368794094646d-09,2.610041112012914d-09,2.457865061808032d-09/ data (p(i),i=586,600)/ & 2.314335852578571d-09,2.178975696160569d-09,2.051332410772609d-09, & 1.930978094185322d-09,1.817507863099436d-09,1.710538655567010d-09, & 1.609708093434255d-09,1.514673401922662d-09,1.425110383596567d-09, & 1.340712444091877d-09,1.261189667101099d-09,1.186267936225734d-09, & 1.115688101417170d-09,1.049205187833157d-09,9.865876450377014d-10/ data (p(i),i=601,615)/ & 9.276166345691163d-10,8.720853539929702d-10,8.197983956451329d-10, & 7.705711383542473d-10,7.242291705137658d-10,6.806077429504144d-10, & 6.395512501096636d-10,6.009127381488439d-10,5.645534385958076d-10, & 5.303423262948813d-10,4.981557004231254d-10,4.678767874181614d-10, & 4.393953647146701d-10,4.126074042396769d-10,3.874147346675664d-10/ data (p(i),i=616,630)/ & 3.637247214840693d-10,3.414499639547374d-10,3.205080081373417d-10, & 3.008210751196840d-10,2.823158037043269d-10,2.649230067999402d-10, & 2.485774408153008d-10,2.332175873867525d-10,2.187854468029041d-10, & 2.052263425218952d-10,1.924887362065496d-10,1.805240527313422d-10, & 1.692865146423052d-10,1.587329855769880d-10,1.488228221762322d-10/ data (p(i),i=631,645)/ & 1.395177340430679d-10,1.307816513264241d-10,1.225805995286331d-10, & 1.148825811560307d-10,1.076574638512162d-10,1.008768746639295d-10, & 9.451410013495054d-11,8.854399188407803d-11,8.294287740902139d-11, & 7.768847581709834d-11,7.275981822590379d-11,6.813717258273434d-11, & 6.380197266544754d-11,5.973675103973090d-11,5.592507575942690d-11/ data (p(i),i=646,660)/ & 5.235149060764004d-11,4.900145868690975d-11,4.586130917672491d-11, & 4.291818708617980d-11,4.016000583859125d-11,3.757540253348849d-11, & 3.515369573951724d-11,3.288484567954244d-11,3.075941667656468d-11, & 2.876854173604333d-11,2.690388914682007d-11,2.515763098911859d-11, & 2.352241344404167d-11,2.199132880464250d-11,2.055788909399524d-11/ data (p(i),i=661,675)/ & 1.921600120077499d-11,1.795994344767311d-11,1.678434351254311d-11, & 1.568415762649946d-11,1.465465097730285d-11,1.369137925025025d-11, & 1.279017124247953d-11,1.194711249009093d-11,1.115852985079933d-11, & 1.042097698796524d-11,9.731220704826869d-12,9.086228080565395d-12, & 8.483154362502072d-12,7.919331571248458d-12,7.392257778017862d-12/ data (p(i),i=676,690)/ & 6.899587015569721d-12,6.439119786395897d-12,6.008794133785011d-12, & 5.606677243315662d-12,5.230957544144618d-12,4.879937281169313d-12, & 4.552025530768030d-12,4.245731634354438d-12,3.959659025435898d-12, & 3.692499427235594d-12,3.443027399236997d-12,3.210095212234586d-12, & 2.992628032635060d-12,2.789619397847642d-12,2.600126965638169d-12/ data (p(i),i=691,705)/ & 2.423268521298989d-12,2.258218227411716d-12,2.104203101851851d-12, & 1.960499710509255d-12,1.826431061976962d-12,1.701363692195681d-12, & 1.584704927736075d-12,1.475900317055532d-12,1.374431219685123d-12, & 1.279812543885835d-12,1.191590623864497d-12,1.109341228159143d-12, & 1.032667691294276d-12,9.611991612689379d-13,8.945889558769916d-13/ data (p(i),i=706,720)/ & 8.325130212702670d-13,7.746684865636528d-13,7.207723086467529d-13, & 6.705600017118691d-13,6.237844463331584d-13,5.802147732383270d-13, & 5.396353172029230d-13,5.018446367696487d-13,4.666545957513176d-13, & 4.338895027178080d-13,4.033853048947575d-13,3.749888331162318d-13, & 3.485570946752429d-13,3.239566111061341d-13,3.010627981117446d-13/ data (p(i),i=721,735)/ & 2.797593850166439d-13,2.599378712863654d-13,2.414970178016698d-13, & 2.243423707173588d-13,2.083858158672078d-13,1.935451618009657d-13, & 1.797437496562159d-13,1.669100881779288d-13,1.549775123019207d-13, & 1.438838638157592d-13,1.335711927020462d-13,1.239854778550311d-13, & 1.150763659422917d-13,1.067969272592302d-13,9.910342749547497d-14/ data (p(i),i=736,750)/ & 9.195511439940067d-14,8.531401838998034d-14,7.914476622443238d-14, & 7.341440688571708d-14,6.809224890620008d-14,6.314970839286098d-14, & 5.856016706548484d-14,5.429883966255224d-14,5.034265011012925d-14, & 4.667011588719071d-14,4.326124005658145d-14,4.009741046439491d-14, & 3.716130564205506d-14,3.443680697493761d-14,3.190891672910920d-14/ data (p(i),i=751,765)/ & 2.956368155375915d-14,2.738812110130123d-14,2.537016142999860d-14, & 2.349857287541125d-14,2.176291209708596d-14,2.015346802575318d-14, & 1.866121145397440d-14,1.727774802974088d-14,1.599527442805322d-14, & 1.480653749004804d-14,1.370479613286938d-14,1.268378584624278d-14, & 1.173768560367169d-14,1.086108702736900d-14,1.004896565652634d-14/ data (p(i),i=766,768)/ & 9.296654178339954d-15,8.599817490408714d-15,7.954429471721532d-15/ ! conversion statement functions float(i) = dble(i) exp(x) = dexp(x) sqrt(x) = dsqrt(x) alog(x) = dlog(x) !------------------------------ start executable statements ------------------- if( q .lt. half ) then qr = q zsign = -one else qr = one - q zsign = one endif ! find closest p(i) to qr if( qr .gt. p(1) ) then r = half - qr r2 = qr - p(1) i = 0 else if( qr .lt. p(768) ) then if( qr .eq. zero ) then dphinv = zsign*big else t = sqrt( alog(one/(qr*qr)) ) zo = t - (a1 + a2*t)/(one + b1*t + b2*t*t) dphinv = zsign*zo endif return else ! find i such that p(i) < qr < p(i+1) if( qr .lt. p(512) ) then if( qr .gt. p(513) ) then i = 512 go to 20 endif i = 640 js = 2 else i = 256 js = 1 endif do 10 j = js, 8 if( qr .gt. p(i) ) then i = i - itwo(j) else i = i + itwo(j) endif 10 continue if( qr .gt. p(i) ) i = i - 1 20 r = p(i) - qr r2 = qr - p(i+1) endif ! now choose p(i) or p(i+1) if( r2 .lt. r ) then r = -r2 i = i + 1 endif ! inverse interpolation zo = pt01*float(i) zo2 = zo*zo zo3 = zo*zo2 t = r*rt2pi*exp(half*zo2) t2 = t*t t3 = t*t2 t4 = t2*t2 zo = zsign*(zo + t + half*zo*t2 + (two*zo2 + one)*t3*sixth & + t4*(f1*zo3+f2*zo+f3)) ! iterate once r = q - dphi(zo) zo2 = zo*zo zo3 = zo*zo2 t = r*rt2pi*exp(half*zo2) t2 = t*t t3 = t*t2 t4 = t2*t2 dphinv = zo + t + half*zo*t2 + (two*zo2 + one)*t3*sixth & + t4*(f1*zo3+f2*zo+f3) return end