!  First order frequency response

!  This is paired with "Second order frequency response."
!  It illustrates the equation governing e.g. the temperature response, within 
!  an insulated container, to a sinusoidally varying ambient temperature.

!  The omega slider at the right controls the circular frequency of the ambient temperature. 
!  The k slider sets the coupling constant k, which is a measure of the degree of insulation.

!  The large time series plane shows the input signal in cyan and the periodic system response in yellow. 
!  The t slider sets the time t, and the corresponding values of the input signal and system response are marked with cyan and yellow diamonds.
!  The left graph plane further illustrates the relationship between the input signal and the system response at the selected time.
!  Time can be animated using the [>>>] button, or stepped incrementally using the [<] and [>] buttons.

!  The minimal period P, and the time lag t_0, are displayed below the k slider.
!  The time lag is also indicated by a red bar on the main graphing window,
!  connecting the  t = 0  line to the first maximum of the system response.
!  
!  A yellow horizontal line indicating the amplitude appears when the cursor is in the main graphing window.
!  
!  The "Bode and Nyquist Plots" button opens up three other graph planes on the right side of the screen.

!  The top plane shows the amplitude response curve, and the amplitude, A, marked by a yellow line segment and diamond,
!  that corresponds to the selected circular frequency, omega.

!  The second plane shows the phase shift curve, with the phase shift, -phi, marked by a green line segment and diamond,
!  that corresponds to the selected circular frequency omega.

!  Rolling the cursor over either one of those windows causes a horizontal rule and a read out value to appear.

!  The lower plane shows the trajectory of the complex number  k/p(i omega)  as omega varies from 0 to +infinity,  
!  with the value corresponding to the selected circular frequency marked by a yellow diamond.

!  The line segment, or modulus, joining the origin to that complex number is the amplitude A and is marked in yellow.
!  The angle, or argument, of that complex number is the phase gain -phi and is marked in green.
!  The characteristic polynomial of the operator,  p(s) , is displayed below the complex plane. 

!  Note: These are not quite true Bode and Nyquist plots. A Bode plot  graphs  log(A)  vs  log(omega)  or  - phi  vs  log(omega).  
!  A Nyquist plot displays  k/p(i omega)  as  omega  ranges from  -infinity  to  +infinity, 
!  so it has a portion above the real axis which is symmetric with what is drawn here.


! notes on "gain" and frequency response and amplitude?


!! File: FrequencyResponseOrderOne
!! May 26, 2003 Hubert Hohn

PUBLIC PCFlag,Mac5Flag,M68KFlag,UnixFlag,xmax,ymax
PUBLIC toolLft,toolRgt,toolBas,toolTop,toolhdr,toolHgt,toolWid  ! tool boundaries
PUBLIC winLft,winRgt,winBas,winTop,winHgt,winWid      ! window
PUBLIC workLft,workRgt,workBas,workTop,workMidx  ! work area
PUBLIC backclr,black,drkgry,drkmid,midgry,litmid,litgry,white
PUBLIC red,yellow,green,cyan,blue,magenta,pink,colorscheme
PUBLIC planeclr,gridclr,rimclr,axisclr,axislabelclr,titleclr,rightsclr
PUBLIC numberlineclr,slotdrkclr,slotlgtclr,slideclr
PUBLIC largefonts, title$, SLUmode

DECLARE DEF QuitWithin, InfoWithin

LET toolHgt= 560
LET toolWid= 780
LET window$= "The d'Arbeloff Interactive Math Project"
LET colorscheme= 0

LET title$ = "First Order Frequency Response"

SUB ThisProgram
    CALL Tide
    CLEAR
END SUB


!! ---------------------------------------------------------
!! ------ Start TB4 Mac Header and Subs ------

!LET M68KFlag = 1
!LIBRARY "MacTools*", "HHLib.trc"
!ASK PIXELS winWid,winHgt
!LET winLft= 0
!LET winTop= 0
!LET winRgt= winWid-1
!LET winBas= winHgt-1
!SET WINDOW 0,winRgt,winBas,0
!CALL Palette
!CLEAR
!
!CALL ToolPanel
!CALL ThisProgram
!
!END
!EXTERNAL
!
!MODULE Mac4Parts
!    SUB SetTextFont(font,size,style$)
!        CALL MacTextFont(font)
!        CALL MacTextSize(size)
!        CALL MacTextFace(style$)
!    END SUB
!
!    SUB StringWidth(sw$,sl)
!        DECLARE DEF MacStringWidth
!        LET sl= MacStringWidth(sw$)
!    END SUB
!
!    SUB SetLineWeight(wgt)
!        CALL MacPenSize(wgt,wgt)
!    END SUB
!
!    SUB BoxDisk(Lft,Rgt,Bas,Top)
!        CALL MacPaintOval(Lft,Rgt,Bas,Top)
!    END SUB
!END MODULE

!! --- End TB4 Mac Header and Subs ---
!!---------------------------------------------------------
!!--- Start TB5 Cross-Platform header and subs ---

LIBRARY "c:\TB Gold 51a\TBLibs\TrueCtrl.trc"     ! windows
LIBRARY "c:\TB Gold 51a\TBLibs\HHLib.trc"

!LIBRARY ":TBLibs:TrueCtrl.trc"                   ! macintosh
!LIBRARY "HHLib.trc"

PUBLIC WinID
DECLARE PUBLIC OBJM_SET,OBJM_SYSINFO
LET winHgt= toolHgt
LET winWid= toolWid
DIM values(1)

CALL TC_Init

CALL Object(OBJM_SYSINFO,WinID,"MACHINE",system$,values())
IF system$="MAC" then
   LET Mac5Flag= 1
ELSE IF system$="WIN32" then
   LET PCFlag  = 1
END IF

CALL TC_SetUnitsToPixels                         ! 5.1 and up needs this
CALL TC_GetScreenSize(scrnLft,scrnRgt,scrnBas,scrnTop)
LET winLft= int((scrnRgt-scrnLft-winWid)/2)
LET winRgt= winLft+winWid-1
LET winTop= int((scrnBas-scrnTop-winHgt)/2) + 10
LET winBas= winTop+winHgt-1

CALL TC_Win_Create (WinID,"TITLE",winLft,winRgt,winBas,winTop)
LET values(1)= 2
CALL Object(OBJM_SET, WinID, "TYPE", "", values())

IF PCFlag=1 then                                 ! kill dithering
   LET  values(1)= 1
   CALL Object(OBJM_SET, WinID, "SOLID MIX", "", values())
END IF

LET  values(1)= 0
CALL TC_SetRect(WinID,winLft,winRgt,winBas,winTop)
CALL TC_Win_SetTitle(WinID,window$)
CALL TC_Show(WinID)

SET MODE "COLORSTANDARD"

ASK PIXELS winWid,winHgt                         ! must follow set mode
LET winLft= 0
LET winTop= 0
LET winRgt= winWid-1
LET winBas= winHgt-1
SET WINDOW 0,winRgt,winBas,0

CALL Palette

IF PCFlag=1 then
   LET  values(1)= 0                             ! now force solid colors
   CALL Object(OBJM_SET, WinID, "SOLID MIX", "", values())
   CALL TC_Win_RealizePalette(WinID)             ! some PCs need this

   CALL TC_Win_SetFont(WinID,"arial",9,"plain")
   CALL StringWidth("0",sw)
   IF sw>7 then LET largefonts=1 else LET largefonts=0
END IF

CALL TC_Win_Switch(WinID)

CALL ToolPanel
CALL ThisProgram

CALL SetTextFont(1,12,"bold")                    ! now shut down and clean up
LET quit$= "click the mouse or press a key to close..."
CALL PlotTextCJ(workmidx,(workbas+worktop)/2,quit$,yellow)
CALL TC_CleanUp

END
EXTERNAL

MODULE TB5Parts
    SUB StringWidth(sw$,sl)
        DECLARE PUBLIC WinID
        LET sl= StrWidth(WinID,sw$)
    END SUB

    SUB SetLineWeight(wgt)
        DECLARE PUBLIC OBJM_SET
        DECLARE PUBLIC WinID
        DIM values(1)
        LET values(1)= wgt
        CALL Object(OBJM_SET,WinID, "WIDTH", "", values())
    END SUB

    SUB SetTextFont(font,size,style$)
        DECLARE PUBLIC WinID,Mac5Flag,PCFlag,largefonts
        IF Mac5Flag=1 then
           SELECT CASE font
           CASE 4
                LET font$= "Courier"
           CASE 16
                LET font$= "Times"
           CASE else
                LET font$= "Geneva"
           END SELECT
        ELSE IF PCFlag=1 then
           IF largefonts=1 then
              IF size<12 then
                 LET size=  6
              ELSE IF size=14 then
                 LET size= 10
              ELSE IF size=18 then
                 LET size= 12
              ELSE IF size=24 then
                 LET size= 14
              ELSE IF size=12 then
                 LET size=  8
              ELSE
                 LET size= round(72/96 * size * .8)
              END IF
           ELSE
              IF size<12 then
                 LET size=  7
              ELSE IF size=14 then
                 LET size= 12
              ELSE IF size=12 then
                 LET size=  9
              ELSE IF size=18 then
                 LET size= 14
              ELSE IF size=24 then
                 LET size= 18
              ELSE
                 LET size= round(72/96 * size)
              END IF
           END IF
           SELECT CASE font
           CASE 4
                LET font$= "Courier New"
           CASE 16
                LET font$= "Times New Roman"
           CASE else
                LET font$= "Verdana"
           END SELECT
        END IF
        IF style$= "normal" then LET style$= "plain"
        CALL TC_Win_SetFont(WinID,font$,size,style$)
    END SUB

    SUB BoxDisk(Lft,Rgt,Bas,Top)
        BOX DISK Lft,Rgt,Bas,Top
    END SUB
END MODULE

!  --- End TB5 Cross-platform header and subs ---
!! ---------------------------------------------------------

!! --- Start Unix Header and Subs ---

!LIBRARY "HHLib.trc"
!LET UnixFlag= 1
!ASK PIXELS winWid,winHgt
!LET winLft= 0
!LET winTop= 0
!LET winRgt= winWid-1
!LET winBas= winHgt-1
!SET WINDOW 0,winRgt,winBas,0
!CALL Palette
!
!CALL ToolPanel
!CALL ThisProgram
!
!END
!EXTERNAL
!
!MODULE UnixParts
!    SHARE CharWidth
!
!    SUB SetTextFont(font,size,style$)
!        LET font$= "-adobe-courier-"
!        IF style$= "normal" then
!           LET style$= "medium-r-normal--"
!        ELSE
!           LET style$= "bold-r-normal--"
!        END IF
!        IF size=9 then
!           LET size$= str$(10)
!        ELSE
!           LET size$= str$(size)
!        END IF
!        LET test= SetFont(font$&style$&size$&"*")
!
!        IF size=9 then
!           LET CharWidth=  6
!        ELSE IF size=12 then                     ! numeric output - axis labels
!           LET CharWidth=  7
!        ELSE IF size=14 then                     ! rare
!           LET CharWidth=  8
!        ELSE IF size=18 then                     ! rare
!           LET CharWidth= 10
!        END IF
!    END SUB
!
!    SUB StringWidth(sw$,sl)                      ! string width in pixels
!        ! LET sl= StrWidth(sw$)
!        LET chars= len(sw$)
!        LET sl   = chars*CharWidth
!    END SUB
!
!    SUB SetLineWeight(wgt)
!        ! CALL PenSize(wgt,wgt)
!    END SUB
!
!    SUB BoxDisk(Lft,Rgt,Bas,Top)
!        CALL Fill_Circle(Lft,Rgt,Bas,Top)
!    END SUB
!END MODULE

!! ------ End of TB Unix Header and Subs ------
! -----------------------------------------------------------


! ***

!  Tides. It's tempting to leave phi as it is when k becomes zero, since
!  the limit as k ---> 0 is indeed pi/2. But in keeping with our rule to
!  drop undefined values, I think really the phi indicator should go off
!  when k = 0. We could use >>,<,> keys here rather than "animate."
!  !
!  Here are the alternatives. 
!  H signifies the complex number k/p(i omega) shown in the "Nyquist plot."
!  u means "undefined," its diamond should vanish - rollover should be disabled
!  
!  inf means infinity
!                          t0     P          A       phi     H
!  
!  k = 0, omega = 0        u      inf        u       u       u
!  
!  k = 0, omega not 0      u      1/omega    0       u       0
!  
!  k not 0, omega = 0      u      inf        1       0       1

SUB Tide
    DECLARE PUBLIC black,drkgry,drkmid,midgry,litmid,litgry,white
    DECLARE PUBLIC red,yellow,green,cyan,blue,magenta,colorscheme
    DECLARE PUBLIC axislabelclr,axisclr,slideclr,true,false
    DECLARE PUBLIC workLft,workRgt,workBas,workTop,workMidx

    DECLARE DEF quitWithin, infoWithin

    ! --- help screen array ---

    DIM info1$(1:12)
    MAT READ info1$
    DATA "First Order Frequency Response"
    DATA ""
    DATA "This illustration is paired with Second Order Frequency Response."
    DATA "It illustrates the equation governing e.g. the temperature response, within an insulated container, to a sinusoidally varying ambient temperature."
    DATA "The omega slider at the right controls the circular frequency of the ambient temperature. "
    DATA "The k slider sets the coupling constant k, which is a measure of the degree of insulation."
    DATA "The large time series plane shows the input signal in cyan and the periodic system response in yellow. "
    DATA "The t slider sets the time t, and the corresponding values of the input signal and system response are marked with cyan and yellow diamonds."
    DATA "The left graph plane further illustrates the relationship between the input signal and the system response at the selected time."
    DATA "Time can also be animated using the [>>>] button, or stepped incrementally using the [<] and [>] buttons."
    DATA "The minimal period P, and the time lag t_0, are displayed below the k slider. The time lag is also indicated by a red bar on the main graphing window, connecting the t = 0 line to the first maximum of the system response."
    DATA "A yellow horizontal line indicating the amplitude appears when the cursor is in the main graphing window."

    DIM info2$(1:12)
    MAT READ info2$
    DATA "First Order Frequency Response"
    DATA ""
    DATA "Click the 'Bode and Nyquist Plots' checkbox to open or close three other graph planes along the right side of the illustration."
    DATA "The top plane shows the amplitude response curve, and the amplitude, A, marked by a yellow line segment and diamond, that corresponds to the selected circular frequency, omega."
    DATA "The second plane shows the phase shift curve, with the phase shift or gain, -phi, marked by a green line segment and diamond, that corresponds to the selected circular frequency omega."
    DATA "Rolling the cursor over either one of those windows causes a horizontal rule and the value to appear."
    DATA "The lower plane shows the trajectory of the complex number k/p(i omega) as omega varies from 0 to +infinity, with the value corresponding to the selected circular frequency marked by a yellow diamond."
    DATA "The line segment, or modulus, joining the origin to that complex number is the amplitude A and is marked in yellow."
    DATA "The angle, or argument, of that complex number is the phase gain -phi and is marked in green."
    DATA "The characteristic polynomial of the operator, p(s), is displayed below the complex plane. "
    DATA "Note: These are not quite true Bode and Nyquist plots."
    DATA "A Bode plot graphs log(A) vs log(omega) or - phi vs log(omega)."
    DATA "A Nyquist plot graphs k/p(i omega) as omega ranges from -infinity to +infinity, so it has a portion above the real axis that is symmetric with what is drawn here."


    ! --- color variables ---

    LET phiclr = green
    LET ocnClr = cyan
    LET bayClr = yellow
    LET ampClr = yellow


    ! --- Utility functions ---

    DECLARE DEF clamp,roundn,e
    DEF RadToDeg(rad)= rad*180/pi


    ! ---- Functions ----

    ! phi  = arctan(w/k)
    ! A    = cos(phi)

    ! ocean= b*cos(w*t)
    ! bay  = A*b*cos(w*t-phi)

    ! ocean= cos(w*t)         ! b = 1
    ! bay  = A*cos(w*t-phi)   ! b = 1

    ! as k gets large...
    !   phi approaches 0
    !   A approaches 1
    !   output approaches synchrony with input
    !   bay approaches synchrony with ocean

    !    1. p(s) is the characteristic polynomial. 
    !    2. p(s) = 0 is the characteristic equation. 
    !    3. The roots of p(s) = 0 are the characteristic roots or eigenvalues.


    LET period= 12                               ! 12 hour cycle
    DEF pi2   = 2*pi

    DEF fPhi(k,w)
        IF k>0 then
           LET fPhi= atn(w/k)
        ELSE
           LET fPhi= pi/2
        END IF
    END DEF

    DEF fAp(phi)= cos(phi)                       ! amplitude polar 
    DEF fAc(k,w)= k/sqr(w^2 + k^2)               ! amplitude cartesian

    DEF fP(w)     = pi2/w                        ! P is the period of the tide
    DEF fT0(phi,w)= phi/w

    DEF fOcean(w,t)  = cos(w*t)
    DEF fBay(w,t,phi)= A*cos(w*t-phi)



    ! ----- Design and layout -----

    ! ---------- Graphing plane parameters and methods ----------

    ! ------ w1 plane - time series for ocean and bay heights -------
    ! --- w1 plane data ---

    DECLARE PUBLIC w1Lft,w1Rgt,w1Bas,w1Top,w1Midx,w1Midy
    DECLARE PUBLIC w1fLft,w1fRgt,w1fBas,w1fTop,w1x0,w1y0
    DECLARE PUBLIC w1xFirst, w1xSTik, w1xLTik, w1xLabel, w1xGridstep
    DECLARE PUBLIC w1yFirst, w1ySTik, w1yLTik, w1yLabel, w1yGridstep
    DECLARE PUBLIC w1wWid,w1wHgt,w1fWid,w1fHgt
    DECLARE PUBLIC w1fxRatio,w1fyRatio,w1wxRatio,w1wyRatio,w1Aspect
    DECLARE PUBLIC w1xPiFlag, w1xMult, w1yPiFlag, w1yMult

    LET w1Flag   = 1
    LET w1xPiFlag= 0
    LET w1xMult  = 1
    LET w1yPiFlag= 0
    LET w1yMult  = 1

    LET w1Lft = workLft+220                      ! window bounds
    LET w1Rgt = w1Lft+312
    LET w1Top = workTop+30
    LET w1Bas = w1Top+240

    LET w1fLft=  0                               ! function bounds
    LET w1fRgt= 24
    LET w1fTop=  1.5
    LET w1fBas= -1.5

    LET w1xAx$= "t"                              ! axis labels
    LET w1yAx$= "x"

    LET w1xGridstep= 0                           ! horizontal grid intervals
    LET w1yGridstep= 0                           ! vertical grid intervals

    LET w1xSTik = 0                              ! horizontal axis Tik marks
    LET w1xLTik = 6
    LET w1xLabel= 0
    LET w1xFirst= w1fLft

    LET w1ySTik = 0                              ! vertical axis Tik marks
    LET w1yLTik = 0.5
    LET w1yLabel= 0.5
    LET w1yFirst= w1fBas

    ! --- w1 Plane methods ---

    DECLARE DEF w1Fncx,w1Fncy,w1Wndx,w1Wndy      ! window/function transforms
    DECLARE DEF w1wWithin, w1Within

    CALL w1Variables

    SUB w1Init
        CALL w1DrawPlane(1,1,1)                  ! x axis, y axis, zeroaxes
        CALL SetTextFont(1,12,"bold")

        CALL PlotTextLJ(w1Rgt+8,w1y0+3,w1Xax$,litgry)      ! axis labels
        CALL PlotTextCJ(w1x0,w1Top-10,w1Yax$,yellow)

        LET bas = w1Top-10

        LET txt$= "x' + kx = k cos(wt)"          ! find width for centering
        CALL StringWidth(txt$,sw)
        LET lft = w1Midx-sw/2

        LET txt$= "x' + kx = k cos("
        CALL PlotTextLJ(lft,bas,txt$,yellow)
        CALL StringWidth(txt$,sw)
        LET lft = lft + sw + 2

        CALL SwapOmega(lft,bas,"w","w",yellow)
        LET lft = lft + 12

        LET txt$= "t)"
        CALL PlotTextLJ(lft,bas,txt$,yellow)

        CALL w1KeepGridLayer

        CALL OceanGraph
        CALL w1KeepGraphLayer
    END SUB

    SUB w1ShowGraphLayer2
        BOX SHOW w1graphLayer2$ at w1Lft-5,w1Bas+5
    END SUB

    SUB w1KeepGraphLayer2
        BOX KEEP w1Lft-5,w1Rgt+5,w1Bas+5,w1Top-5 in w1graphLayer2$
    END SUB

    SUB w1KeepTemp
        BOX KEEP w1Lft,w1Rgt,w1Bas,w1Top in w1temp$
    END SUB

    SUB w1ShowTemp
        BOX SHOW w1temp$ at w1Lft,w1Bas
    END SUB


    ! ------- w2 plane - water level cross section ------------------
    ! ------- w2 plane data ---

    DECLARE PUBLIC w2Lft,w2Rgt,w2Bas,w2Top,w2Midx,w2Midy
    DECLARE PUBLIC w2fLft,w2fRgt,w2fBas,w2fTop,w2x0,w2y0
    DECLARE PUBLIC w2xFirst, w2xSTik, w2xLTik, w2xLabel, w2xGridstep
    DECLARE PUBLIC w2yFirst, w2ySTik, w2yLTik, w2yLabel, w2yGridstep
    DECLARE PUBLIC w2wWid,w2wHgt,w2fWid,w2fHgt
    DECLARE PUBLIC w2fxRatio,w2fyRatio,w2wxRatio,w2wyRatio,w2Aspect
    DECLARE PUBLIC w2xPiFlag, w2xMult, w2yPiFlag, w2yMult

    LET w2Flag   = 1
    LET w2xPiFlag= 0
    LET w2xMult  = 1
    LET w2yPiFlag= 0
    LET w2yMult  = 1

    LET w2Lft = worklft + 50                     ! pixel bounds
    LET w2Rgt = w2Lft + 85
    LET w2Top = w1Top
    LET w2Bas = w1Bas

    LET w2fLft=  0                               ! function bounds * pi
    LET w2fRgt=  3
    LET w2fTop=  1.5
    LET w2fBas= -1.5

    LET w2xAx$= "t"                              ! axis labels
    LET w2yAx$= "I"

    LET w2xGridstep= 0                           ! horizontal grid intervals
    LET w2yGridstep= 0                           ! vertical grid intervals

    LET w2xSTik = 0                              ! horizontal axis Tik marks
    LET w2xLTik = 0
    LET w2xLabel= 0
    LET w2xFirst= w2fLft

    LET w2ySTik = 0                              ! vertical axis Tik marks
    LET w2yLTik = 0.5
    LET w2yLabel= 0.5
    LET w2yFirst= w2fBas

    ! --- w2 plane methods ---

    DECLARE DEF w2Fncx, w2Fncy, w2Wndx, w2Wndy   ! window/function transforms
    DECLARE DEF w2wWithin, w2Within

    CALL w2Variables

    LET w2x1= w2Wndx(1)
    LET w2x2= w2Wndx(2)

    SUB w2Init
        CALL w2DrawPlane(1,1,0)                  ! x axis, y axis, zeroaxes

        LET basln = w2Bas + 33
        CALL SetTextFont(1,12,"bold")
        CALL PlotTextLJ(w2Lft,basln   ,"Input Signal",ocnclr)
        CALL PlotTextLJ(w2Lft,basln+20,"Coupling",white)
        CALL PlotTextLJ(w2Lft,basln+40,"Response",bayclr)

        ! draw height range mark

        LET wx2= w2Wndx( 1)
        LET wy3= w2Wndy( 1)
        LET wy4= w2Wndy(-1)
        CALL PlotLine( wx2,wy3, wx2,wy4, ocnClr)
        CALL PlotLine( w2Lft+1,w2y0, w2Rgt-1,w2y0, drkmid)

        CALL w2KeepGridLayer
    END SUB

    SUB w2KeepTemp
        BOX KEEP w2Lft,w2Rgt,w2Bas,w2Top in w2temp$
    END SUB

    SUB w2ShowTemp
        BOX SHOW w2temp$ at w2Lft,w2Bas
    END SUB

    ! --------- w3 Plane - middle omega/phi graph ----------------
    ! ----- w3 plane data -----

    DECLARE PUBLIC w3Lft,w3Rgt,w3Bas,w3Top,w3Midx,w3Midy
    DECLARE PUBLIC w3fLft,w3fRgt,w3fBas,w3fTop,w3x0,w3y0
    DECLARE PUBLIC w3xFirst, w3xSTik, w3xLTik, w3xLabel, w3xGridstep
    DECLARE PUBLIC w3yFirst, w3ySTik, w3yLTik, w3yLabel, w3yGridstep
    DECLARE PUBLIC w3wWid,w3wHgt,w3fWid,w3fHgt
    DECLARE PUBLIC w3fxRatio,w3fyRatio,w3wxRatio,w3wyRatio,w3Aspect
    DECLARE PUBLIC w3xPiFlag, w3xMult, w3yPiFlag, w3yMult

    LET w3Flag   = 1
    LET w3xPiFlag= 0
    LET w3xMult  = 1
    LET w3yPiFlag= 1
    LET w3yMult  = pi

    LET w3Lft = w1Rgt + 100                      ! pixel bounds
    LET w3Rgt = w3Lft + 80
    LET w3Top = w1Midy + 22
    LET w3Bas = w3Top + 80

    LET w3fLft=  0                               ! function bounds
    LET w3fRgt=  1
    LET w3fTop=  0
    LET w3fBas= -0.5

    LET w3xAx$= "w"                              ! axis labels
    LET w3yAx$= ""

    LET w3xGridstep= 0                           ! grid line intervals
    LET w3yGridstep= 0

    LET w3xSTik = 0                              ! axis Tik marks
    LET w3xLTik = 0.5
    LET w3xLabel= 0.5
    LET w3xFirst= w3fLft

    LET w3ySTik  = 0
    LET w3yLTik  = 0.25
    LET w3yLabel = 0.5
    LET w3yFirst = w3fBas

    ! --- w3 plane methods ---

    DECLARE DEF w3Fncx,w3Fncy,w3Wndx,w3Wndy      ! window/function transforms
    DECLARE DEF w3wWithin, w3Within

    CALL w3Variables

    SUB w3Init
        CALL w3DrawPlane(-1,1,1)                 ! x axis, y axis, zeroaxes

        CALL SetTextFont(1,12,"bold")
        CALL PlotTextRJ(w3x0-5,w3Bas+13,"-",phiclr)
        CALL DrawPhi12(w3x0-3,w3Bas+13,phiclr)

        CALL SwapOmega(w3Rgt+8,w3y0+3,"w","w",axislabelclr)

        CALL w3KeepGridLayer
        CALL w3KeepGraphLayer
    END SUB

    SUB w3KeepTemp
        BOX KEEP w3Lft,w3Rgt,w3Bas,w3Top in w3temp$
    END SUB

    SUB w3ShowTemp
        BOX SHOW w3temp$ at w3Lft,w3Bas
    END SUB


    ! ------------- w4 plane - Nyquist real/imaginary ---------------
    ! --- w4 plane data ---

    DECLARE PUBLIC w4Lft,w4Rgt,w4Bas,w4Top,w4Midx,w4Midy
    DECLARE PUBLIC w4fLft,w4fRgt,w4fBas,w4fTop,w4x0,w4y0
    DECLARE PUBLIC w4xFirst, w4xSTik, w4xLTik, w4xLabel, w4xGridstep
    DECLARE PUBLIC w4yFirst, w4ySTik, w4yLTik, w4yLabel, w4yGridstep
    DECLARE PUBLIC w4wWid,w4wHgt,w4fWid,w4fHgt
    DECLARE PUBLIC w4fxRatio,w4fyRatio,w4wxRatio,w4wyRatio,w4Aspect
    DECLARE PUBLIC w4xPiFlag, w4xMult, w4yPiFlag, w4yMult

    LET w4Flag   = 1
    LET w4xPiFlag= 0
    LET w4xMult  = 1
    LET w4yPiFlag= 0
    LET w4yMult  = 1

    LET w4Lft = w3Rgt - 112                      ! pixel bounds
    LET w4Rgt = w4Lft + 112                      ! 160
    LET w4Top = workBas - 145
    LET w4Bas = w4Top + 64

    LET w4fLft= -0.2                             ! function bounds
    LET w4fRgt=  1.2
    LET w4fTop=  0.2
    LET w4fBas= -0.6

    LET w4xAx$= "Re"                             ! axis labels
    LET w4yAx$= "Im"

    LET w4xGridstep= 0                           ! grid line intervals
    LET w4yGridstep= 0

    LET w4xSTik = 0                              ! axis Tik marks
    LET w4xLTik = 0.5
    LET w4xLabel= 0.5
    LET w4xFirst= 0

    LET w4ySTik =  0
    LET w4yLTik =  0.5
    LET w4yLabel=  0.5
    LET w4yFirst= -0.5

    ! --- w4 plane methods ---

    DECLARE DEF w4Fncx,w4Fncy,w4Wndx,w4Wndy      ! window/function transforms
    DECLARE DEF w4wWithin, w4Within

    CALL w4Variables

    LET w4xScale= w4wWid/w4fWid

    SUB w4Init
        CALL w4DrawPlane(-1,1,1)                 ! x axis, y axis, zeroaxes

        CALL SetTextFont(1,12,"bold")
        CALL PlotTextLJ(w4Rgt+8,w4y0+3,w4xAx$,axislabelclr)     ! axis labels
        CALL PlotTextCJ(w4x0,w4Top-17,w4yAx$,axislabelclr)

        LET txt$ = "k / p(iw)"
        CALL StringWidth(txt$,sw)
        CALL SwapOmega(w4Midx-sw/2,w4Bas+20,txt$,"w",white)

        CALL PlotTextCJ(w4Midx,w4Bas+40,"p(s) = s + k",white)   ! axis labels

        CALL w4KeepGridLayer
    END SUB

    SUB w4KeepTemp
        BOX KEEP w4Lft,w4Rgt,w4Bas,w4Top in w4temp$
    END SUB

    SUB w4ShowTemp
        BOX SHOW w4temp$ at w4Lft,w4Bas
    END SUB


    ! ------------ w5 plane - omega/amplitude graph ----------------
    ! --- w5 plane data ---

    DECLARE PUBLIC w5Lft,w5Rgt,w5Bas,w5Top,w5Midx,w5Midy
    DECLARE PUBLIC w5fLft,w5fRgt,w5fBas,w5fTop,w5x0,w5y0
    DECLARE PUBLIC w5xFirst, w5xSTik, w5xLTik, w5xLabel, w5xGridstep
    DECLARE PUBLIC w5yFirst, w5ySTik, w5yLTik, w5yLabel, w5yGridstep
    DECLARE PUBLIC w5wWid,w5wHgt,w5fWid,w5fHgt
    DECLARE PUBLIC w5fxRatio,w5fyRatio,w5wxRatio,w5wyRatio,w5Aspect
    DECLARE PUBLIC w5xPiFlag, w5xMult, w5yPiFlag, w5yMult

    LET w5Flag   = 1
    LET w5xPiFlag= 0
    LET w5xMult  = 1
    LET w5yPiFlag= 0
    LET w5yMult  = 1

    LET w5Lft = w3Lft                            ! pixel bounds
    LET w5Rgt = w5Lft + 80
    LET w5Top = w1Top + 40
    LET w5Bas = w5Top + 80

    LET w5fLft=  0                               ! function bounds
    LET w5fRgt=  1
    LET w5fTop=  1
    LET w5fBas=  0

    LET w5xAx$= "w"                              ! axis labels
    LET w5yAx$= "A"

    LET w5xGridstep= 0                           ! grid line intervals
    LET w5yGridstep= 0

    LET w5xSTik = 0                              ! axis Tik marks
    LET w5xLTik = 0.5
    LET w5xLabel= 0
    LET w5xFirst= w5fLft

    LET w5ySTik = 0
    LET w5yLTik = 0.5
    LET w5yLabel= 0.5
    LET w5yFirst= w5fBas

    ! --- plane 5 methods ---

    DECLARE DEF w5Fncx,w5Fncy,w5Wndx,w5Wndy      ! window/function transforms
    DECLARE DEF w5wWithin, w5Within

    CALL w5Variables

    SUB w5Init
        CALL w5DrawPlane(1,1,1)                  ! x axis, y axis, zeroaxes

        CALL SetTextFont(1,12,"bold")
        CALL SwapOmega(w5Rgt+8,w5y0+3,"w","w",axislabelclr)
        CALL SuperSubScriptCJ(w5x0,w5Top-10,w5yAx$,ampclr)
        CALL w5KeepGridLayer
    END SUB

    SUB w5KeepTemp
        BOX KEEP w5Lft,w5Rgt,w5Bas,w5Top in w5temp$
    END SUB

    SUB w5ShowTemp
        BOX SHOW w5temp$ at w5Lft,w5Bas
    END SUB


    ! ---------- Slider parameters and methods --------
    ! --------------- horizontal sliders --------------

    ! --- h1 horizontal slider - t value ---

    DECLARE PUBLIC h1axis,h1wLft,h1wRgt,h1wBas,h1wTop,h1fLft,h1fRgt
    DECLARE PUBLIC h1name$,h1form$,h1clr,h1First,h1STik,h1LTik,h1Label
    DECLARE PUBLIC h1PiAxis,h1Mult,h1fMin,h1fMax

    DECLARE DEF h1Within                         ! window/function transforms

    LET h1PiAxis= 0
    LET h1Mult  = 1

    LET h1clr   = slideclr
    LET h1name$ = "t"
    LET h1form$ = "--%.#"
    LET h1Places= 1
    LET h1Click = 1

    LET h1axis  = w1Bas + 22
    LET h1wLft  = w1Lft
    LET h1wRgt  = w1Rgt

    LET h1fLft  =  0
    LET h1fRgt  = 24

    LET h1STik  = 1                              ! short tick marks
    LET h1LTik  = 6                              ! long tick marks
    LET h1Label = 6                              ! labels
    LET h1First = h1fLft                         ! first tick mark

    CALL h1SliderVariables

    SUB h1SetZero
        LET t= 0
        CALL h1mark(t)
    END SUB

    DECLARE DEF h1AnimWithin, h1LStpWithin, h1RStpWithin
    DECLARE DEF h1AnimStopWithin

    LET h1AnimStep= 0.2                          ! two pixels

    SUB h1Init
        LET t= 0
        CALL h1DrawSlider(h1name$,t)
        CALL h1AnimMoveStep(1,0,1)               ! animstate,movestate,stepstate
    END SUB


    ! --- h2 Slider - omega value ---

    DECLARE PUBLIC h2axis,h2wLft,h2wRgt,h2wBas,h2wTop,h2fLft,h2fRgt
    DECLARE PUBLIC h2name$,h2form$,h2clr,h2First,h2STik,h2LTik,h2Label
    DECLARE PUBLIC h2PiAxis,h2Mult,h2fMin,h2fMax

    DECLARE DEF h2Within                         ! window/function transforms

    LET h2PiAxis= 0
    LET h2Mult  = 1

    LET h2clr   = slideclr
    LET h2name$ = ""
    LET h2form$ = "%.##"
    LET h2Places= 2
    LET h2Click = 0.25

    LET h2axis  = h1axis
    LET h2wLft  = w5Lft
    LET h2wRgt  = w5Rgt

    LET h2fLft  = 0
    LET h2fRgt  = 1

    LET h2STik  = 0.25
    LET h2LTik  = 0.5
    LET h2Label = 0.5
    LET h2First = h2fLft

    CALL h2SliderVariables

    SUB h2Init
        CALL h2DrawSlider(h2name$,w)
        CALL SwapOmega(h2wLft-16,h2wBas-3,"w","w",white)
    END SUB


    ! --- h3 slider ---

    DECLARE PUBLIC h3axis,h3wLft,h3wRgt,h3wBas,h3wTop,h3fLft,h3fRgt
    DECLARE PUBLIC h3name$,h3form$,h3clr,h3First,h3STik,h3LTik,h3Label
    DECLARE PUBLIC h3PiAxis,h3Mult,h3fMin,h3fMax

    DECLARE DEF h3Within                         ! window/function transforms

    LET h3PiAxis= 0
    LET h3Mult  = 1

    LET h3clr   = slideclr
    LET h3name$ = "k"
    LET h3form$ = "%.##"
    LET h3Places= 2

    LET h3axis  = w1Bas + 110
    LET h3wLft  = h1wLft
    LET h3wRgt  = h3wLft + 200

    LET h3fLft  = 0
    LET h3fRgt  = 2

    LET h3STik  = 0.25
    LET h3LTik  = 0.5
    LET h3Label = 0.5
    LET h3First = h3fLft
    LET h3click = 0.25

    CALL h3SliderVariables

    SUB h3Init
        CALL h3DrawSlider(h3name$,k)
    END SUB


    ! --- text output boxes ---

    ! --- text rectangle 1 - period equation and value ---
    ! Period = 2pi/omega

    LET t1BasLn = w1Bas + 172
    LET t1Lft   = w1Lft + 15
    LET t1Rgt   = w1Rgt - 70
    LET t1Bas   = t1BasLn +  5
    LET t1Top   = t1BasLn - 15
    LET t1Label$= "P = "
    LET t1Clr   = cyan

    SUB t1Label
        CALL SetTextFont(1,12,"bold")
        CALL SuperSubScriptRJ(t1Lft,t1BasLn,t1Label$,t1Clr)
    END SUB

    SUB t1Set
        CALL SetTextFont(1,12,"bold")
        LET lft= t1Lft

        CALL SwapPi(lft,t1BasLn,"2w","w",t1Clr)
        CALL StringWidth("2w",sw)
        LET lft = t1Lft + sw + 3

        CALL DivSign(lft,t1BasLn,t1Clr)
        LET lft = lft + 8

        CALL SwapOmega(lft,t1BasLn,"w","w",t1Clr)
        CALL StringWidth("w",sw)

        LET lft = lft + sw + 6
        CALL PlotTextLJ(lft,t1BasLn,"=",t1Clr)

        LET t1vLft = lft + 12
        LET t1vRgt = t1vLft + 100
    END SUB

    SUB t1Val
        CALL t1ClearVal
        CALL SetTextFont(1,12,"bold")
        IF w=0 then
           CALL DrawInf12(t1vLft,t1BasLn,t1Clr)
        ELSE
           LET t$= trim$(using$("---%.##",fP(w)))
           CALL PlotTextLJ(t1vLft,t1BasLn,t$,t1Clr)
        END IF
    END SUB

    SUB t1Clear
        BOX CLEAR t1Lft-2,t1Rgt,t1Bas,t1Top
    END SUB

    SUB t1ClearVal
        BOX CLEAR t1vLft,t1vRgt,t1Bas,t1Top
    END SUB

    SUB t1Init
        CALL t1Label
        CALL t1Set
        CALL t1Val
    END SUB


    ! --- text rectangle 2 - t_0 phase lag equation and value ---

    ! timelag/phaseshift = ( phi/pi2)*Period = phi/omega

    LET t2BasLn = t1BasLn + 24
    LET t2Lft   = t1Lft
    LET t2Rgt   = t1Rgt
    LET t2Bas   = t2BasLn +  5
    LET t2Top   = t2BasLn - 15
    LET t2Label$= "t_[0] = "

    SUB t2Init
        CALL t2Label
        CALL t2Set
        CALL t2Val
    END SUB

    SUB t2Label
        CALL SetTextFont(1,12,"bold")
        CALL SuperSubScriptRJ(t2Lft,t2BasLn,t2Label$,red)
    END SUB

    SUB t2Set
        CALL SetTextFont(1,12,"bold")
        ! t_0 = (phi/2 pi)P = phi/omega

        CALL StringWidth("(0",sw)
        CALL SwapPhi(t2Lft,t2BasLn,"(0","0",red)
        LET lft = t2Lft + sw + 4

        CALL DivSign(lft,t2BasLn,red)
        LET lft = lft + 8

        CALL SwapPi(lft,t2BasLn,"2p)P","p",red)
        CALL StringWidth("2p)P",sw)
        LET lft = lft + sw + 8

        CALL PlotTextLJ(lft,t2BasLn,"=",red)
        LET lft= lft + 12

        CALL SwapPhi(lft,t2BasLn,"0","0",red)
        CALL StringWidth("0",sw)
        LET lft= lft + sw + 4

        CALL DivSign(lft,t2BasLn,red)
        LET lft= lft + 8

        CALL SwapOmega(lft,t2BasLn,"w","w",red)
        LET lft= lft+16

        CALL PlotTextLJ(lft,t2BasLn,"=",red)

        LET t2vLft = lft + 12
        LET t2vRgt = t2vLft + 40
    END SUB

    SUB t2Val
        CALL t2ClearVal
        IF w>0 and k>0 then
           CALL SetTextFont(1,12,"bold")
           LET t$ = trim$(using$("---%.##",fT0(phi,w)))
           CALL PlotTextLJ(t2vLft,t2BasLn,t$,red)
        END IF
    END SUB

    SUB t2Clear
        BOX CLEAR t2Lft-2,t2Rgt,t2Bas,t2Top
    END SUB

    SUB t2ClearVal
        BOX CLEAR t2vLft-2,t2vRgt,t2Bas,t2Top
    END SUB


    ! --- checkbox - Bode and Nyquist switch ---

    DECLARE PUBLIC cb1Lft,cb1Rgt,cb1Bas,cb1Top,cb1Txt$,cb1Clr,cb1State
    DECLARE DEF cb1Within

    LET cb1Lft = w4Lft-22
    LET cb1Bas = workBas-3
    LET cb1Txt$= "Bode and Nyquist Plots"
    LET cb1Clr = litgry

    CALL cb1Variables

    ! --- end of design and layout ---


    ! --- initialize default parameters ---

    LET w,oldw    = 0.523599
    LET k,oldk    = 0.2
    LET t,oldt    = 0
    LET cb1State  = 0
    LET NyState   = cb1State

    LET phi,oldphi= fPhi(k,w)
    LET A,oldA    = fAp(phi)
    LET timeLag   = fT0(phi,w)

    ! --- Draw the screen ---

    CALL InitScreen
    CALL SetTimer

    SUB InitScreen
        BOX CLEAR worklft,workrgt,workbas,worktop

        CALL w1Init
        CALL w2Init

        CALL t1Init
        CALL t2Init

        CALL BayGraph
        CALL w2BayRangeMark
        CALL TsliderAction(t)

        CALL h1Init
        CALL h2Init
        CALL h3Init
        CALL cb1Init
        CALL InitBodeNyquist(NyState)
    END SUB

    SUB InitBodeNyquist(NyState)
        IF NyState=1 then                        ! add planes to interface
           CALL w5Init                           ! Amplitude
           CALL w3Init                           ! Phase
           CALL w4Init                           ! Nyquist
           CALL w3Curve
           CALL w3Value
           CALL w4Graph
           CALL w4Value
           CALL w5Graph
           CALL w5Value
        ELSE                                     ! remove planes from interface
           BOX CLEAR w3Lft-38,w3Rgt+24,w3Bas+18,w3Top-18
           BOX CLEAR w4Lft-30,w4Rgt+30,w4Bas+44,w4Top-30
           BOX CLEAR w5Lft-26,w5Rgt+24,w5Bas+ 8,w5Top-24
        END IF
    END SUB


    ! ----------------- Event manager -----------------

    DO
       LET ClearFlag1,ClearFlag3= 0
       IF ms<>2 then
          CALL StoreGraphsForClear
          DO
             GET MOUSE: mx,my,ms
             IF w1Within(mx,my)=true or (w5Within(mx,my)=true and NyState=true) then
                CALL w1AmplitudeRollover
             ELSE IF w3Within(mx,my)=true and NyState=true then
                CALL w3PhiRollover
             ELSE IF ClearFlag1=true then
                CALL ClearAmpRollover
             ELSE IF ClearFlag3=true then
                CALL ClearPhiRollover
             END IF
          LOOP until ms=2
       END IF

       ! clean up as needed
       IF ClearFlag1=true then
          CALL ClearAmpRollover
       END IF
       IF ClearFlag3=true then
          CALL ClearPhiRollover
       END IF

       IF h1Within(mx,my)=true then              ! time
          IF my<h1axis then
             CALL h1MouseClick
          ELSE
             CALL h1MouseDrag
          END IF
       ELSE IF h1AnimWithin(mx,my)=true then
          CALL h1AnimButtonUp(ms)
          CALL Animate
       ELSE IF h1LstpWithin(mx,my)=true then
          CALL h1LstpButtonUp(ms)
          CALL h1LeftStep(t,h1AnimStep)
          CALL h1Action
       ELSE IF h1RstpWithin(mx,my)=true then
          CALL h1RstpButtonUp(ms)
          CALL h1RightStep(t,h1AnimStep)
          CALL h1Action

       ELSE IF h2Within(mx,my)=true then         ! omega
          CALL h1SetZero
          IF my<h2axis then
             CALL h2MouseClick
          ELSE
             CALL h2MouseDrag
          END IF
          CALL TsliderAction(t)
       ELSE IF h3Within(mx,my)=true then         ! k
          CALL h1SetZero
          CALL OceanGraph                        ! ocean graph in background
          IF my<h3axis then
             CALL h3MouseClick
          ELSE
             CALL h3MouseDrag
          END IF
          CALL TsliderAction(t)

       ELSE IF cb1Within(mx,my)=true then        ! bode plot switch
          CALL cb1MouseUp(ms)
          LET NyState= cb1State
          CALL InitBodeNyquist(NyState)

       ELSE IF infoWithin(mx,my,ms)=true then    ! info button
          CALL infoButtonUp(ms)
          CALL InfoPage(Info1$)
          CALL InfoPage(Info2$)
          CALL InitScreen
       ELSE IF quitWithin(mx,my,ms)=true then    ! quit button
          CALL quitButtonUp(ms)
          EXIT SUB
       ELSE
          CALL MouseUp(mx,my,ms)
       END IF
    LOOP

    ! --- rollover methods ---

    SUB w1AmplitudeRollover
        IF clearflag3=1 then CALL ClearPhiRollover
        IF ClearFlag1=0 then                     ! don't flicker if on
           LET ra= round(a,5)
           LET wy= w1Wndy(ra)

           SET COLOR yellow
           IF k<>0 or w<>0 then
              PLOT w1Lft,wy; w1Rgt,wy
              PLOT w2Lft,wy; w2Rgt,wy
              IF NyState=1 then PLOT w5Lft,wy; w5Rgt,wy
              CALL PlotTextLJ(w2Rgt+10,wy+4,using$("%.##",ra),yellow)
           END IF
           LET ClearFlag1= 1
        END IF
    END SUB

    SUB w3PhiRollover
        IF clearflag1=1 then CALL ClearAmpRollover
        IF k>0 then                              ! undefined for k=0
           IF ClearFlag3=0 then                  ! don't flicker if on
              LET wy= w3Wndy(-phi)
              CALL PlotLine( w3Lft,wy, w3Rgt,wy, phiclr)
              LET n = round(RadToDeg(-phi))
              CALL PlotDegreesLJ(w3Rgt+10,w3Midy+4,n,phiclr)
              LET ClearFlag3= 1
           END IF
        END IF
    END SUB

    SUB StoreGraphsForClear
        CALL w1KeepTemp
        CALL w2KeepTemp
        IF NyState=1 then                        ! Bode and Nyquest on?
           CALL w5KeepTemp
           CALL w3KeepTemp
        END IF
    END SUB

    SUB ClearAmpRollover
        CALL w1ShowTemp
        CALL w2ShowTemp
        BOX CLEAR w2Rgt+8,w2Rgt+45,w2Bas+5,w2Top-5
        IF NyState=1 then CALL w5ShowTemp
        LET ClearFlag1= 0
    END SUB

    SUB ClearPhiRollover
        CALL w3ShowTemp
        BOX CLEAR w3Rgt+8,workrgt,w3Midy+10,w3Midy-10
        LET ClearFlag3= 0
    END SUB

    ! --- mouse events ---
    ! --- slider 1 - t value ---

    SUB h1MouseClick
        CALL h1GetClickVal(ms,h1Click,t)
        LET oldt= -999                           ! force redraw
        CALL h1Action
    END SUB
    SUB h1MouseDrag
        DO
           CALL h1GetDragVal(ms,h1Places,t)
           CALL h1Action
        LOOP until ms=3
    END SUB

    SUB h1Action
        IF t<>oldt then
           CALL TsliderAction(t)                 ! update heights and time marker
           LET oldt= t
        END IF
    END SUB

    SUB TsliderAction(t)                         ! response to change in t
        ! heights in plane 2 and time marker in plane 1

        LET ocean= fOcean(w,t)
        LET bay  = fBay(w,t,phi)
        LET wyO  = round(w1Wndy(ocean))
        LET wyB  = round(w1Wndy(bay))
        IF abs(wyO-wyB)<=2 then LET wyO= wyB     ! = if within a pixel

        ! time graphs - current time indicator
        ! connects points on the two sinusoidal curves in plane 1

        CALL w1ShowGraphLayer2

        LET wx= w1Wndx(t)
        IF k>0 and w>0 then                      ! points on ocean and bay
           CALL PlotLine( wx,wyO, wx,wyB, white)      ! draw the line
           CALL PlotDiamondClr(wx,wyO,ocnclr)
           CALL PlotDiamondClr(wx,wyB,bayclr)
        ELSE                                     ! ocean only
           CALL PlotDiamondClr(wx,wyO,ocnclr)
        END IF

        ! cross section - ocean and bay heights in plane 2

        CALL w2Redraw
    END SUB

    SUB w2Redraw
        IF w>0 then                              ! tidal oscillation
           CALL w2ShowGraphLayer
           CALL PlotLine( w2Lft+1,wyO, w2x1   ,wyO, ocnclr)     ! ocean
           CALL PlotLine( w2x1   ,wyO, w2x2   ,wyB, white)      ! pass
           CALL PlotLine( w2x2   ,wyB, w2Rgt-1,wyB, bayclr)     ! bay
        ELSE IF w=0 then                         ! no tidal oscillation
           IF k>0 then
              CALL w2Clear
              CALL PlotLine( w2Lft+1,wyO, w2x1   ,wyO, ocnclr)  ! ocean
              CALL PlotLine( w2x1   ,wyO, w2x2   ,wyB, white)   ! pass
              CALL PlotLine( w2x2   ,wyB, w2Rgt-1,wyB, bayclr)  ! bay
           ELSE                                  ! both zero - no bay
              CALL w2Clear
              CALL PlotLine( w2Lft+1,wyO, w2x1   ,wyO, ocnclr)  ! ocean
           END IF
        END IF
    END SUB

    ! --- t slider animation ---

    SUB Animate
        LET start= round(w1Wndx(t))              ! start at time t
        IF start>=w1Rgt then LET start= w1Lft    ! start over at t=0?
        FOR wx= start to w1Rgt
            LET t= w1Fncx(wx)                    ! get t
            CALL TsliderAction(t)                ! update left plane
            CALL h1mark(t)                       ! update slider handle
            GET MOUSE: mx,my,ms                  ! check for mouse event
            IF ms=2 then                         ! mousedown?
               CALL h1AnimStopButtonUp(ms)       ! get mouse up
               EXIT FOR
            END IF
            CALL Delay(1/30)
        NEXT wx
        CALL h1StopButtonClear
    END SUB


    ! --- h2 slider - w (omega) value ---

    SUB h2MouseClick
        CALL h2GetClickVal(ms,h2Click,w)
        CALL h2Action
    END SUB
    SUB h2MouseDrag
        DO
           CALL h2GetDragVal(ms,h2Places,w)
           CALL h2Action
        LOOP until ms=3
    END SUB

    SUB h2Action
        IF w<>oldw then
           LET phi = fPhi(k,w)
           LET A   = fAp(phi)
           IF NyState=1 then                     ! Bode and Nyquest on?
              CALL w5Value
              CALL w4Value
              CALL w3Value
           END IF
           CALL OceanGraph                       ! draw the ocean height curve
           CALL BayGraph                         ! draw the bay height curve
           CALL w2BayRangeMark                   ! plane 2 bay height range mark

           CALL t1Val
           CALL t2Val

           LET oldw= w
        END IF
    END SUB


    ! --- h3 slider - k value ---

    SUB h3MouseClick
        CALL h3GetClickVal(ms,h3Click,k)
        CALL h3Action
    END SUB
    SUB h3MouseDrag
        DO
           CALL h3GetDragVal(ms,h3Places,k)
           CALL h3Action
        LOOP until ms=3
    END SUB

    SUB h3Action
        IF k<>oldk then
           LET phi= fPhi(k,w)
           LET A  = fAp(phi)

           CALL BayGraph                         ! draw the bay height curve
           CALL w2BayRangeMark                   ! plane 2 bay height range mark

           IF NyState=1 then                     ! Bode and Nyquest on?
              CALL w3Curve
              CALL w3Value
              CALL w5Graph
              CALL w5Value
              CALL w4Value
           END IF
           CALL t2Val

           LET oldk= k
        END IF
    END SUB

    ! --- plane 2 cross section - water levels ---

    SUB w2BayRangeMark                           ! ocean and bay cross section
        IF w>0 then
           CALL w2ShowGridLayer
        ELSE
           CALL w2Clear
        END IF
        IF w>0 then                              ! show range line for bay height
           LET wx2= w2Wndx( 2)
           LET wy1= w2Wndy( A)
           LET wy2= w2Wndy(-A)
           CALL PlotLine(wx2,wy1, wx2,wy2, bayclr)
        END IF
        CALL w2KeepGraphLayer
    END SUB

    ! --- plane 3 - phi from omega ---

    SUB w3Curve                                  ! phi graph
        CALL w3ShowGridLayer
        IF k>0 then
           SET COLOR litmid
           FOR wx= w3Lft to w3Rgt
               LET x = w3Fncx(wx)                ! omega
               LET y = -atn(x/k)                 ! phi
               LET wy= w3Wndy(y)
               IF wy>w3Top and wy<w3Bas then
                  PLOT wx,wy;
               ELSE
                  PLOT
               END IF
           NEXT wx
           PLOT
        END IF
        CALL w3KeepGraphLayer
    END SUB

    SUB w3Value                                  ! phi point and line
        CALL w3ShowGraphLayer
        IF k>0 then
           CALL w3MathToPixels(w,-fPhi(k,w),wx,wy)
           CALL PlotLine(wx,wy,wx,w3Top+1,phiclr)
           CALL PlotDiamondClr(wx,wy,phiclr)
        END IF
    END SUB

    ! ---- plane 4 - complex plane - Nyquist plot ----

    SUB w4Graph                                  ! Nyquist curve
        CALL w4ShowGridLayer
        IF k>0 and w>0 then
           SET COLOR litmid
           LET r,xc= 0.5
           LET stp = 1/(w4xScale*r)              ! angle per pixel
           FOR ang= 0 to pi step stp
               CALL PolarToCartesian(r,ang,x,y)  ! A,phi to x,y
               CALL w4MathToPixels(xc+x,-y,wx,wy)
               PLOT wx,wy;
           NEXT ang
           PLOT
        END IF
        CALL w4KeepGraphLayer
    END SUB

    SUB w4Value
        CALL w4ShowGraphLayer
        CALL PolarToCartesian(A,-phi,x,y)

        ! set the line and point

        CALL w4MathToPixels(x,y,wx,wy)
        CALL PlotLine(w4x0,w4y0, wx,wy,ampclr)
        CALL PlotDiamondClr(wx,wy,ampclr)

        ! draw the green arc

        IF k>0 or w>0 then
           SET COLOR phiclr
           LET arcrad= sqr(x*x+y*y)/2
           IF arcrad>0 then
              LET stp= 1/(w4xScale*arcrad)       ! angle per pixel
              FOR ang= 0 to phi step stp
                  CALL PolarToCartesian(arcrad,ang,x,y)
                  PLOT w4Wndx(x),w4Wndy(-y);
              NEXT ang
              PLOT
           END IF
        ELSE
           CALL w4ShowGridLayer
        END IF
    END SUB

    ! ---- w5 plane - amplitude graph ---

    SUB w5Graph
        ! A = k/sqrt(omega^2 + k^2)
        SET COLOR litmid
        CALL w5ShowGridLayer
        IF w>0 or k>0 then
           FOR wx= w5lft+1 to w5rgt
               LET x = w5Fncx(wx)                ! omega
               LET y = fAc(k,x)                  ! amplitude
               LET wy= w5Wndy(y)
               PLOT wx,wy;
           NEXT wx
           PLOT
        END IF
        CALL w5KeepGraphLayer
    END SUB

    SUB w5Value
        CALL w5ShowGraphLayer
        IF w>0 or k>0 then
           LET y = fAc(k,w)                      ! amplitude
           LET wx= w5Wndx(w)
           LET wy= w5Wndy(y)
           CALL PlotDiamondClr(wx,wy,ampclr)
           CALL PlotLine( wx,wy, wx,w5Bas-1, ampclr)
        END IF
    END SUB

    ! --- w1 plane graphing methods - sinusoidal ocean and bay curves ---

    SUB OceanGraph
        CALL w1ShowGridLayer
        IF w>0 then                              ! show period of ocean tides
           LET period= pi2/w                     ! period value
           LET wpx   = w1Wndx(period)
           IF wpx<w1Rgt then                     ! if period is in window
              CALL PlotLine(wpx,w1Top+1, wpx,w1Bas-1, drkmid)
           END IF
        END IF

        SET COLOR ocnclr                         ! draw the ocean height curve
        FOR wx= w1Lft to w1Rgt
            LET ot   = w1Fncx(wx)                ! find t
            LET ocean= fOcean(w,ot)              ! get ocean height
            PLOT wx,w1Wndy(ocean);
        NEXT wx
        PLOT
        CALL w1KeepGraphLayer
    END SUB

    SUB BayGraph
        CALL w1ShowGraphLayer

        IF k>0 and w>0 then                      ! if k AND omega are non-zero
           LET lag= fT0(phi,w)                   ! display the lag or phase shift
           LET wx = w1Wndx(lag)
           IF wx<w1Rgt then                      ! in graph window?
              LET wy= w1Wndy(-1.25)              ! vertical position of lag line
              CALL PlotLine( w1x0,wy, wx,wy, red)     ! horizontal red
              LET wy= w1Wndy(fBay(w,lag,phi))
              CALL PlotLine( wx,wy, wx,w1Bas-1, drkmid)    ! vertical gray
           END IF
        END IF

        IF k>0 or w>0 then                       ! if k OR omega is non-zero
           SET COLOR bayclr                      ! draw the bay curve
           FOR wx= w1Lft to w1Rgt
               LET bt = w1Fncx(wx)               ! find t
               LET bay= fBay(w,bt,phi)           ! get bay height
               PLOT wx,w1Wndy(bay);
           NEXT wx
           PLOT
        END IF
        CALL w1KeepGraphLayer2
    END SUB
END SUB

! --- end of tide code -------------------