git.maemo.org Git - gnuplot/blob - demo/multiplt.dem

   1 #   GNUPLOT v3.6 beta multiplot script file
   2 #
   3 # Second Order System Characteristics
   4 #
   5 #       D**2 + 2*zeta*wn*D + (wn**2)y = (wn**2)*x
   6 #
   7 #   x           input variable
   8 #   y           output variable
   9 #   w           frequency ratio (w/wn)
  10 #   wn          natural frequency
  11 #   wd          damped natural frequency
  12 #   zeta        damping ratio
  13 #   mag(w)      amplitude response
  14 #   phi(w)      phase response
  15 #   wdwn        damped natural frequency ratio
  16 #   wnt         normalized time
  17 #
  18 #   Plots:
  19 #   Frequency domain    magnitude response
  20 #                           phase response
  21 #
  22 #   Time domain         unit step response
  23 #                       unit impulse response
  24 #
  25 #
  26 # Created by: W. D. Kirby email: wdkirby@ix.netcom.com
  27 # Date: 1/18/96
  28 # Released to the public domain with no warranty of any kind
  29 #
  30 reset
  31 set style function lines
  32 set size 1.0, 1.0
  33 set origin 0.0, 0.0
  34 set multiplot
  35 set size 0.5,0.5
  36 set origin 0.0,0.5
  37 set grid
  38 unset key
  39 set angles radians
  40 set samples 250
  41 #  Plot Magnitude Response
  42 set title "Second Order System Transfer Function - Magnitude"
  43 mag(w) = -10*log10( (1-w**2)**2 + 4*(zeta*w)**2)
  44 set dummy w
  45 set logscale x
  46 set xlabel "Frequency (w/wn)"
  47 set ylabel "Magnitude (dB)" offset 1,0
  48 set label 1 "Damping =.1,.2,.3,.4,.5,.707,1.0,2.0" at .14,17
  49 set xrange [.1:10]
  50 set yrange [-40:20]
  51 plot \
  52   zeta=.1,mag(w), \
  53   zeta=.2,mag(w), \
  54   zeta=.3,mag(w), \
  55   zeta=.4,mag(w), \
  56   zeta=.5,mag(w), \
  57   zeta=.707,mag(w), \
  58   zeta=1.0,mag(w), \
  59   zeta=2.0,mag(w),-6
  60 #  Plot Phase Response
  61 set size 0.5,0.5
  62 set origin 0.0,0.0
  63 set title "Second Order System Transfer Function - Phase"
  64 set label 1 ""
  65 set ylabel "Phase (deg)" offset 1,0
  66 set ytics -180, 30, 0
  67 set yrange [-180:0]
  68 tmp(w) = (-180/pi)*atan( 2*zeta*w/(1-w**2) )
  69 # Fix for atan function wrap problem
  70 tmp1(w)= w<1?tmp(w):(tmp(w)-180)
  71 phi(w)=zeta==1?(-2*(180/pi)*atan(w)):tmp1(w)
  72 plot \
  73   zeta=.1,phi(w), \
  74   zeta=.2,phi(w), \
  75   zeta=.3,phi(w), \
  76   zeta=.4,phi(w), \
  77   zeta=.5,phi(w), \
  78   zeta=.707,phi(w), \
  79   zeta=1,phi(w), \
  80   zeta=2.0,phi(w), \
  81   -90
  82 #  Plot Step Response
  83 set size 0.5,0.5
  84 set origin 0.5,0.5
  85 set dummy wnt
  86 unset logscale x
  87 set title "Second Order System - Unit Step Response"
  88 set ylabel "Amplitude y(wnt)" offset 1,0
  89 set xlabel "Normalized Time (wnt)"
  90 set xrange [0:20]
  91 set xtics 0,5,20
  92 set yrange [0:2.0]
  93 set ytics 0, .5, 2.0
  94 set mytics 5
  95 set mxtics 10
  96 wdwn(zeta)=sqrt(1-zeta**2)
  97 shift(zeta) = atan(wdwn(zeta)/zeta)
  98 alpha(zeta)=zeta>1?sqrt(zeta**2-1.0):0
  99 tau1(zeta)=1/(zeta-alpha(zeta))
 100 tau2(zeta)=1/(zeta+alpha(zeta))
 101 c1(zeta)=(zeta + alpha(zeta))/(2*alpha(zeta))
 102 c2(zeta)=c1(zeta)-1
 103 y1(wnt)=zeta==1?1 - exp(-wnt)*(wnt + 1):0
 104 y2(wnt)=zeta<1?(1 - (exp(-zeta*wnt)/wdwn(zeta))*sin(wdwn(zeta)*wnt + shift(zeta))):y1(wnt)
 105 y(wnt)=zeta>1?1-c1(zeta)*exp(-wnt/tau1(zeta))+c2(zeta)*exp(-wnt/tau2(zeta)):y2(wnt)
 106 plot \
 107   zeta=.1,y(wnt), \
 108   zeta=.2,y(wnt), \
 109   zeta=.3,y(wnt), \
 110   zeta=.4,y(wnt), \
 111   zeta=.5,y(wnt), \
 112   zeta=.707,y(wnt), \
 113   zeta=1,y(wnt), \
 114   zeta=2,y(wnt)
 115 #
 116 #  Plot Impulse Response
 117 set origin .5,0.
 118 set title "Second Order System - Unit Impulse Response"
 119 y(wnt)=exp(-zeta*wnt) * sin(wdwn(zeta)*wnt) / wdwn(zeta)
 120 set yrange [-1. :1.]
 121 set ytics -1,.5,1.
 122 plot \
 123   zeta=.1,y(wnt), \
 124   zeta=.2,y(wnt), \
 125   zeta=.3,y(wnt), \
 126   zeta=.4,y(wnt), \
 127   zeta=.5,y(wnt), \
 128   zeta=.707,y(wnt), \
 129   zeta=1,y(wnt), \
 130   zeta=2,y(wnt)
 131 unset multiplot
 132 pause -1 "Hit return to continue"
 133 #
 134 #  Clean up: reset parameter defaults
 135 #
 136 reset