Dimpulse.java
/*
* $Id: Dimpulse.java,v 1.29 2008/03/24 23:45:43 koga Exp $
*
* Copyright (C) 2004 Koga Laboratory. All rights reserved.
*/
package org.mklab.tool.control;
import java.io.IOException;
import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;
import org.mklab.nfc.matrix.DoubleMatrix;
import org.mklab.nfc.matrix.DoubleRationalPolynomialMatrix;
import org.mklab.nfc.scalar.DoubleRationalPolynomial;
import org.mklab.tool.graph.gnuplot.Canvas;
import org.mklab.tool.graph.gnuplot.Gnuplot;
/**
* 離散時間線形システムのインパルス応答を求めるクラスです。
*
* <p>Impulse response of discrete-time linear systems
*
* @author koga
* @version $Revision: 1.29 $
* @see org.mklab.tool.control.Impulse
* @see org.mklab.tool.control.Dstep
*/
public class Dimpulse {
/**
* 離散時間システム
*
* <pre><code>
*
* x[n+1] = Ax[n] + Bu[n] y[n] = Cx[n] + Du[n]
*
* </code></pre>
*
* の <code>iu</code> 番目の入力に単位インパルスを加えたときの応答を計算します。
*
* <p><code>iu = 0</code> のとき、
*
* <pre><code>
*
* [[Y for 1st input] [[X for 1st input] [Y for 2nd input] [X for 2nd input] [...............] [...............] [Y for m'th input]] と [X for m'th input]].
*
* </code></pre>
*
* を返します。
*
* @param A A行列
* @param B B行列
* @param C C行列
* @param D D行列
* @param inputNumber 入力番号
* @param seriesSize 応答の個数
* @return {YY, XX} (出力と状態の応答) response
*/
public static List<DoubleMatrix> dimpulse(DoubleMatrix A, DoubleMatrix B, DoubleMatrix C, DoubleMatrix D, int inputNumber, int seriesSize) {
String message;
if ((message = Abcdchk.abcdchk(A, B, C, D)).length() > 0) {
throw new IllegalArgumentException(message);
}
DoubleMatrix XX, YY;
if (inputNumber == 0) {
int inputSize = B.getColumnSize();
YY = A.createZero(C.getRowSize() * inputSize, seriesSize);
XX = A.createZero(A.getRowSize() * inputSize, seriesSize);
for (int i = 1; i <= inputSize; i++) {
List<DoubleMatrix> yx = Dlsim.dlsim(A, B.getColumnVector(i), C, D.getColumnVector(i), A.createOnes(1).appendRight(A.createZero(1, seriesSize - 1)));
DoubleMatrix Y = yx.get(0);
DoubleMatrix X = yx.get(1);
YY.setSubMatrix(i, 1, Y, Y);
XX.setSubMatrix(i, 1, X, X);
}
} else {
List<DoubleMatrix> yx = Dlsim.dlsim(A, B.getColumnVector(inputNumber), C, D.getColumnVector(inputNumber), A.createOnes(1).appendRight(A.createZero(1, seriesSize - 1)));
YY = yx.get(0);
XX = yx.get(1);
}
return new ArrayList<>(Arrays.asList(new DoubleMatrix[] {YY, XX}));
//return new MatxList(new Object[] {YY, XX});
}
/**
* 単位インパルス応答を計算する
*
* @param g 伝達関数
* @param seriesSize 応答の個数
* @return {YY, XX} (出力と状態の応答) response
*/
public static List<DoubleMatrix> dimpulse(DoubleRationalPolynomial g, int seriesSize) {
List<DoubleMatrix> abcd = Tfn2ss.tfn2ss(g);
DoubleMatrix A = abcd.get(0);
DoubleMatrix B = abcd.get(1);
DoubleMatrix C = abcd.get(2);
DoubleMatrix D = abcd.get(3);
return dimpulse(A, B, C, D, 1, seriesSize);
}
/**
* 単位インパルス応答を計算する
*
* @param G 伝達関数行列
* @param inputNumber 入力番号
* @param seriesSize 応答の個数
* @return {YY, XX} (出力と状態の応答) response
*/
public static List<DoubleMatrix> dimpulse(DoubleRationalPolynomialMatrix G, int inputNumber, int seriesSize) {
if (inputNumber == 0) {
List<DoubleMatrix> abcd = Tfm2ss.tfm2ss(G);
DoubleMatrix A = abcd.get(0);
DoubleMatrix B = abcd.get(1);
DoubleMatrix C = abcd.get(2);
DoubleMatrix D = abcd.get(3);
return dimpulse(A, B, C, D, inputNumber, seriesSize);
}
List<DoubleMatrix> abcd = Tfm2ss.tfm2ss(G, inputNumber);
DoubleMatrix A = abcd.get(0);
DoubleMatrix B = abcd.get(1);
DoubleMatrix C = abcd.get(2);
DoubleMatrix D = abcd.get(3);
return dimpulse(A, B, C, D, 1, seriesSize);
}
/**
* 単位インパルス応答を計算する
*
* @param num 伝達関数の分子多項式の係数
* @param den 伝達関数の分母多項式の係数
* @param seriesSize 応答の個数
* @return {YY, XX} (出力と状態の応答) response
*/
public static List<DoubleMatrix> dimpulse(DoubleMatrix num, DoubleMatrix den, int seriesSize) {
List<DoubleMatrix> abcd = Tf2ss.tf2ss(num, den);
DoubleMatrix A = abcd.get(0);
DoubleMatrix B = abcd.get(1);
DoubleMatrix C = abcd.get(2);
DoubleMatrix D = abcd.get(3);
return dimpulse(A, B, C, D, 1, seriesSize);
}
/**
* @param g 伝達関数
* @param seriesSize 入力番号
* @return Gnuplot
* @throws IOException gnuplotプロセスを起動できない場合
*/
public static Gnuplot plot(DoubleRationalPolynomial g, int seriesSize) throws IOException {
final Gnuplot gp = new Gnuplot();
return plot(gp, g, seriesSize);
}
/**
* 単位インパルス応答をプロットする
*
* @param gnuplot gnuplot
* @param g 伝達関数
* @param seriesSize 応答の個数
* @return Gnuplot
*/
@SuppressWarnings("nls")
public static Gnuplot plot(Gnuplot gnuplot, DoubleRationalPolynomial g, int seriesSize) {
List<DoubleMatrix> yx = dimpulse(g, seriesSize);
DoubleMatrix Y = yx.get(0);
Canvas canvas = gnuplot.createCanvas();
canvas.setGridVisible(true);
canvas.setXLabel("k");
canvas.setYLabel("y");
canvas.plot(DoubleMatrix.series(1, seriesSize, 1), Y, new String[] {"y"});
return gnuplot;
}
/**
* @param numerator 分子多項式の係数
* @param denominator 分母多項式の係数
* @param seriesSize 応答の個数
* @return Gnuplot
* @throws IOException gnuplotプロセスを起動できない場合
*/
public static Gnuplot plot(DoubleMatrix numerator, DoubleMatrix denominator, int seriesSize) throws IOException {
Gnuplot gp = new Gnuplot();
return plot(gp, numerator, denominator, seriesSize);
}
/**
* 単位インパルス応答をプロットする
*
* @param gnuplot gnuplot
* @param numerator 分子多項式の係数
* @param denominator 分母多項式の係数
* @param seriesSize 応答の個数
* @return Gnuplot
*/
public static Gnuplot plot(Gnuplot gnuplot, DoubleMatrix numerator, DoubleMatrix denominator, int seriesSize) {
DoubleRationalPolynomial g = Tf2tfn.tf2tfn(numerator, denominator);
return plot(gnuplot, g, seriesSize);
}
/**
* @param G 伝達関数行列
* @param inputNumber 入力番号
* @param seriesSize 応答の個数
* @return Gnuplot
* @throws IOException gnuplotプロセスを起動できない場合
*/
public static Gnuplot plot(DoubleRationalPolynomialMatrix G, int inputNumber, int seriesSize) throws IOException {
Gnuplot gp = new Gnuplot();
return plot(gp, G, inputNumber, seriesSize);
}
/**
* 単位インパルス応答をプロットする
*
* @param gnuplot gnuplot
* @param G 伝達関数行列
* @param inputNumber 入力番号
* @param seriesSize 応答の個数
* @return Gnuplot
*/
@SuppressWarnings("nls")
public static Gnuplot plot(Gnuplot gnuplot, DoubleRationalPolynomialMatrix G, int inputNumber, int seriesSize) {
List<DoubleMatrix> yx = dimpulse(G, inputNumber, seriesSize);
DoubleMatrix Y = yx.get(0);
int outputSize = G.getRowSize();
String[] yy = new String[outputSize];
for (int i = 1; i <= outputSize; i++) {
if (outputSize == 1) {
yy[i - 1] = "y";
} else {
yy[i - 1] = "y" + i;
}
}
Canvas canvas = gnuplot.createCanvas();
canvas.setGridVisible(true);
canvas.setXLabel("k");
canvas.setYLabel("y");
canvas.plot(DoubleMatrix.series(1, seriesSize, 1), Y, yy);
return gnuplot;
}
/**
* @param A A行列
* @param B B行列
* @param C C行列
* @param D D行列
* @param inputNumber 入力番号
* @param seriesSize 応答の個数
* @return Gnuplot
* @throws IOException gnuplotプロセスを起動できない場合
*/
public static Gnuplot plot(DoubleMatrix A, DoubleMatrix B, DoubleMatrix C, DoubleMatrix D, int inputNumber, int seriesSize) throws IOException {
Gnuplot gp = new Gnuplot();
return plot(gp, A, B, C, D, inputNumber, seriesSize);
}
/**
* 単位インパルス応答をプロットする
*
* @param gnuplot gnuplot
* @param A A行列
* @param B B行列
* @param C C行列
* @param D D行列
* @param inputNumber 入力番号
* @param seriesSize 応答の個数
* @return Gnuplot
*/
@SuppressWarnings("nls")
public static Gnuplot plot(Gnuplot gnuplot, DoubleMatrix A, DoubleMatrix B, DoubleMatrix C, DoubleMatrix D, int inputNumber, int seriesSize) {
List<DoubleMatrix> yx = dimpulse(A, B, C, D, inputNumber, seriesSize);
DoubleMatrix Y = yx.get(0);
DoubleMatrix X = yx.get(1);
int outputSize = C.getRowSize();
String[] yy = new String[outputSize];
for (int i = 1; i <= outputSize; i++) {
if (outputSize == 1) {
yy[i - 1] = "y";
} else {
yy[i - 1] = "y" + i;
}
}
gnuplot.createCanvas(2, 1);
Canvas canvas1 = gnuplot.getCanvas(0, 0);
canvas1.setGridVisible(true);
canvas1.setXLabel("k");
canvas1.setYLabel("y");
canvas1.plot(DoubleMatrix.series(1, seriesSize, 1), Y, yy);
int stateSize = A.getRowSize();
String[] xx = new String[stateSize];
for (int i = 1; i <= stateSize; i++) {
xx[i - 1] = "x" + i;
}
Canvas canvas2 = gnuplot.getCanvas(1, 0);
canvas2.setGridVisible(true);
canvas2.setXLabel("k");
canvas2.setYLabel("x");
canvas2.plot(DoubleMatrix.series(1, seriesSize, 1), X, xx);
return gnuplot;
}
}