Example 2 -- interpolated periodic function and derivative


Please note that to compile the graphics portion of the code below, you must link with a double precision version of NCAR Graphics, or else you need to convert the code back to single precision before you call the graphics routines.

      PROGRAM FTEX02
C
C  Example of CURVP1DP, CURVPIDP.
C
C  Define dimensions, declare arrays.
C
      PARAMETER (IDIM=10,IOUT=201)
      DOUBLE PRECISION X
      DOUBLE PRECISION Y
      DOUBLE PRECISION YP
      DOUBLE PRECISION TEMP
      DOUBLE PRECISION XO
      DOUBLE PRECISION YO
      DOUBLE PRECISION YI
      DOUBLE PRECISION PERIOD
      DOUBLE PRECISION SIGMA
      DOUBLE PRECISION XR
      DOUBLE PRECISION XL
      DOUBLE PRECISION XINC
      DOUBLE PRECISION CURVP2DP
      DOUBLE PRECISION CURVPIDP
      DIMENSION X(IDIM),Y(IDIM),YP(IDIM),TEMP(IDIM,2)
      DIMENSION XO(IOUT),YO(IOUT),YI(IOUT)
C
C Specify the input data.
C
      DATA X/0.000D0,0.210D0,0.360D0,0.540D0,1.000D0,1.500D0,1.970D0,
     +     2.300D0,2.500D0,2.700D0/
      DATA Y/0.000D0,2.600D0,3.000D0,2.500D0,0.000D0,-1.000D0,0.000D0,
     +     0.800D0,0.920D0,0.700D0/
C
C  Call CURVP1DP setup.
C
      PERIOD = 3.D0
      SIGMA = 1.D0
      CALL CURVP1DP(IDIM,X,Y,PERIOD,YP,TEMP,SIGMA,IERR)
C
C  Call CURVP2DP and calculate the interpolated values and the integrals.
C
      XR = 5.D0
      XL = -1.D0
      XINC = (XR-XL)/ (IOUT-1)
      DO 10 I = 1,IOUT
          XO(I) = XL + (I-1)*XINC
          YO(I) = CURVP2DP(XO(I),IDIM,X,Y,PERIOD,YP,SIGMA)
          YI(I) = CURVPIDP(0.D0,XO(I),IDIM,X,Y,PERIOD,YP,SIGMA)
   10 CONTINUE
C
C  Draw a plot of the interpolated functions and mark the original
C  points.
C
      CALL DRWFT2(XL,XR,IDIM,X,Y,IOUT,XO,YO,YI)
C
      STOP
      END
      SUBROUTINE DRWFT2(XL,XR,II,X,Y,IO,XO,YO,YI)
C
C  Define error file, Fortran unit number, and workstation type,
C  and workstation ID.
C
      PARAMETER (IERRF=6,LUNIT=2,IWTYPE=1,IWKID=1)
      DOUBLE PRECISION XL
      DOUBLE PRECISION XR
      DOUBLE PRECISION X
      DOUBLE PRECISION Y
      DOUBLE PRECISION XO
      DOUBLE PRECISION YO
      DOUBLE PRECISION YI
      DOUBLE PRECISION YPOS_TOP
      DOUBLE PRECISION YB
      DOUBLE PRECISION YT
C
      DATA YPOS_TOP/0.85D0/
C
C  Open GKS, open and activate a workstation.
C
      CALL GOPKS(IERRF,ISZDM)
      CALL GOPWK(IWKID,LUNIT,IWTYPE)
      CALL GACWK(IWKID)
C
C  Define a color table.
C
      CALL GSCR(IWKID,0,1.0D0,1.0D0,1.0D0)
      CALL GSCR(IWKID,1,0.0D0,0.0D0,0.0D0)
      CALL GSCR(IWKID,2,1.0D0,0.0D0,0.0D0)
      CALL GSCR(IWKID,3,0.0D0,1.0D0,0.0D0)
      CALL GSCR(IWKID,4,0.0D0,0.0D0,1.0D0)
      CALL GSCLIP(0)
C
C  Plot main title.
C
      CALL PLCHHQ(0.50D0,0.95D0,':F25:Demo for CURVP1DP, CURVPIDP',
     +            0.032D0,0.D0,0.D0)
C
C  Graph the interpolated function values and mark the original
C  input data points.
C
      YB = -2.0D0
      YT = 3.0D0
      CALL BKGFT2(XL,XR,YPOS_TOP,'Function',0.42D0,YB,YT)
      CALL GRIDAL(6,5,5,1,1,1,10,XL,YB)
C
C  Mark the original data points.
C
      CALL GSMKSC(2.D0)
      CALL GSPMCI(4)
      CALL GPM(II,X,Y)
C
C  Graph the interpolated function values.
C
      CALL GPL(IO,XO,YO)
C
C  Graph the integral.
C
      YB = -1.0D0
      YT = 4.0D0
      CALL BKGFT2(XL,XR,YPOS_TOP-0.47D0,'Integral (from X = 0.)',0.2D0,
     +            YB,YT)
      CALL GRIDAL(6,5,5,1,1,1,10,XL,YB)
      CALL GPL(IO,XO,YI)
      CALL GSPLCI(1)
C
C  Indicate the period.
C
      CALL DRWPRD(0.D0,3.D0,6.5D0)
C
      CALL FRAME
C
      CALL GDAWK(IWKID)
      CALL GCLWK(IWKID)
      CALL GCLKS
C
      RETURN
      END
      SUBROUTINE BKGFT2(XL,XR,YPOS,LABEL,XLP,YB,YT)
      DOUBLE PRECISION XL
      DOUBLE PRECISION XR
      DOUBLE PRECISION YPOS
      DOUBLE PRECISION XLP
      DOUBLE PRECISION YB
      DOUBLE PRECISION YT
      DOUBLE PRECISION XX
      DOUBLE PRECISION YY
      DIMENSION XX(2),YY(2)
      CHARACTER*(*) LABEL
C
      CALL SET(0.D0,1.D0,0.D0,1.D0,0.D0,1.D0,0.D0,1.D0,1)
      CALL PCSETI('FN',21)
      CALL PLCHHQ(XLP,YPOS-0.03D0,LABEL,0.025D0,0.D0,-1.D0)
      CALL SET(0.13D0,0.93D0,YPOS-0.25D0,YPOS,XL,XR,YB,YT,1)
      XX(1) = XL
      XX(2) = XR
      YY(1) = 0.D0
      YY(2) = 0.D0
      CALL GSPLCI(2)
      CALL GPL(2,XX,YY)
      CALL GSPLCI(1)

      CALL GASETI('LTY',1)
      CALL PCSETI('FN',21)
      CALL GASETR('XLS',0.02D0)
      CALL GASETC('XLF','(F4.1)')
      CALL GASETR('YLS',0.02D0)
      CALL GASETC('YLF','(F5.2)')
      CALL GASETR('XMJ',0.02D0)
      CALL GASETR('YMJ',0.02D0)
C
      RETURN
      END
      SUBROUTINE DRWPRD(XL,XR,Y)
      DOUBLE PRECISION XL
      DOUBLE PRECISION XR
      DOUBLE PRECISION Y
      DOUBLE PRECISION XX
      DOUBLE PRECISION YY
      DOUBLE PRECISION YOFF
      DOUBLE PRECISION XMID
      DOUBLE PRECISION XB
      DOUBLE PRECISION XE
      DOUBLE PRECISION CFUX
      DOUBLE PRECISION YI
C
C  Draws a bounding indicator for the period of the function.
C
      DIMENSION XX(2),YY(2)
C
      YOFF = 0.4D0
      CALL GSPLCI(2)
      XMID = 0.5D0* (XR-XL)
      CALL PLCHHQ(XMID,Y,':F25:Period',0.02D0,0.D0,0.D0)
C
C  Vertical lines at the period limits.
C
      XX(1) = XL
      XX(2) = XL
      YY(1) = Y + YOFF
      YY(2) = Y - YOFF
      CALL GPL(2,XX,YY)
      XX(1) = XR
      XX(2) = XR
      CALL GPL(2,XX,YY)
C
C  Horizontal lines between label and vertical lines.
C
      CALL PCSETI('TE',1)
      CALL PCGETR('XB',XB)
      CALL PCGETR('XE',XE)
      XX(1) = XL
      XX(2) = CFUX(XB) - 0.09D0
      YY(1) = Y
      YY(2) = Y
      CALL GPL(2,XX,YY)
C
C  Left arrow.
C
      XX(1) = XL
      YY(1) = Y
      YI = 0.5D0*YOFF
      YY(2) = Y + YI
      XX(2) = XL + YI
      CALL GPL(2,XX,YY)
      YY(2) = Y - YI
      CALL GPL(2,XX,YY)
C
      XX(1) = XR
      XX(2) = CFUX(XE) + 0.09D0
      YY(1) = Y
      YY(2) = Y
      CALL GPL(2,XX,YY)
C
C  Right arrow.
C
      XX(1) = XR
      YY(1) = Y
      YY(2) = Y + YI
      XX(2) = XR - YI
      CALL GPL(2,XX,YY)
      YY(2) = Y - YI
      CALL GPL(2,XX,YY)
C
      RETURN
      END

home | contents | defs | params | procedures | exmpls