examples/pid_repsonse.py

   1 #!/usr/bin/env python
   2 # Copyright (C) 2011-2012 W. Trevor King <wking@tremily.us>
   3 #
   4 # This file is part of pypid.
   5 #
   6 # pypid is free software: you can redistribute it and/or modify it under the
   7 # terms of the GNU General Public License as published by the Free Software
   8 # Foundation, either version 3 of the License, or (at your option) any later
   9 # version.
  10 #
  11 # pypid is distributed in the hope that it will be useful, but WITHOUT ANY
  12 # WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR
  13 # A PARTICULAR PURPOSE.  See the GNU General Public License for more details.
  14 #
  15 # You should have received a copy of the GNU General Public License along with
  16 # pypid.  If not, see <http://www.gnu.org/licenses/>.
  17
  18
  19 from argparse import ArgumentParser
  20 from sys import stdout
  21 from time import sleep
  22
  23 try:
  24     from matplotlib import pyplot
  25     from numpy import loadtxt
  26 except (ImportError,RuntimeError), e:
  27     pyplot = None
  28     loadtxt = None
  29     plot_import_error = e
  30
  31 from pypid.backend.test import TestBackend
  32 from pypid.rules import ziegler_nichols_step_response
  33
  34
  35 parser = ArgumentParser(description='Simulate a step response.')
  36 parser.add_argument(
  37     '-K', '--process-gain', metavar='GAIN', type=float, default=1,
  38     help='process gain (PV-units over MV-units)',)
  39 parser.add_argument(
  40     '-L', '--dead-time', metavar='TIME', type=float, default=1,
  41     help='system dead time (lag)')
  42 parser.add_argument(
  43     '-T', '--decay-time', metavar='TIME', type=float, default=1,
  44     help='exponential decay timescale')
  45 parser.add_argument(
  46     '-P', '--proportional', metavar='GAIN', type=float, default=None,
  47     help='process gain (output units over input units)',)
  48 parser.add_argument(
  49     '-I', '--integral', metavar='TIME', type=float, default=None,
  50     help='intergral gain timescale')
  51 parser.add_argument(
  52     '-D', '--derivative', metavar='TIME', type=float, default=None,
  53     help='derivative gain timescale')
  54 parser.add_argument(
  55     '-M', '--max-mv', metavar='MV', type=float, default=100.,
  56     help='maximum manipulated variable')
  57 parser.add_argument(
  58     '-A', '--tuning-algorithm', metavar='TUNER', default=None,
  59     choices=['ZN'], help='step tuning algorithm')
  60 parser.add_argument(
  61     '-m', '--mode', metavar='MODE', default='PID',
  62     choices=['P', 'PI', 'PID'], help='controller mode')
  63 parser.add_argument(
  64     '-t', '--time', metavar='TIME', type=float, default=10.,
  65     help='simulation time')
  66 parser.add_argument(
  67     '-o', '--output', default='-', help='output log file')
  68 parser.add_argument(
  69     '-p', '--plot', action='store_true', default=False,
  70     help='plot the repsonse')
  71
  72 args = parser.parse_args()
  73
  74 if args.plot and not (pyplot and loadtxt) :
  75     raise plot_import_error
  76
  77 if args.output == '-':
  78     log_stream = stdout
  79     if args.plot:
  80         raise ValueError('can only plot when outputing to a file')
  81 else:
  82     log_stream = open(args.output, 'w')
  83
  84 K = args.process_gain
  85 L = args.dead_time
  86 T = args.decay_time
  87
  88 p,i,d = (0, float('inf'), 0)
  89 if args.tuning_algorithm == 'ZN':
  90     p,i,d = ziegler_nichols_step_response(
  91         process_gain=K, dead_time=L, decay_time=T, mode=args.mode)
  92 else:
  93     if args.proportional:
  94         p = args.proportional
  95     if args.integral:
  96         i = args.integral
  97     if args.derivative:
  98         d = args.derivative
  99
 100 b = TestBackend(
 101     process_gain=K, dead_time=L, decay_time=T, max_mv=args.max_mv,
 102     log_stream=log_stream)
 103 try:
 104     b.set_up_gains(proportional=p, integral=i, derivative=d)
 105     b.set_down_gains(proportional=p, integral=i, derivative=d)
 106     b.set_setpoint(1.0)
 107     sleep(args.time)
 108 finally:
 109     b.cleanup();
 110     if args.output != '-':
 111         log_stream.close()
 112
 113 if args.plot:
 114     header = open(args.output, 'r').readline()
 115     label = header.strip('#\n').split('\t')
 116     data = loadtxt(args.output)
 117     times = data[:,0] - data[0,0]
 118     pyplot.hold(True)
 119     subplot = 1
 120     for i in range(1, len(label)):
 121         if i in [1, 4, 6]:
 122             if i:
 123                 pyplot.legend(loc='best')  # add legend to previous subplot
 124             pyplot.subplot(3, 1, subplot)
 125             subplot += 1
 126         pyplot.plot(times, data[:,i], '.', label=label[i])
 127     pyplot.legend(loc='best')
 128     pyplot.show()