English
Source / code snippets

AKIMA-Interpolation with Kubischen Splines

 

p.specht

One of US-Volkswirt Hiroshi Akima 1978 vorgestelltes take action uses to the usual "linearen Interpolation" between Stützstellen Polynome third Grades.
In its Berechnung from whom too unregelmäßig distributed - Stützstellen (Messwertetabelle) weg too The each benachbarten Stützstellen one as well as additional one gewogener rose one.
Window Title mkstr$(" ",26)+upper$("AKIMA-Interpolation through Kubischer Splines [H. Akima, 1978]")
'  Portiert to XProfan11 (CL) CopyLeft 2014-08 by P. woodpecker, Vienna (Austria).
'{ No however geartete warranty. Use on alleinige menace the Anwenders!
'
' Übersetzungsbasis: Turbo Pascal manuals 1974; for Detailklärungen watts uses:
' F90/95: https://www.tacc.utexas.edu/documents/13601/162125/fortran_class.pdf
' F95: https://www.mrao.cam.ac.uk/~rachael/compphys/SelfStudyF95.pdf
'
' Pascal-DemoSource: https://jean-pierre.moreau.pagesperso-orange.fr/Pascal/akima_pas.txt
' Orig. 2D-Algorithmus: Hiroshi Akima,u.s. Department of Commerce, NTIA/ITS, F90-Version of 1995/08,
' published in "Hiroshi Akima, 'A Method of Bivariate Interpolation and Smooth Surface Fitting
' for Irregularly Distributed Data Points", ACM Transactions on Math. software Vol.4 No.2, June 1978
' pp. 148-159. copyright 1978, Association for Computing Machinery, Inc.
' Fortran-90: TOMS 22.yr,No3 (Sept. 1996) p.357ff,  https://calgo.acm.org/760.gz (Abfrg. 2014-08)
' Test Data: https://cran.r-project.org/web/packages/akima/akima.pdf
'
'********************************************************
'*   Program to demonstrate the Akima spline fitting    *
'*       of Function SIN(X) in double precision         *
'* ---------------------------------------------------- *
'*   Reference: BASIC Scientific Subroutines, Vol. II   *
'*   By F.R. Ruckdeschel, BYTE/McGRAWW-HILL, 1981 [1].  *
'*                                                      *
'*             Pascal Version by J.-P. Moreau, Paris.   *
'*                       (www.jpmoreau.fr)              *
'* ---------------------------------------------------- *
'* SAMPLE RUN:                                          *
'*                                                      *
'* Akima spline fitting of SIN(X):                      *
'*                                                      *
'*   X   SIN(X) HANDBOOK   AKIMA INTERPOLATION   ERROR  *
'* ---------------------------------------------------- *
'* 0.00    0.0000000          0.0000000       0.0000000 *
'* 0.05    0.0499792          0.0500402      -0.0000610 *
'* 0.10    0.0998334          0.0998435      -0.0000101 *
'* 0.15    0.1494381          0.1494310       0.0000072 *
'* 0.20    0.1986693          0.1986459       0.0000235 *
'* 0.25    0.2474040          0.2474157      -0.0000118 *
'* 0.30    0.2955202          0.2955218      -0.0000016 *
'* 0.35    0.3428978          0.3428916       0.0000062 *
'* 0.40    0.3894183          0.3894265      -0.0000081 *
'* 0.45    0.4349655          0.4349655       0.0000000 *
'* 0.50    0.4794255          0.4794204       0.0000051 *
'* 0.55    0.5226872          0.5226894      -0.0000021 *
'* 0.60    0.5646425          0.5646493      -0.0000068 *
'* 0.65    0.6051864          0.6051821       0.0000043 *
'* 0.70    0.6442177          0.6442141       0.0000035 *
'* 0.75    0.6816388          0.6816405      -0.0000017 *
'* ---------------------------------------------------- *
'*                                                      *
'*******************************************************}
'}

PROC Interpol_Akima

    '********************************************************
    '*       Akima Spline Fitting Subroutine                *
    '* -----------------------------------------------------*
    '* The input table is X(i), Y(i), where Y(i) is the     *
    '* dependant variable. The interpolation point is xx,   *
    '* which is assumed to be in the interval of the table  *
    '* with at least one table value to the left, and three *
    '* to the right. The interpolated returned value is yy. *
    '* n is returned as on error checked (n=0 implies error). *
    '* It is means assumed that the X(i) are in ascending    *
    '* order!                                               *
    '********************************************************
    parameters i&
    n&=1
    ' {special case xx=0}

    if xx!=0: yy!=0:return:endif

        ' {Check to see if interpolation point is correct}

        if (xx!<X![1]) or (xx!>=X![iv&-3]): n&=0:return:endif

            X![0]=2*X![1]-X![2]
            ' {Calculate Akima coefficients, a and b}

            whileloop iv&-1:i&=&Loop

                '   {Invoice values i to i+2}
                XM![i&+2]=(Y![i&+1]-Y![i&])/(X![i&+1]-X![i&])

            endwhile

            XM![iv&+2]=2*XM![iv&+1]-XM![iv&]
            XM![iv&+3]=2*XM![iv&+2]-XM![iv&+1]
            XM![2]=2*XM![3]-XM![4]
            XM![1]=2*XM![2]-XM![3]

            whileloop iv&:i&=&Loop

                a!=ABS(XM![i&+3]-XM![i&+2])
                b!=ABS(XM![i&+1]-XM![i&])

                if (a!+b!)<>0

                    Z![i&]=(a!*XM![i&+1]+b!*XM![i&+2])/(a!+b!)

                else

                    Z![i&]=(XM![i&+2]+XM![i&+1])/2

                endif

            endwhile

            ' {Find relevant table interval}
            i&=0

            repeat

                inc i&

            until xx!<=X![i&]

            dec i&
            ' {Begin interpolation}
            b!=X![i&+1]-X![i&]
            a!=xx!-X![i&]
            yy!=Y![i&]+Z![i&]*a!+(3*XM![i&+2]-2*Z![i&]-Z![i&+1])*a!*a!/b!
            yy!=yy!+(Z![i&]+Z![i&+1]-2*XM![i&+2])*a!*a!*a!/(b!*b!)

        ENDPROC

        '{ Main Program AKIMA-SINTEST
        Windowstyle 24:font 0:set("decimals",8)
        CLS'Window 0,0-%maxx,%maxy
        var fmt$=" 0.0000000;-0.0000000"
        var SIZE& = 14
        declare X![SIZE&+1],Y![SIZE&+1],XM![SIZE&+4],Z![SIZE&+1]
        declare a!,b!,xx!,yy!, i&,iv&,n&
        iv&=14'{Number of points in table}
        CLS
        Font 2:Print "\n";tab(10);"  Hiroshi Akima's Interpolation with kubischen Splines"
        Font 0:print tab(10);" "+mkstr$("-",53)
        print tab(10);" Base are 14 unregelmäßig distributed Stützwerte of/ one "
        print tab(10);" Eigenbau-Sinusfunktion. Errechnet becomes the absolute   "
        Print tab(10);" Error (as difference to internen Sin(x)-function).  "
        print tab(10);" in the Programmtext are The Vergleichswerte to find!   \n"
        ' -----------------------------------------------------------------
        ' Input y=Sine(x) table
        ' -----------------------------------------------------------------
        ' Sine table values from  Handbook of mathematical functions
        ' by M. Abramowitz and I.A. Stegun, NBS, june 1964
        ' -----------------------------------------------------------------
        X![1]=0.000:Y![1]=0.00000000  :  X![2]=0.125:Y![2]=0.12467473
        X![3]=0.217:Y![3]=0.21530095  :  X![4]=0.299:Y![4]=0.29456472
        X![5]=0.376:Y![5]=0.36720285  :  X![6]=0.450:Y![6]=0.43496553
        X![7]=0.520:Y![7]=0.49688014  :  X![8]=0.589:Y![8]=0.55552980
        X![9]=0.656:Y![9]=0.60995199  :  X![10]=0.721:Y![10]=0.66013615
        X![11]=0.7853981634           : Y![11]=0.7071067812
        X![12]=0.849:Y![12]=0.75062005:  X![13]=0.911:Y![13]=0.79011709
        X![14]=0.972:Y![14]=0.82601466
        ' -----------------------------------------------------------------
        print "\n";mkstr$("-",77)
        font 2
        print tab(3+3);"X";tab(%pos+15);"SIN(X)";\
        tab(%pos+3);"AKIMA-Stützwertinterpolation";tab(%pos+2);"ABS.FEHLER"
        font 0
        print mkstr$("=",77)
        ' {main loop}
        xx!=0

        whileloop 16:i&=&Loop

            Interpol_Akima(i&)
            print tab(4);xx!;tab(%pos+6);SIN(xx!);tab(%pos+9);yy!;tab(%pos+11);stature$(fmt$,SIN(xx!) - yy!)
            casenot i& mod 4:print
            xx! = xx! + 0.05

        endwhile

        ' {print footer}
        print mkstr$("-",77)
        beep
        waitinput
        END
        '} of file akima.prf
 
XProfan 11
Computer: Gerät, daß es in Mikrosekunden erlaubt, 50.000 Fehler zu machen, zB 'daß' statt 'das'...
05/11/21  
 



Zum Quelltext


Topictitle, max. 100 characters.
 

Systemprofile:

no Systemprofil laid out. [anlegen]

XProfan:

 Posting  Font  Smilies  ▼ 

Please register circa a Posting To verfassen.
 

Topic-Options

838 Views

Themeninformationen

this Topic has 1 subscriber:

p.specht (1x)


Admins  |  AGB  |  Applications  |  Authors  |  Chat  |  Privacy Policy  |  Download  |  Entrance  |  Help  |  Merchantportal  |  Imprint  |  Mart  |  Interfaces  |  SDK  |  Services  |  Games  |  Search  |  Support

One proposition all XProfan, The there's!


My XProfan
Private Messages
Own Storage Forum
Topics-Remember-List
Own Posts
Own Topics
Clipboard
Log off
 Deutsch English Français Español Italia
Translations

Privacy Policy


we use Cookies only as Session-Cookies because of the technical necessity and with us there no Cookies of Drittanbietern.

If you here on our Website click or navigate, stimmst You ours registration of Information in our Cookies on XProfan.Net To.

further Information To our Cookies and moreover, How You The control above keep, find You in ours nachfolgenden Datenschutzerklärung.


all rightDatenschutzerklärung
i want none Cookie