| 静岡理工科大学 | 菅沼ホーム | 目次 | 索引 |
周波数１ ゲイン１（または，位相１） 周波数２ ゲイン２（または，位相２） ・・・・・・・
/********************************/
/* 伝達関数のゲインと位相の計算 */
/* coded by Y.Suganuma */
/********************************/
#include <stdio.h>
#include <stdlib.h>
#include <iostream>
using namespace std;
/****************************/
/* クラスComplex */
/* coded by Y.Suganuma */
/****************************/
class Complex {
double r; // 実部
double i; // 虚部
public:
Complex (double a = 0.0, double b = 0.0); // コンストラクタ
friend double abs (Complex a); // 絶対値
friend double angle (Complex a); // 角度
friend Complex pow (Complex a, int n); // べき乗
friend Complex operator +(Complex a, Complex b); // ＋のオーバーロード
friend Complex operator -(Complex a, Complex b); // ーのオーバーロード
friend Complex operator *(Complex a, Complex b); // ＊のオーバーロード
friend Complex operator /(Complex a, Complex b); // ／のオーバーロード
friend ostream& operator << (ostream&, Complex); // ＜＜のオーバーロード
};
/******************/
/* コンストラクタ */
/* a : 実部 */
/* b : 虚部 */
/******************/
Complex::Complex(double a, double b)
{
r = a;
i = b;
}
/******************/
/* 複素数の絶対値 */
/******************/
double abs(Complex a)
{
double x;
x = sqrt(a.r * a.r + a.i * a.i);
return x;
}
/****************/
/* 複素数の角度 */
/****************/
double angle(Complex a)
{
double x;
if (a.r == 0.0 && a.i == 0.0)
x = 0.0;
else
x = atan2(a.i, a.r);
return x;
}
/****************/
/* 複素数のｎ乗 */
/****************/
Complex pow(Complex a, int n)
{
int i1;
Complex c(1, 0);
if (n >= 0) {
for (i1 = 0; i1 < n; i1++)
c = c * a;
}
else {
for (i1 = 0; i1 < -n; i1++)
c = c * a;
c = 1.0 / c;
}
return c;
}
/**********************/
/* ＋のオーバーロード */
/**********************/
Complex operator +(Complex a, Complex b)
{
Complex c;
c.r = a.r + b.r;
c.i = a.i + b.i;
return c;
}
/**********************/
/* ーのオーバーロード */
/**********************/
Complex operator -(Complex a, Complex b)
{
Complex c;
c.r = a.r - b.r;
c.i = a.i - b.i;
return c;
}
/**********************/
/* ＊のオーバーロード */
/**********************/
Complex operator *(Complex a, Complex b)
{
Complex c;
c.r = a.r * b.r - a.i * b.i;
c.i = a.i * b.r + a.r * b.i;
return c;
}
/**********************/
/* ／のオーバーロード */
/**********************/
Complex operator /(Complex a, Complex b)
{
double x;
Complex c;
x = 1.0 / (b.r * b.r + b.i * b.i);
c.r = (a.r * b.r + a.i * b.i) * x;
c.i = (a.i * b.r - a.r * b.i) * x;
return c;
}
/************************/
/* ＜＜のオーバーロード */
/************************/
ostream& operator << (ostream& stream, Complex a)
{
stream << "(" << a.r << ", " << a.i << ")";
return stream;
}
Complex G_s(double, int, int, double *, int, int, double *);
Complex value(Complex, int, int, double *);
/*******************/
/* main プログラム */
/*******************/
int main()
{
double f, ff, fl, fu, h, *bo, *si, gain, phase = 0.0, x, uc;
int i1, k, k1, ms, mb, m, n;
char o_fg[100], o_fp[100];
Complex g;
FILE *og, *op;
uc = 90.0 / asin(1.0);
/*
データの入力
*/
// 出力ファイル名の入力とファイルのオープン
scanf("%*s %s %*s %s", o_fg, o_fp);
og = fopen(o_fg, "w");
op = fopen(o_fp, "w");
// 伝達関数データ
scanf("%*s %lf %lf %*s %d", &fl, &fu, &k);
// 分子
scanf("%*s %d %*s %d %*s", &ms, &m);
// 因数分解されていない
if (ms == 0) {
k1 = m + 1;
si = new double [k1];
for (i1 = 0; i1 < k1; i1++)
scanf("%lf", &si[i1]);
}
// 因数分解されている
else {
k1 = (m == 0) ? 1 : 2 * m;
si = new double [k1];
for (i1 = 0; i1 < k1; i1++)
scanf("%lf", &si[i1]);
}
// 分母
scanf("%*s %d %*s %d %*s", &mb, &n);
// 因数分解されていない
if (mb == 0) {
k1 = n + 1;
bo = new double [k1];
for (i1 = 0; i1 < k1; i1++)
scanf("%lf", &bo[i1]);
}
// 因数分解されている
else {
k1 = (n == 0) ? 1 : 2 * n;
bo = new double [k1];
for (i1 = 0; i1 < k1; i1++)
scanf("%lf", &bo[i1]);
}
/*
ゲインと位相の計算
*/
h = (log10(fu) - log10(fl)) / k;
ff = log10(fl);
for (i1 = 0; i1 <= k; i1++) {
// 周波数の対数を元に戻す
f = pow(10.0, ff);
// 値の計算
g = G_s(f, ms, m, si, mb, n, bo);
// ゲインと位相の計算
gain = 20.0 * log10(abs(g));
x = angle(g) * uc;
while (fabs(phase-x) > 180.0) {
if (x-phase > 180.0)
x -= 360.0;
else {
if (x-phase < -180.0)
x += 360.0;
}
}
phase = x;
// 出力
fprintf(og, "%f %f\n", f, gain);
fprintf(op, "%f %f\n", f, phase);
// 次の周波数
ff += h;
}
return 0;
}
/********************************************************************/
/* 伝達関数のｓにｊωを代入した値の計算 */
/* ff : ω（周波数） */
/* ms : 分子の表現方法 */
/* =0 : 多項式（si[0]+si[1]*s+si[1]*s^2+･･･） */
/* =1 : 因数分解（(si[0]+si[1]*s)*(si[2]+si[3]*s)*･･･） */
/* m : 分子の次数 */
/* si : 分子多項式の係数 */
/* mb : 分母の表現方法 */
/* =0 : 多項式（bo[0]+bo[1]*s+bo[1]*s^2+･･･） */
/* =1 : 因数分解（(bo[0]+bo[1]*s)*(bo[2]+bo[3]*s)*･･･） */
/* n : 分母の次数 */
/* bo : 分母多項式の係数 */
/* return : 結果 */
/********************************************************************/
Complex G_s(double ff, int ms, int m, double *si, int mb, int n, double *bo)
{
Complex f, x, y;
// 周波数を複素数に変換
f = Complex (0.0, ff);
// 分子
x = value(f, ms, m, si);
// 分母
y = value(f, mb, n, bo);
return x / y;
}
/****************************************************************/
/* 多項式のｓにｊωを代入した値の計算 */
/* f : ｊω（周波数，複素数） */
/* ms : 多項式の表現方法 */
/* =0 : 多項式（a[0]+a[1]*s+a[1]*s^2+･･･） */
/* =1 : 因数分解（(a[0]+a[1]*s)*(a[2]+a[3]*s)*･･･） */
/* n : 多項式の次数 */
/* a : 多項式の係数 */
/* return : 結果 */
/****************************************************************/
Complex value(Complex f, int ms, int n, double *a)
{
int i1, k1;
Complex u0(0.0, 0.0), u1(1.0, 0.0), x;
// 因数分解されていない
if (ms == 0) {
x = u0;
for (i1 = 0; i1 <= n; i1++)
x = x + a[i1] * pow(f, i1);
}
// 因数分解されている
else {
if (n == 0)
x = a[0] * u1;
else {
x = u1;
k1 = 2 * n;
for (i1 = 0; i1 < k1; i1 += 2)
x = x * (a[i1] + a[i1+1] * f);
}
}
return x;
}
/*
------------------------data1.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 0 次数 0 係数（次数が低い順） 1.0
分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0
------------------------data2.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 1 次数 0 係数（次数が低い順） 1.0
分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0
*/
/********************************/
/* 伝達関数のゲインと位相の計算 */
/* coded by Y.Suganuma */
/********************************/
import java.io.*;
import java.util.StringTokenizer;
public class Test {
/****************/
/* main program */
/****************/
public static void main(String args[]) throws IOException
{
double f, ff, fl, fu, h, gain, phase = 0.0, x, uc;
double bo[] = new double [1];
double si[] = new double [1];
int i1, k, k1, ms, mb, m, n;
Complex g;
String str;
StringTokenizer token;
BufferedReader in = new BufferedReader(new InputStreamReader(System.in));
uc = 90.0 / Math.asin(1.0);
/*
データの入力
*/
// 出力ファイル名とファイルのオープン
str = in.readLine();
token = new StringTokenizer(str, " ");
token.nextToken();
PrintStream g_o = new PrintStream(new FileOutputStream(token.nextToken()));
token.nextToken();
PrintStream p_o = new PrintStream(new FileOutputStream(token.nextToken()));
// 伝達関数データ
str = in.readLine();
token = new StringTokenizer(str, " ");
token.nextToken();
fl = Double.parseDouble(token.nextToken());
fu = Double.parseDouble(token.nextToken());
token.nextToken();
k = Integer.parseInt(token.nextToken());
// 分子
str = in.readLine();
token = new StringTokenizer(str, " ");
token.nextToken();
ms = Integer.parseInt(token.nextToken());
token.nextToken();
m = Integer.parseInt(token.nextToken());
token.nextToken();
// 因数分解されていない
if (ms == 0) {
k1 = m + 1;
si = new double [k1];
for (i1 = 0; i1 < k1; i1++)
si[i1] = Double.parseDouble(token.nextToken());
}
// 因数分解されている
else {
k1 = (m == 0) ? 1 : 2 * m;
si = new double [k1];
for (i1 = 0; i1 < k1; i1++)
si[i1] = Double.parseDouble(token.nextToken());
}
// 分母
str = in.readLine();
token = new StringTokenizer(str, " ");
token.nextToken();
mb = Integer.parseInt(token.nextToken());
token.nextToken();
n = Integer.parseInt(token.nextToken());
token.nextToken();
// 因数分解されていない
if (mb == 0) {
k1 = n + 1;
bo = new double [k1];
for (i1 = 0; i1 < k1; i1++)
bo[i1] = Double.parseDouble(token.nextToken());
}
// 因数分解されている
else {
k1 = (n == 0) ? 1 : 2 * n;
bo = new double [k1];
for (i1 = 0; i1 < k1; i1++)
bo[i1] = Double.parseDouble(token.nextToken());
}
/*
ゲインと位相の計算
*/
h = (log10(fu) - log10(fl)) / k;
ff = log10(fl);
for (i1 = 0; i1 <= k; i1++) {
// 周波数の対数を元に戻す
f = Math.pow(10.0, ff);
// 値の計算
g = G_s(f, ms, m, si, mb, n, bo);
// ゲインと位相の計算
gain = 20.0 * log10(Complex.abs(g));
x = Complex.angle(g) * uc;
while (Math.abs(phase-x) > 180.0) {
if (x-phase > 180.0)
x -= 360.0;
else {
if (x-phase < -180.0)
x += 360.0;
}
}
phase = x;
// 出力
g_o.println(f + " " + gain);
p_o.println(f + " " + phase);
// 次の周波数
ff += h;
}
g_o.close();
p_o.close();
}
/************/
/* log10(x) */
/************/
static double log10(double x)
{
return Math.log(x) / Math.log(10.0);
}
/********************************************************************/
/* 伝達関数のｓにｊωを代入した値の計算 */
/* ff : ω（周波数） */
/* ms : 分子の表現方法 */
/* =0 : 多項式（si[0]+si[1]*s+si[1]*s^2+･･･） */
/* =1 : 因数分解（(si[0]+si[1]*s)*(si[2]+si[3]*s)*･･･） */
/* m : 分子の次数 */
/* si : 分子多項式の係数 */
/* mb : 分母の表現方法 */
/* =0 : 多項式（bo[0]+bo[1]*s+bo[1]*s^2+･･･） */
/* =1 : 因数分解（(bo[0]+bo[1]*s)*(bo[2]+bo[3]*s)*･･･） */
/* n : 分母の次数 */
/* bo : 分母多項式の係数 */
/* return : 結果 */
/********************************************************************/
static Complex G_s(double ff, int ms, int m, double si[], int mb, int n, double bo[])
{
Complex f, x, y;
// 周波数を複素数に変換
f = new Complex (0.0, ff);
// 分子
x = value(f, ms, m, si);
// 分母
y = value(f, mb, n, bo);
return Complex.dev(x, y);
}
/****************************************************************/
/* 多項式のｓにｊωを代入した値の計算 */
/* f : ｊω（周波数，複素数） */
/* ms : 多項式の表現方法 */
/* =0 : 多項式（a[0]+a[1]*s+a[1]*s^2+･･･） */
/* =1 : 因数分解（(a[0]+a[1]*s)*(a[2]+a[3]*s)*･･･） */
/* n : 多項式の次数 */
/* a : 多項式の係数 */
/* return : 結果 */
/****************************************************************/
static Complex value(Complex f, int ms, int n, double a[])
{
int i1, k1;
Complex x;
// 因数分解されていない
if (ms == 0) {
x = new Complex (0.0, 0.0);
for (i1 = 0; i1 <= n; i1++)
x = Complex.add(x, Complex.mul(new Complex(a[i1]), Complex.pow(f, i1)));
}
// 因数分解されている
else {
if (n == 0)
x = new Complex (a[0]);
else {
x = new Complex (1.0, 0.0);
k1 = 2 * n;
for (i1 = 0; i1 < k1; i1 += 2)
x = Complex.mul(x, Complex.add(new Complex(a[i1]), Complex.mul(new Complex(a[i1+1]), f)));
}
}
return x;
}
}
/****************************/
/* クラスComplex */
/* coded by Y.Suganuma */
/****************************/
class Complex {
private double r; // 実部
private double i; // 虚部
/******************/
/* コンストラクタ */
/* a : 実部 */
/* b : 虚部 */
/******************/
Complex(double a, double b)
{
r = a;
i = b;
}
/******************/
/* コンストラクタ */
/* a : 実部 */
/******************/
Complex(double a)
{
r = a;
i = 0.0;
}
/******************/
/* コンストラクタ */
/******************/
Complex()
{
r = 0.0;
i = 0.0;
}
/********/
/* 出力 */
/********/
void output(PrintStream out, Complex a)
{
out.print("(" + a.r + ", " + a.i + ")");
}
/******************/
/* 複素数の絶対値 */
/******************/
static double abs(Complex a)
{
double x;
x = Math.sqrt(a.r * a.r + a.i * a.i);
return x;
}
/****************/
/* 複素数の角度 */
/****************/
static double angle(Complex a)
{
double x;
if (a.r == 0.0 && a.i == 0.0)
x = 0.0;
else
x = Math.atan2(a.i, a.r);
return x;
}
/****************/
/* 複素数のｎ乗 */
/****************/
static Complex pow(Complex a, int n)
{
int i1;
Complex c = new Complex (1);
if (n >= 0) {
for (i1 = 0; i1 < n; i1++)
c = Complex.mul(c, a);
}
else {
for (i1 = 0; i1 < -n; i1++)
c = Complex.mul(c, a);
c = Complex.dev(new Complex(1.0), c);
}
return c;
}
/****************/
/* 複素数の加算 */
/****************/
static Complex add(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r + b.r;
c.i = a.i + b.i;
return c;
}
/****************/
/* 複素数の減算 */
/****************/
static Complex sub(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r - b.r;
c.i = a.i - b.i;
return c;
}
/****************/
/* 複素数の乗算 */
/****************/
static Complex mul(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r * b.r - a.i * b.i;
c.i = a.i * b.r + a.r * b.i;
return c;
}
/****************/
/* 複素数の除算 */
/****************/
static Complex dev(Complex a, Complex b)
{
double x;
Complex c = new Complex();
x = 1.0 / (b.r * b.r + b.i * b.i);
c.r = (a.r * b.r + a.i * b.i) * x;
c.i = (a.i * b.r - a.r * b.i) * x;
return c;
}
}
/*
------------------------data1----------------------------
ゲイン gain 位相 phase
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 0 次数 0 係数（次数が低い順） 1.0
分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0
------------------------data2----------------------------
ゲイン gain 位相 phase
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 1 次数 0 係数（次数が低い順） 1.0
分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0
*/
/*****************************************/
/* 伝達関数のゲインと位相の計算（Swing） */
/* coded by Y.Suganuma */
/*****************************************/
import java.io.*;
import java.awt.*;
import java.awt.event.*;
import javax.swing.*;
import javax.swing.event.*;
import java.text.*;
public class BodeApplet_s extends JApplet implements ActionListener, DocumentListener
{
double fl, fu, si[][], bo[][];
int k, m[], n[], n_g = -1;
JTextField fl_t, fu_t, k_t, n_g_t, ex[];
JTextArea ta;
JButton bt;
Coef bunsi[], bunbo[];
Container cp;
JScrollPane sp;
Font f = new Font("TimesRoman", Font.BOLD, 20);
/************/
/* 初期設定 */
/************/
public void init()
{
// レイアウト，背景色，フォント
cp = getContentPane();
cp.setLayout(new BorderLayout(5, 5));
cp.setBackground(new Color(225, 255, 225));
// 上のパネル
JPanel pn1 = new JPanel();
pn1.setBackground(new Color(225, 255, 225));
pn1.setLayout(new GridLayout(2, 1, 5, 5));
cp.add(pn1, BorderLayout.NORTH);
JPanel pn11 = new JPanel();
pn11.setBackground(new Color(225, 255, 225));
pn1.add(pn11);
JLabel lb1 = new JLabel("式の数（グラフの数）");
lb1.setFont(f);
pn11.add(lb1);
n_g_t = new JTextField(3);
n_g_t.setFont(f);
n_g_t.getDocument().putProperty("name", "no");
n_g_t.getDocument().addDocumentListener(this);
pn11.add(n_g_t);
pn11.add(new JLabel(" "));
bt = new JButton("実行");
bt.setBackground(Color.pink);
bt.setFont(f);
bt.addActionListener(this);
pn11.add(bt);
JPanel pn12 = new JPanel();
pn12.setBackground(new Color(225, 255, 225));
pn1.add(pn12);
JLabel lb2 = new JLabel("周波数下限");
lb2.setFont(f);
pn12.add(lb2);
fl = 0.01;
fl_t = new JTextField("0.01", 3);
fl_t.setFont(f);
fl_t.getDocument().putProperty("name", "lower");
fl_t.getDocument().addDocumentListener(this);
pn12.add(fl_t);
JLabel lb3 = new JLabel(" 周波数上限");
lb3.setFont(f);
pn12.add(lb3);
fu = 100.0;
fu_t = new JTextField("100", 3);
fu_t.setFont(f);
fu_t.getDocument().putProperty("name", "upper");
fu_t.getDocument().addDocumentListener(this);
pn12.add(fu_t);
JLabel lb4 = new JLabel(" データ数");
lb4.setFont(f);
pn12.add(lb4);
k = 100;
k_t = new JTextField("100", 3);
k_t.setFont(f);
k_t.getDocument().putProperty("name", "data");
k_t.getDocument().addDocumentListener(this);
pn12.add(k_t);
// 中央のパネル
sp = new JScrollPane();
cp.add(sp, BorderLayout.CENTER);
// 下のパネル
JPanel pn6 = new JPanel();
pn6.setBackground(new Color(225, 255, 225));
cp.add(pn6, BorderLayout.SOUTH);
ta = new JTextArea(5, 35);
ta.setFont(f);
JScrollPane scroll2 = new JScrollPane(ta);
pn6.add(scroll2);
}
/************/
/* log10(x) */
/************/
static double log10(double x)
{
return Math.log(x) / Math.log(10.0);
}
/****************************************/
/* 伝達関数のｓにｊωを代入した値の計算 */
/* ff : ω（周波数） */
/* m : 分子の次数 */
/* si : 分子多項式の係数 */
/* n : 分母の次数 */
/* bo : 分母多項式の係数 */
/* return : 結果 */
/****************************************/
static Complex G_s(double ff, int m, double si[], int n, double bo[])
{
Complex f, x, y;
// 周波数を複素数に変換
f = new Complex (0.0, ff);
// 分子
x = value(f, m, si);
// 分母
y = value(f, n, bo);
return Complex.dev(x, y);
}
/**************************************/
/* 多項式のｓにｊωを代入した値の計算 */
/* f : ｊω（周波数，複素数） */
/* n : 多項式の次数 */
/* a : 多項式の係数 */
/* return : 結果 */
/**************************************/
static Complex value(Complex f, int n, double a[])
{
int i1, k1;
Complex x;
x = new Complex (0.0, 0.0);
for (i1 = 0; i1 <= n; i1++)
x = Complex.add(x, Complex.mul(new Complex(a[i1]), Complex.pow(f, i1)));
return x;
}
/******************************/
/* 上，左，下，右の余白の設定 */
/******************************/
public Insets getInsets()
{
return new Insets(10, 10, 10, 10);
}
/******************************/
/* ボタンが押されたときの処理 */
/******************************/
public void actionPerformed(ActionEvent e)
{
double x, h, f, ff, g_min, g_max, p_min, p_max, uc = 90.0 / Math.asin(1.0);
int i1, i2, k1, p, sw = 0;
String g_title[] = null;
Complex g;
ta.setForeground(Color.black);
ta.setText("");
// エラーの判定
if (n_g < 0) {
sw = 1;
ta.setForeground(Color.red);
ta.append(" 式の数を入力してください\n");
}
else {
g_title = new String [n_g];
for (i1 = 0; i1 < n_g; i1++) {
if (ex[i1].getText().length() <= 0) {
sw = 1;
ta.setForeground(Color.red);
ta.append(" " + (i1+1) + " 番目の式の説明(凡例)を入力してください\n");
}
else
g_title[i1] = ex[i1].getText();
}
for (i1 = 0; i1 < n_g; i1++) {
if (m[i1] < 0) {
sw = 1;
ta.setForeground(Color.red);
ta.append(" " + (i1+1) + " 番目の式の分子の次数を入力してください\n");
}
else {
for (i2 = 0; i2 <= m[i1]; i2++) {
if (bunsi[i1].bun_t[i2].getText().length() <= 0) {
sw = 1;
ta.setForeground(Color.red);
ta.append(" " + (i1+1) + " 番目の式の分子の s の " + i2 + " 次の係数を入力してください\n");
}
else
si[i1][i2] = Double.parseDouble(bunsi[i1].bun_t[i2].getText());
}
}
}
for (i1 = 0; i1 < n_g; i1++) {
if (n[i1] < 0) {
sw = 1;
ta.setForeground(Color.red);
ta.append(" " + (i1+1) + " 番目の式の分母の次数を入力してください\n");
}
else {
for (i2 = 0; i2 <= n[i1]; i2++) {
if (bunbo[i1].bun_t[i2].getText().length() <= 0) {
sw = 1;
ta.setForeground(Color.red);
ta.append(" " + (i1+1) + " 番目の式の分母の s の " + i2 + " 次の係数を入力してください\n");
}
else
bo[i1][i2] = Double.parseDouble(bunbo[i1].bun_t[i2].getText());
}
}
}
}
// ゲインと位相の計算
if (sw == 0) {
// 下限と上限の再設定
if (fl < 1.0) {
x = 0.1;
while (fl < x-1.0e-10)
x /= 10.0;
fl = x;
}
else {
x = 1.0;
while (fl > x-1.0e-10)
x *= 10.0;
fl = x;
}
if (fu < 1.0) {
x = 0.1;
while (fu < x+1.0e-10)
x /= 10.0;
fu = 10.0 * x;
}
else {
x = 1.0;
while (fu > x+1.0e-10)
x *= 10.0;
fu = x;
}
// 初期設定
double freq1[][] = new double [n_g][k+1];
double freq2[][] = new double [n_g][k+1];
double gain[][] = new double [n_g][k+1];
double phase[][] = new double [n_g][k+1];
h = (log10(fu) - log10(fl)) / k;
g_min = 0.0;
g_max = 0.0;
p_min = 0.0;
p_max = 0.0;
ta.setText(" 角周波数 ゲイン(dB) 位相(度)\n");
for (i1 = 0; i1 < n_g; i1++) {
ff = log10(fl);
ta.append(g_title[i1] + "\n");
for (i2 = 0; i2 <= k; i2++) {
// 周波数の対数を元に戻す
f = Math.pow(10.0, ff);
freq1[i1][i2] = f;
freq2[i1][i2] = f;
ta.append(" " + f);
// 値の計算
g = G_s(f, m[i1], si[i1], n[i1], bo[i1]);
// ゲインと位相の計算
gain[i1][i2] = 20.0 * log10(Complex.abs(g));
ta.append(" " + gain[i1][i2]);
if (i1 == 0 && i2 == 0) {
g_min = gain[i1][i2];
g_max = gain[i1][i2];
}
else {
if (gain[i1][i2] > g_max)
g_max = gain[i1][i2];
else if (gain[i1][i2] < g_min)
g_min = gain[i1][i2];
}
x = Complex.angle(g) * uc;
if (i2 > 0) {
while (Math.abs(phase[i1][i2-1]-x) > 180.0) {
if (x-phase[i1][i2-1] > 180.0)
x -= 360.0;
else {
if (x-phase[i1][i2-1] < -180.0)
x += 360.0;
}
}
}
phase[i1][i2] = x;
ta.append(" " + x + "\n");
if (i1 == 0 && i2 == 0) {
p_min = phase[i1][i2];
p_max = phase[i1][i2];
}
else {
if (phase[i1][i2] > p_max)
p_max = phase[i1][i2];
else if (phase[i1][i2] < p_min)
p_min = phase[i1][i2];
}
// 次の周波数
ff += h;
}
}
// グラフの描画
// グラフ，ｘ軸，及び，ｙ軸のタイトル
String title1[] = new String [3];
title1[0] = "ゲイン線図";
title1[1] = "角周波数";
title1[2] = "ゲイン（dB）";
// ゲイン線図
// ｘ軸目盛り
double x_scale1[] = new double[3];
x_scale1[0] = fl; // 最小値
x_scale1[1] = fu; // 最大値
x_scale1[2] = 1.0; // 最大値
// ｙ軸目盛り
double y_scale1[] = new double[3];
if ((g_max-g_min) < 1.0e-5) {
if (g_min > 0.0) {
k1 = (int)Math.round(g_min / 20);
y_scale1[0] = 20.0 * k1 - 20.0; // 最小値
y_scale1[1] = 20.0 * k1 + 20.0; // 最大値
}
else {
k1 = (int)Math.round(-g_min / 20);
y_scale1[0] = -20.0 * k1 - 20.0; // 最小値
y_scale1[1] = -20.0 * k1 + 20.0; // 最大値
}
}
else {
if (g_min > 0.0) {
k1 = (int)(g_min / 20);
y_scale1[0] = 20.0 * k1; // 最小値
}
else {
k1 = (int)(-g_min / 20);
if (Math.abs(-20.0*k1-g_min) > 1.0e-5)
k1++;
y_scale1[0] = -20.0 * k1; // 最小値
}
if (g_max > 0.0) {
k1 = (int)(g_max / 20);
if (Math.abs(20.0*k1-g_max) > 1.0e-5)
k1++;
y_scale1[1] = 20.0 * k1; // 最大値
}
else {
k1 = (int)(-g_max / 20);
y_scale1[1] = -20.0 * k1; // 最大値
}
}
y_scale1[2] = 20.0; // 刻み幅
// ｘ軸の小数点以下桁数
p = 0;
if (fl < 1.0) {
p = 1;
x = 0.1;
while (fl < x-1.0e-10) {
x /= 10.0;
p++;
}
}
Bode gp1 = new Bode(title1, g_title, x_scale1, p, y_scale1, 0, freq1, gain, true, true);
// 位相線図
// グラフ，ｘ軸，及び，ｙ軸のタイトル
String title2[] = new String [3];
title2[0] = "位相線図";
title2[1] = "角周波数";
title2[2] = "位相（度）";
// ｘ軸目盛り
double x_scale2[] = new double[3];
x_scale2[0] = fl; // 最小値
x_scale2[1] = fu; // 最大値
x_scale2[2] = 1.0; // 最大値
// ｙ軸目盛り
double y_scale2[] = new double[3];
if ((p_max-p_min) < 1.0e-5) {
if (p_min > 0.0) {
k1 = (int)Math.round(p_min / 90);
y_scale2[0] = 90.0 * k1 - 90.0; // 最小値
y_scale2[1] = 90.0 * k1 + 90.0; // 最大値
}
else {
k1 = (int)Math.round(-p_min / 90);
y_scale2[0] = -90.0 * k1 - 90.0; // 最小値
y_scale2[1] = -90.0 * k1 + 90.0; // 最大値
}
}
else {
if (p_min > 0.0) {
k1 = (int)(p_min / 90);
y_scale2[0] = 90.0 * k1; // 最小値
}
else {
k1 = (int)(-p_min / 90);
if (Math.abs(-90.0*k1-p_min) > 1.0e-5)
k1++;
y_scale2[0] = -90.0 * k1; // 最小値
}
if (p_max > 0.0) {
k1 = (int)(p_max / 90);
if (Math.abs(90.0*k1-p_max) > 1.0e-5)
k1++;
y_scale2[1] = 90.0 * k1; // 最大値
}
else {
k1 = (int)(-p_max / 90);
y_scale2[1] = -90.0 * k1; // 最大値
}
}
y_scale2[2] = 90.0; // 刻み幅
Bode gp2 = new Bode(title2, g_title, x_scale2, p, y_scale2, 0, freq2, phase, true, true);
}
}
/************************************/
/* パラメータが変更されたときの処理 */
/************************************/
public void changedUpdate(DocumentEvent e) {}
public void removeUpdate(DocumentEvent e) {}
public void insertUpdate(DocumentEvent e)
{
int i1;
String key = e.getDocument().getProperty("name").toString();
try {
if (key.equals("no")) {
n_g = Integer.parseInt(n_g_t.getText());
ta.setForeground(Color.black);
ta.setText("");
if (n_g <= 0) {
ta.setForeground(Color.red);
ta.setText(" 0 より大きい数値を入力してください");
}
else {
remove(sp);
si = new double [n_g][];
bo = new double [n_g][];
m = new int [n_g];
n = new int [n_g];
for (i1 = 0; i1 < n_g; i1++) {
m[i1] = -1;
n[i1] = -1;
}
JPanel pn2 = new JPanel();
pn2.setBackground(new Color(225, 255, 225));
pn2.setLayout(new GridLayout(n_g, 1, 5, 5));
sp = new JScrollPane(pn2);
JScrollPane sp1[] = new JScrollPane [n_g];
JScrollPane sp2[] = new JScrollPane [n_g];
JPanel pn3[] = new JPanel [n_g];
JPanel pn4[] = new JPanel [n_g];
JPanel pn5[] = new JPanel [n_g];
JLabel lb[] = new JLabel [n_g];
ex = new JTextField [n_g];
bunsi = new Coef [n_g];
bunbo = new Coef [n_g];
for (i1 = 0; i1 < n_g; i1++) {
pn3[i1] = new JPanel();
pn3[i1].setBackground(new Color(255, 255, 225));
pn3[i1].setLayout(new BorderLayout(5, 5));
pn2.add(pn3[i1]);
pn4[i1] = new JPanel();
pn4[i1].setBackground(Color.white);
lb[i1] = new JLabel((i1+1) + "番目の式の説明：");
lb[i1].setFont(f);
pn4[i1].add(lb[i1]);
ex[i1] = new JTextField(10);
ex[i1].setFont(f);
pn4[i1].add(ex[i1]);
pn3[i1].add(pn4[i1], BorderLayout.NORTH);
pn5[i1] = new JPanel();
pn5[i1].setBackground(Color.white);
pn5[i1].setLayout(new GridLayout(1, 2, 5, 5));
bunsi[i1] = new Coef(0, this, i1);
sp1[i1] = new JScrollPane(bunsi[i1]);
pn5[i1].add(sp1[i1]);
bunbo[i1] = new Coef(1, this, i1);
sp2[i1] = new JScrollPane(bunbo[i1]);
pn5[i1].add(sp2[i1]);
pn3[i1].add(pn5[i1], BorderLayout.CENTER);
}
cp.add(sp, BorderLayout.CENTER);
getParent().validate();
}
}
else if (key.equals("lower")) {
fl = Double.parseDouble(fl_t.getText());
fu = Double.parseDouble(fu_t.getText());
ta.setForeground(Color.black);
ta.setText("");
if (fl <= 0.0) {
ta.setForeground(Color.red);
ta.setText(" 0 より大きい数値を入力してください");
}
else if (fu <= fl) {
ta.setForeground(Color.red);
ta.setText(" 上限より小さい数値を入力してください");
}
}
else if (key.equals("upper")) {
fl = Double.parseDouble(fl_t.getText());
fu = Double.parseDouble(fu_t.getText());
ta.setForeground(Color.black);
ta.setText("");
if (fu <= 0.0) {
ta.setForeground(Color.red);
ta.setText(" 0 より大きい数値を入力してください");
}
else if (fu <= fl) {
ta.setForeground(Color.red);
ta.setText(" 下限より大きい数値を入力してください");
}
}
else if (key.equals("data")) {
k = Integer.parseInt(k_t.getText());
ta.setForeground(Color.black);
ta.setText("");
if (k <= 0) {
ta.setForeground(Color.red);
ta.setText(" 0 より大きい数値を入力してください");
}
}
}
catch (NumberFormatException em) {}
}
}
/****************************/
/* ボード線図の描画 */
/* coded by Y.Suganuma */
/****************************/
class Bode extends JFrame {
Draw_bode pn;
/*********************************************************/
/* コンストラクタ */
/* title_i : グラフ，ｘ軸，及び，ｙ軸のタイトル */
/* g_title_i : 凡例 */
/* x_scale_i : データの最小値，最大値，目盛幅（ｙ） */
/* place_x_i : 小数点以下の桁数（ｘ軸） */
/* y_scale_i : データの最小値，最大値，目盛幅（ｙ） */
/* place_y_i : 小数点以下の桁数（ｙ軸） */
/* data_x_i : グラフのデータ（ｘ軸） */
/* data_y_i : グラフのデータ（ｙ軸） */
/* d_t_i : タイトル表示の有無 */
/* d_g_i : 凡例表示の有無 */
/*********************************************************/
Bode(String title_i[], String g_title_i[], double x_scale_i[],
int place_x_i, double y_scale_i[], int place_y_i,
double data_x_i[][], double data_y_i[][], boolean d_t_i,
boolean d_g_i)
{
// JFrameクラスのコンストラクタの呼び出し
super("ボード線図");
// Windowサイズと表示位置を設定
int width = 900, height = 600; // Windowの大きさ（初期サイズ）
setSize(width, height);
Toolkit tool = getToolkit();
Dimension d = tool.getScreenSize();
setLocation(d.width / 2 - width / 2, d.height / 2 - height / 2);
// 描画パネル
Container cp = getContentPane();
pn = new Draw_bode(title_i, g_title_i, x_scale_i, place_x_i, y_scale_i, place_y_i, data_x_i, data_y_i, d_t_i, d_g_i, this);
cp.add(pn);
// ウィンドウを表示
setVisible(true);
// イベントアダプタ
addWindowListener(new WinEnd());
addComponentListener(new ComponentResize());
}
/**********************/
/* Windowのサイズ変化 */
/**********************/
class ComponentResize extends ComponentAdapter
{
public void componentResized(ComponentEvent e)
{
pn.repaint();
}
}
/************/
/* 終了処理 */
/************/
class WinEnd extends WindowAdapter
{
public void windowClosing(WindowEvent e) {
setVisible(false);
}
}
}
class Draw_bode extends JPanel {
String title[]; // グラフのタイトル
String g_title[]; // 凡例（グラフの内容）
double xx_scale[]; // ｙ軸目盛り
double x_scale[]; // 元のｙ軸目盛り
double y_scale[]; // ｙ軸目盛り
double data_x[][]; // 元のデータ
double data_xx[][], data_y[][]; // データ
boolean d_t; // タイトル表示の有無
boolean d_g; // 凡例表示の有無
boolean log_c = false; // 対数に変換したか否か
int place_x; // 小数点以下の桁数（ｘ軸）
int place_y; // 小数点以下の桁数（ｙ軸）
int width = 900, height = 600; // Windowの大きさ（初期サイズ）
int bx1, bx2, by1, by2; // 表示切り替えボタンの位置
Bode bd;
String change = " 色 "; // 表示切り替えボタン
float line_w = 1.0f; // 折れ線グラフ等の線の太さ
Color cl[] = {Color.black, Color.magenta, Color.blue, Color.orange, Color.cyan,
Color.pink, Color.green, Color.yellow, Color.darkGray, Color.red}; // グラフの色
int n_g; // グラフの数
/******************/
/* コンストラクタ */
/******************/
Draw_bode(String title_i[], String g_title_i[], double x_scale_i[],
int place_x_i, double y_scale_i[], int place_y_i,
double data_x_i[][], double data_y_i[][], boolean d_t_i,
boolean d_g_i, Bode bd1) {
// 背景色
setBackground(Color.white);
// テーブルデータの保存
title = title_i;
g_title = g_title_i;
x_scale = x_scale_i;
place_x = place_x_i;
y_scale = y_scale_i;
place_y = place_y_i;
data_x = data_x_i;
data_y = data_y_i;
d_t = d_t_i;
d_g = d_g_i;
bd = bd1;
int i1, i2;
int n_g = g_title.length;
int n_p = data_x[0].length;
xx_scale = new double [3];
data_xx = new double [n_g][n_p];
xx_scale[0] = x_scale[0];
xx_scale[1] = x_scale[1];
for (i1 = 0; i1 < n_g; i1++) {
for (i2 = 0; i2 < n_p; i2++)
data_xx[i1][i2] = data_x[i1][i2];
}
// イベントアダプタ
addMouseListener(new ClickMouse(this));
}
/********/
/* 描画 */
/********/
public void paintComponent (Graphics g)
{
super.paintComponent(g); // 親クラスの描画(必ず必要)
double r, x1, y1, y2, sp, x_scale_org = 0.0;
int i1, i2, k, k_x, k_y, k1, k2, kx, kx1, ky, ky1, han, len;
int x_l, x_r, y_u, y_d; // 描画領域
int f_size; // フォントサイズ
int n_p; // データの数
String s1;
Font f;
FontMetrics fm;
Graphics2D g2 = (Graphics2D)g;
//
// Windowサイズの取得
//
Insets insets = bd.getInsets();
Dimension d = bd.getSize();
width = d.width - (insets.left + insets.right);
height = d.height - (insets.top + insets.bottom);
x_l = insets.left + 10;
x_r = d.width - insets.right - 10;
y_u = 20;
y_d = d.height - insets.bottom - insets.top;
//
// グラフタイトルの表示
//
r = 0.05; // タイトル領域の割合
f_size = ((y_d - y_u) < (x_r - x_l)) ? (int)((y_d - y_u) * r) : (int)((x_r - x_l) * r);
if (f_size < 5)
f_size = 5;
if (d_t) {
f = new Font("TimesRoman", Font.BOLD, f_size);
g.setFont(f);
fm = g.getFontMetrics(f);
len = fm.stringWidth(title[0]);
g.drawString(title[0], (x_l+x_r)/2-len/2, y_d-f_size/2);
y_d -= f_size;
}
//
// 表示切り替えボタンの設置
//
f_size = (int)(0.8 * f_size);
if (f_size < 5)
f_size = 5;
f = new Font("TimesRoman", Font.PLAIN, f_size);
fm = g.getFontMetrics(f);
g.setFont(f);
g.setColor(Color.yellow);
len = fm.stringWidth(change);
bx1 = x_r - len - 7 * f_size / 10;
by1 = y_u - f_size / 2;
bx2 = bx1 + len + f_size / 2;
by2 = by1 + 6 * f_size / 5;
g.fill3DRect(bx1, by1, len+f_size/2, 6*f_size/5, true);
g.setColor(Color.black);
g.drawString(change, x_r-len-f_size/2, y_u+f_size/2);
//
// 凡例の表示
//
n_g = g_title.length;
if (d_g) {
han = 0;
for (i1 = 0; i1 < n_g; i1++) {
len = fm.stringWidth(g_title[i1]);
if (len > han)
han = len;
}
han += 15;
r = 0.2; // 凡例領域の割合
k1 = (int)((x_r - x_l) * r);
if (han > k1)
han = k1;
kx = x_r - han;
ky = y_u + 3 * f_size / 2;
k = 0;
g2.setStroke(new BasicStroke(7.0f));
for (i1 = 0; i1 < n_g; i1++) {
g.setColor(cl[k]);
g.drawLine(kx, ky, kx+10, ky);
g.setColor(Color.black);
g.drawString(g_title[i1], kx+15, ky+2*f_size/5);
k++;
if (k >= cl.length)
k = 0;
ky += f_size;
}
g2.setStroke(new BasicStroke(1.0f));
x_r -= (han + 10);
}
else
x_r -= (int)(0.03 * (x_r - x_l));
//
// ｘ軸の対数
//
n_p = data_x[0].length;
x_scale_org = x_scale[0];
xx_scale[0] = Math.log(x_scale[0]) / Math.log(10.0);
xx_scale[1] = Math.log(x_scale[1]) / Math.log(10.0);
xx_scale[2] = 1.0;
for (i1 = 0; i1 < n_g; i1++) {
for (i2 = 0; i2 < n_p; i2++)
data_xx[i1][i2] = Math.log(data_x[i1][i2]) / Math.log(10.0);
}
//
// ｘ軸及びｙ軸のタイトルの表示
//
if (title[1].length() > 0 && !title[1].equals("-")) {
len = fm.stringWidth(title[1]);
g.drawString(title[1], (x_l+x_r)/2-len/2, y_d-4*f_size/5);
y_d -= 7 * f_size / 4;
}
else
y_d -= f_size / 2;
if (title[2].length() > 0 && !title[2].equals("-")) {
g.drawString(title[2], x_l, y_u+f_size/2);
y_u += f_size;
}
//
// ｘ軸，ｙ軸，及び，各軸の目盛り
//
f_size = (int)(0.8 * f_size);
if (f_size < 5)
f_size = 5;
f = new Font("TimesRoman", Font.PLAIN, f_size);
fm = g.getFontMetrics(f);
y_d -= 3 * f_size / 2;
k_y = (int)((y_scale[1] - y_scale[0]) / (0.99 * y_scale[2]));
k_x = (int)((xx_scale[1] - xx_scale[0]) / (0.99 * xx_scale[2]));
g.setFont(f);
DecimalFormat df_x, df_y;
df_x = new DecimalFormat("#");
df_y = new DecimalFormat("#");
if (place_x != 0) {
s1 = "0.";
for (i1 = 0; i1 < place_x; i1++)
s1 += "0";
df_x = new DecimalFormat(s1);
}
if (place_y != 0) {
s1 = "#.";
for (i1 = 0; i1 < place_y; i1++)
s1 += "0";
df_y = new DecimalFormat(s1);
}
// ｙ軸
y1 = y_scale[0];
len = 0;
for (i1 = 0; i1 < k_y+1; i1++) {
s1 = df_y.format(y1);
k1 = fm.stringWidth(s1);
if (k1 > len)
len = k1;
y1 += y_scale[2];
}
g.drawLine(x_l+len+5, y_u, x_l+len+5, y_d);
g.drawLine(x_r, y_u, x_r, y_d);
y1 = y_scale[0];
x1 = y_d;
sp = (double)(y_d - y_u) / k_y;
for (i1 = 0; i1 < k_y+1; i1++) {
ky = (int)Math.round(x1);
s1 = df_y.format(y1);
k1 = fm.stringWidth(s1);
g.drawString(s1, x_l+len-k1, ky+f_size/2);
g.drawLine(x_l+len+5, ky, x_r, ky);
y1 += y_scale[2];
x1 -= sp;
}
x_l += (len + 5);
// ｘ軸
x1 = x_scale_org;
y1 = x_l;
sp = (double)(x_r - x_l) / k_x;
for (i1 = 0; i1 < k_x+1; i1++) {
kx = (int)Math.round(y1);
s1 = df_x.format(x1);
k1 = fm.stringWidth(s1);
g.drawString(s1, kx-k1/2, y_d+6*f_size/5);
g.drawLine(kx, y_d, kx, y_u);
if (i1 != k_x) {
g.setColor(Color.darkGray);
for (i2 = 2; i2 <= 9; i2++) {
y2 = Math.log(x1 * i2) / Math.log(10.0);
kx = x_l + (int)Math.round(((x_r - x_l) * (y2 - xx_scale[0]) / (xx_scale[1] - xx_scale[0])));
g.drawLine(kx, y_d, kx, y_u);
}
g.setColor(Color.black);
}
x1 *= 10.0;
y1 += sp;
}
//
// グラフの表示
//
g2.setStroke(new BasicStroke(line_w));
k1 = 0;
for (i1 = 0; i1 < n_g; i1++) {
g.setColor(cl[k1]);
kx1 = 0;
ky1 = 0;
for (i2 = 0; i2 < n_p; i2++) {
kx = x_l + (int)((x_r - x_l) * (data_xx[i1][i2] - xx_scale[0]) / (xx_scale[1] - xx_scale[0]));
ky = y_d - (int)((y_d - y_u) * (data_y[i1][i2] - y_scale[0]) / (y_scale[1] - y_scale[0]));
if (i2 > 0)
g.drawLine(kx1, ky1, kx, ky);
kx1 = kx;
ky1 = ky;
}
k1++;
if (k1 >= cl.length)
k1 = 0;
}
g2.setStroke(new BasicStroke(1.0f));
}
/************************************/
/* マウスがクリックされたときの処理 */
/************************************/
class ClickMouse extends MouseAdapter
{
Draw_bode dd;
ClickMouse(Draw_bode dd1)
{
dd = dd1;
}
public void mouseClicked(MouseEvent e)
{
int xp = e.getX();
int yp = e.getY();
// グラフの色，線の太さ等
if (xp > bx1 && xp < bx2 && yp > by1 && yp < by2) {
Modify md = new Modify(dd.bd, dd);
md.setVisible(true);
}
}
}
}
/****************************/
/* クラスCoef */
/* coded by Y.Suganuma */
/****************************/
class Coef extends JPanel implements DocumentListener
{
int type, no;
BodeApplet_s ba;
JTextField n_t;
JTextField bun_t[];
JLabel lb[];
JPanel pn;
Font f = new Font("TimesRoman", Font.BOLD, 20);
/****************************/
/* コンストラクタ */
/* type_i : =0 : 分子 */
/* =1 : 分母 */
/* ba_i : BodeApplet_s */
/* no_i : 式番号 */
/****************************/
Coef(int type_i, BodeApplet_s ba_i, int no_i)
{
int i1;
type = type_i;
ba = ba_i;
no = no_i;
setBackground(Color.white);
// 次数
setLayout(new BorderLayout(5, 5));
JPanel pn1 = new JPanel();
pn1.setBackground(Color.white);
add(pn1, BorderLayout.NORTH);
String str;
if (type == 0)
str = (no + 1) + ".分子の次数";
else
str = (no + 1) + ".分母の次数";
JLabel jisu = new JLabel(str);
jisu.setFont(f);
pn1.add(jisu);
n_t = new JTextField(2);
n_t.setFont(f);
n_t.getDocument().addDocumentListener(this);
pn1.add(n_t);
// 係数
pn = new JPanel();
add(pn, BorderLayout.CENTER);
}
/************************************/
/* パラメータが変更されたときの処理 */
/************************************/
public void changedUpdate(DocumentEvent e) {}
public void removeUpdate(DocumentEvent e) {}
public void insertUpdate(DocumentEvent e)
{
int i1, n;
ba.ta.setForeground(Color.black);
ba.ta.setText("");
try {
n = Integer.parseInt(n_t.getText());
if (n < 0) {
n_t.setText("0");
n = 0;
}
if (type == 0) {
ba.m[no] = n;
ba.si[no] = new double [n+1];
}
else {
ba.n[no] = n;
ba.bo[no] = new double [n+1];
}
remove(pn);
pn = new JPanel();
pn.setBackground(Color.white);
pn.setLayout(new GridLayout(n+1, 1, 5, 0));
add(pn, BorderLayout.CENTER);
JPanel pnx[] = new JPanel [n+1];
bun_t = new JTextField [n+1];
lb = new JLabel [n+1];
pnx[0] = new JPanel();
pnx[0].setBackground(Color.white);
pn.add(pnx[0]);
lb[0] = new JLabel("定数項");
lb[0].setFont(f);
pnx[0].add(lb[0]);
bun_t[0] = new JTextField(3);
bun_t[0].setFont(f);
pnx[0].add(bun_t[0]);
for (i1 = 1; i1 <= n; i1++) {
pnx[i1] = new JPanel();
pnx[i1].setBackground(Color.white);
pn.add(pnx[i1]);
lb[i1] = new JLabel("s の "+i1+" 次の項");
lb[i1].setFont(f);
pnx[i1].add(lb[i1]);
bun_t[i1] = new JTextField(3);
bun_t[i1].setFont(f);
pnx[i1].add(bun_t[i1]);
}
getParent().validate();
}
catch (NumberFormatException em) {}
}
}
/****************************/
/* クラスComplex */
/* coded by Y.Suganuma */
/****************************/
class Complex {
private double r; // 実部
private double i; // 虚部
/******************/
/* コンストラクタ */
/* a : 実部 */
/* b : 虚部 */
/******************/
Complex(double a, double b)
{
r = a;
i = b;
}
/******************/
/* コンストラクタ */
/* a : 実部 */
/******************/
Complex(double a)
{
r = a;
i = 0.0;
}
/******************/
/* コンストラクタ */
/******************/
Complex()
{
r = 0.0;
i = 0.0;
}
/********/
/* 出力 */
/********/
void output(PrintStream out, Complex a)
{
out.print("(" + a.r + ", " + a.i + ")");
}
/******************/
/* 複素数の絶対値 */
/******************/
static double abs(Complex a)
{
double x;
x = Math.sqrt(a.r * a.r + a.i * a.i);
return x;
}
/****************/
/* 複素数の角度 */
/****************/
static double angle(Complex a)
{
double x;
if (a.r == 0.0 && a.i == 0.0)
x = 0.0;
else
x = Math.atan2(a.i, a.r);
return x;
}
/****************/
/* 複素数のｎ乗 */
/****************/
static Complex pow(Complex a, int n)
{
int i1;
Complex c = new Complex (1);
if (n >= 0) {
for (i1 = 0; i1 < n; i1++)
c = Complex.mul(c, a);
}
else {
for (i1 = 0; i1 < -n; i1++)
c = Complex.mul(c, a);
c = Complex.dev(new Complex(1.0), c);
}
return c;
}
/****************/
/* 複素数の加算 */
/****************/
static Complex add(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r + b.r;
c.i = a.i + b.i;
return c;
}
/****************/
/* 複素数の減算 */
/****************/
static Complex sub(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r - b.r;
c.i = a.i - b.i;
return c;
}
/****************/
/* 複素数の乗算 */
/****************/
static Complex mul(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r * b.r - a.i * b.i;
c.i = a.i * b.r + a.r * b.i;
return c;
}
/****************/
/* 複素数の除算 */
/****************/
static Complex dev(Complex a, Complex b)
{
double x;
Complex c = new Complex();
x = 1.0 / (b.r * b.r + b.i * b.i);
c.r = (a.r * b.r + a.i * b.i) * x;
c.i = (a.i * b.r - a.r * b.i) * x;
return c;
}
}
/****************************/
/* 色及び線の太さの変更 */
/* coded by Y.Suganuma */
/****************************/
class Modify extends JDialog implements ActionListener, TextListener {
Draw_bode dd; // ボード線図
JButton bt_dd;
TextField rgb[];
TextField r[];
TextField g[];
TextField b[];
JTextField tx;
JRadioButton r1, r2;
float line_w = 1.0f; // 折れ線グラフ等の線の太さ
Color cl[]; // グラフの色
int n_g; // グラフの数
int wd; // 線の太さを変更するか
int n;
JPanel jp[];
// ボード線図
Modify(Frame host, Draw_bode dd1)
{
super(host, "色と線の変更", true);
// 初期設定
dd = dd1;
wd = 1;
n_g = dd.n_g;
if (n_g > 10)
n_g = 10;
n = n_g + 2;
line_w = dd.line_w;
cl = new Color[n_g];
for (int i1 = 0; i1 < n_g; i1++)
cl[i1] = dd.cl[i1];
set();
// ボタン
Font f = new Font("TimesRoman", Font.BOLD, 20);
bt_dd = new JButton("OK");
bt_dd.setFont(f);
bt_dd.addActionListener(this);
jp[n-1].add(bt_dd);
}
// 設定
void set()
{
setSize(450, 60*(n));
Container cp = getContentPane();
cp.setBackground(Color.white);
cp.setLayout(new GridLayout(n, 1, 5, 5));
jp = new JPanel[n];
for (int i1 = 0; i1 < n; i1++) {
jp[i1] = new JPanel();
cp.add(jp[i1]);
}
Font f = new Font("TimesRoman", Font.BOLD, 20);
// 色の変更
JLabel lb[][] = new JLabel[n_g][3];
rgb = new TextField[n_g];
r = new TextField[n_g];
g = new TextField[n_g];
b = new TextField[n_g];
for (int i1 = 0; i1 < n_g; i1++) {
rgb[i1] = new TextField(3);
rgb[i1].setFont(f);
rgb[i1].setBackground(new Color(cl[i1].getRed(), cl[i1].getGreen(), cl[i1].getBlue()));
jp[i1].add(rgb[i1]);
lb[i1][0] = new JLabel(" 赤");
lb[i1][0].setFont(f);
jp[i1].add(lb[i1][0]);
r[i1] = new TextField(3);
r[i1].setFont(f);
r[i1].setBackground(Color.white);
r[i1].setText(Integer.toString(cl[i1].getRed()));
r[i1].addTextListener(this);
jp[i1].add(r[i1]);
lb[i1][1] = new JLabel("緑");
lb[i1][1].setFont(f);
jp[i1].add(lb[i1][1]);
g[i1] = new TextField(3);
g[i1].setFont(f);
g[i1].setBackground(Color.white);
g[i1].setText(Integer.toString(cl[i1].getGreen()));
g[i1].addTextListener(this);
jp[i1].add(g[i1]);
lb[i1][2] = new JLabel("青");
lb[i1][2].setFont(f);
jp[i1].add(lb[i1][2]);
b[i1] = new TextField(3);
b[i1].setFont(f);
b[i1].setBackground(Color.white);
b[i1].setText(Integer.toString(cl[i1].getBlue()));
b[i1].addTextListener(this);
jp[i1].add(b[i1]);
}
// 線の変更
if (wd > 0) {
JLabel lb1 = new JLabel("線の太さ：");
lb1.setFont(f);
jp[n_g].add(lb1);
tx = new JTextField(2);
tx.setFont(f);
tx.setBackground(Color.white);
tx.setText(Integer.toString((int)line_w));
jp[n_g].add(tx);
}
}
// TextFieldの内容が変更されたときの処理
public void textValueChanged(TextEvent e)
{
for (int i1 = 0; i1 < n_g; i1++) {
if (e.getSource() == r[i1] || e.getSource() == g[i1] || e.getSource() == b[i1]) {
String str = r[i1].getText();
int rc = str.length()>0 ? Integer.parseInt(str) : 0;
str = g[i1].getText();
int gc = str.length()>0 ? Integer.parseInt(str) : 0;
str = b[i1].getText();
int bc = str.length()>0 ? Integer.parseInt(str) : 0;
rgb[i1].setBackground(new Color(rc, gc, bc));
}
}
}
// 値の設定
public void actionPerformed(ActionEvent e)
{
for (int i1 = 0; i1 < n_g; i1++) {
String str = r[i1].getText();
int rc = str.length()>0 ? Integer.parseInt(str) : 0;
str = g[i1].getText();
int gc = str.length()>0 ? Integer.parseInt(str) : 0;
str = b[i1].getText();
int bc = str.length()>0 ? Integer.parseInt(str) : 0;
dd.cl[i1] = new Color(rc, gc, bc);
}
dd.line_w = Integer.parseInt(tx.getText());
dd.repaint();
setVisible(false);
}
}
<!DOCTYPE HTML>
<HTML>
<HEAD>
<TITLE>ボード線図</TITLE>
<META HTTP-EQUIV="Content-Type" CONTENT="text/html; charset=utf-8">
<SCRIPT TYPE="text/javascript">
max_g = 5;
max_order = 10;
n_g = 1;
function main()
{
// 入力データの有無チェック
n_g = parseInt(document.getElementById("n_g").value);
let err = 0;
for (let i1 = 0; i1 < n_g; i1++) {
let str1 = "s" + i1;
let str2 = document.getElementById(str1).value;
if (str2 == "") {
err = 1;
alert((i1+1) + " 番目の式の分子の次数を入力して下さい");
}
else {
let k1 = parseInt(str2);
for (let i2 = 0; i2 <= k1; i2++) {
str1 = "s" + i1 + i2;
str2 = document.getElementById(str1).value;
if (str2 == "") {
err = 1;
alert((i1+1) + " 番目の式の分子の" + i2 + "次の項を入力して下さい");
}
}
}
str1 = "b" + i1;
str2 = document.getElementById(str1).value;
if (str2 == "") {
err = 1;
alert((i1+1) + " 番目の式の分母の次数を入力して下さい");
}
else {
let k1 = parseInt(str2);
for (let i2 = 0; i2 <= k1; i2++) {
str1 = "b" + i1 + i2;
str2 = document.getElementById(str1).value;
if (str2 == "") {
err = 1;
alert((i1+1) + " 番目の式の分子の" + i2 + "次の項を入力して下さい");
}
}
}
}
if (err == 0) {
// ゲインと位相の計算
let fl = parseFloat(document.getElementById("low").value);
let fu = parseFloat(document.getElementById("up").value);
let k = parseInt(document.getElementById("n_data").value);
// 下限と上限の再設定
if (fl < 1.0) {
let x = 0.1;
while (fl < x-1.0e-10)
x /= 10.0;
fl = x;
}
else {
let x = 1.0;
while (fl > x-1.0e-10)
x *= 10.0;
fl = x;
}
if (fu < 1.0) {
let x = 0.1;
while (fu < x+1.0e-10)
x /= 10.0;
fu = 10.0 * x;
}
else {
let x = 1.0;
while (fu > x+1.0e-10)
x *= 10.0;
fu = x;
}
// 初期設定
let uc = 90.0 / Math.asin(1.0);
let g_title = new Array(); // グラフタイトル
let m = new Array(); // 分子の次数
let n = new Array(); // 分母の次数
let si = new Array(); // 分子の係数
let bo = new Array(); // 分母の係数
let freq1 = new Array(); // ゲイン線図の周波数
let freq2 = new Array(); // 位相線図の周波数
let gain = new Array(); // ゲイン
let phase = new Array(); // 位相
for (let i1 = 0; i1 < n_g; i1++) {
si[i1] = new Array();
bo[i1] = new Array();
freq1[i1] = new Array();
freq2[i1] = new Array();
gain[i1] = new Array();
phase[i1] = new Array();
}
let h = (log10(fu) - log10(fl)) / k;
let g_min = 0.0;
let g_max = 0.0;
let p_min = 0.0;
let p_max = 0.0;
let ta = " 角周波数 ゲイン(dB) 位相(度)\n";
for (let i1 = 0; i1 < n_g; i1++) {
// データの取得
g_title[i1] = document.getElementById("exp"+i1).value;
let str = "s" + i1;
m[i1] = parseInt(document.getElementById(str).value);
for (let i2 = 0; i2 <= m[i1]; i2++) {
str = "s" + i1 + i2;
si[i1][i2] = parseFloat(document.getElementById(str).value);
}
str = "b" + i1;
n[i1] = parseInt(document.getElementById(str).value);
for (let i2 = 0; i2 <= n[i1]; i2++) {
str = "b" + i1 + i2;
bo[i1][i2] = parseFloat(document.getElementById(str).value);
}
// 計算
ff = log10(fl);
ta = ta + g_title[i1] + "\n";
for (let i2 = 0; i2 <= k; i2++) {
// 周波数の対数を元に戻す
let f = Math.pow(10.0, ff);
freq1[i1][i2] = Math.round(100 * f) / 100;
freq2[i1][i2] = Math.round(100 * f) / 100;
ta = ta + " " + f;
// 値の計算
let g = G_s(f, m[i1], si[i1], n[i1], bo[i1]);
// ゲインの計算
gain[i1][i2] = 20.0 * log10(g.abs());
ta = ta + " " + gain[i1][i2];
if (i1 == 0 && i2 == 0) {
g_min = gain[i1][i2];
g_max = gain[i1][i2];
}
else {
if (gain[i1][i2] > g_max)
g_max = gain[i1][i2];
else if (gain[i1][i2] < g_min)
g_min = gain[i1][i2];
}
gain[i1][i2] = Math.round(100 * gain[i1][i2]) / 100;
// 位相の計算
let x = g.angle() * uc;
if (i2 > 0) {
while (Math.abs(phase[i1][i2-1]-x) > 180.0) {
if (x-phase[i1][i2-1] > 180.0)
x -= 360.0;
else {
if (x-phase[i1][i2-1] < -180.0)
x += 360.0;
}
}
}
phase[i1][i2] = Math.round(10 * x) / 10;
ta = ta + " " + x + "\n";
if (i1 == 0 && i2 == 0) {
p_min = phase[i1][i2];
p_max = phase[i1][i2];
}
else {
if (phase[i1][i2] > p_max)
p_max = phase[i1][i2];
else if (phase[i1][i2] < p_min)
p_min = phase[i1][i2];
}
// 次の周波数
ff += h;
}
}
document.getElementById("ans").value = ta;
// グラフの描画
// ゲイン線図
// タイトル
let gp1 = "7,ゲイン線図,角周波数,ゲイン（dB）," + n_g;
for (let i1 = 0; i1 < n_g; i1++)
gp1 = gp1 + "," + g_title[i1];
// ｘ軸目盛り
gp1 = gp1 + "," + fl + "," + fu + ",1.0";
let p = 0;
if (fl < 1.0) {
p = 1;
x = 0.1;
while (fl < x-1.0e-10) {
x /= 10.0;
p++;
}
}
gp1 = gp1 + "," + p;
// ｙ軸目盛り
if ((g_max-g_min) < 1.0e-5) {
if (g_min > 0.0) {
let k1 = Math.round(g_min / 20);
gp1 = gp1 + "," + (20.0 * k1 - 20.0);
gp1 = gp1 + "," + (20.0 * k1 + 20.0);
}
else {
let k1 = Math.round(-g_min / 20);
gp1 = gp1 + "," + (-20.0 * k1 - 20.0);
gp1 = gp1 + "," + (-20.0 * k1 + 20.0);
}
}
else {
if (g_min > 0.0) {
let k1 = Math.floor(g_min / 20);
gp1 = gp1 + "," + (20.0 * k1);
}
else {
let k1 = Math.floor(-g_min / 20);
if (Math.abs(-20.0*k1-g_min) > 1.0e-5)
k1++;
gp1 = gp1 + "," + (-20.0 * k1);
}
if (g_max > 0.0) {
let k1 = Math.floor(g_max / 20);
if (Math.abs(20.0*k1-g_max) > 1.0e-5)
k1++;
gp1 = gp1 + "," + (20.0 * k1);
}
else {
let k1 = Math.floor(-g_max / 20);
gp1 = gp1 + "," + (-20.0 * k1);
}
}
gp1 = gp1 + ",20.0,0";
// データ
gp1 = gp1 + "," + (k + 1);
for (let i1 = 0; i1 < n_g; i1++) {
for (let i2 = 0; i2 <= k; i2++)
gp1 = gp1 + "," + freq1[i1][i2];
for (let i2 = 0; i2 <= k; i2++)
gp1 = gp1 + "," + gain[i1][i2];
}
// 表示
gp1 = gp1 + ",1";
if (n_g == 1)
gp1 = gp1 + ",0";
else
gp1 = gp1 + ",1";
let str = "graph_js.htm?gp=" + gp1;
open(str, "gain", "width=950, height=700");
// 位相線図
// タイトル
let gp2 = "7,位相線図,角周波数,位相（度）," + n_g;
for (let i1 = 0; i1 < n_g; i1++)
gp2 = gp2 + "," + g_title[i1];
// ｘ軸目盛り
gp2 = gp2 + "," + fl + "," + fu + ",1.0," + p;
// ｙ軸目盛り
if ((p_max-p_min) < 1.0e-5) {
if (p_min > 0.0) {
let k1 = Math.round(p_min / 90);
gp2 = gp2 + "," + (90.0 * k1 - 90.0);
gp2 = gp2 + "," + (90.0 * k1 + 90.0);
}
else {
let k1 = Math.round(-p_min / 90);
gp2 = gp2 + "," + (-90.0 * k1 - 90.0);
gp2 = gp2 + "," + (-90.0 * k1 + 90.0);
}
}
else {
if (p_min > 0.0) {
let k1 = Math.floor(p_min / 90);
gp2 = gp2 + "," + (90.0 * k1);
}
else {
let k1 = Math.floor(-p_min / 90);
if (Math.abs(-90.0*k1-p_min) > 1.0e-5)
k1++;
gp2 = gp2 + "," + (-90.0 * k1);
}
if (p_max > 0.0) {
let k1 = Math.floor(p_max / 90);
if (Math.abs(90.0*k1-p_max) > 1.0e-5)
k1++;
gp2 = gp2 + "," + (90.0 * k1);
}
else {
let k1 = Math.floor(-p_max / 90);
gp2 = gp2 + "," + (-90.0 * k1);
}
}
gp2 = gp2 + ",90.0,0";
// データ
gp2 = gp2 + "," + (k + 1);
for (let i1 = 0; i1 < n_g; i1++) {
for (let i2 = 0; i2 <= k; i2++)
gp2 = gp2 + "," + freq2[i1][i2];
for (let i2 = 0; i2 <= k; i2++)
gp2 = gp2 + "," + phase[i1][i2];
}
// 表示
gp2 = gp2 + ",1";
if (n_g == 1)
gp2 = gp2 + ",0";
else
gp2 = gp2 + ",1";
str = "graph_js.htm?gp=" + gp2;
open(str, "phase", "width=950, height=700");
}
}
/************/
/* log10(x) */
/************/
function log10(x)
{
return Math.log(x) / Math.log(10.0);
}
/****************************************/
/* 伝達関数のｓにｊωを代入した値の計算 */
/* ff : ω（周波数） */
/* m : 分子の次数 */
/* si : 分子多項式の係数 */
/* n : 分母の次数 */
/* bo : 分母多項式の係数 */
/* return : 結果 */
/****************************************/
function G_s(ff, m, si, n, bo)
{
let f;
let x;
let y;
// 周波数を複素数に変換
f = new Complex (0.0, ff);
// 分子
x = value(f, m, si);
// 分母
y = value(f, n, bo);
return x.dev(y);
}
/**************************************/
/* 多項式のｓにｊωを代入した値の計算 */
/* f : ｊω（周波数，複素数） */
/* n : 多項式の次数 */
/* a : 多項式の係数 */
/* return : 結果 */
/**************************************/
function value(f, n, a)
{
let i1;
let k1;
let z = new Complex (0.0, 0.0);
for (i1 = 0; i1 <= n; i1++) {
let x = new Complex(a[i1], 0.0);
let y = f.pow(i1);
x = x.mul(y);
z = z.add(x);
}
return z;
}
/************************/
/* Complex オブジェクト */
/************************/
function Complex(rr, ii)
{
this.r = rr; // 実部
this.i = ii; // 虚部
}
/******************/
/* 複素数の絶対値 */
/******************/
Complex.prototype.abs = function()
{
return Math.sqrt(this.r * this.r + this.i * this.i);
}
/****************/
/* 複素数の角度 */
/****************/
Complex.prototype.angle = function()
{
let x;
if (this.r == 0.0 && this.i == 0.0)
x = 0.0;
else
x = Math.atan2(this.i, this.r);
return x;
}
/****************/
/* 複素数のｎ乗 */
/****************/
Complex.prototype.pow = function(n)
{
let i1;
let c = new Complex(1.0, 0.0);
if (n >= 0) {
for (i1 = 0; i1 < n; i1++)
c = c.mul(this);
}
else {
for (i1 = 0; i1 < -n; i1++)
c = c.mul(this);
b = new Complex(1.0, 0.0);
c = b.dev(c);
}
return c;
}
/****************/
/* 複素数の加算 */
/****************/
Complex.prototype.add = function(b)
{
let c = new Complex(0.0, 0.0);
c.r = this.r + b.r;
c.i = this.i + b.i;
return c;
}
/****************/
/* 複素数の減算 */
/****************/
Complex.prototype.sub = function(b)
{
let c = new Complex(0.0, 0.0);
c.r = this.r - b.r;
c.i = this.i - b.i;
return c;
}
/****************/
/* 複素数の乗算 */
/****************/
Complex.prototype.mul = function(b)
{
let c = new Complex(0.0, 0.0);
c.r = this.r * b.r - this.i * b.i;
c.i = this.i * b.r + this.r * b.i;
return c;
}
/****************/
/* 複素数の除算 */
/****************/
Complex.prototype.dev = function(b)
{
let c = new Complex(0.0, 0.0);
let x = 1.0 / (b.r * b.r + b.i * b.i);
c.r = (this.r * b.r + this.i * b.i) * x;
c.i = (this.i * b.r - this.r * b.i) * x;
return c;
}
/*****************************************/
/* 式の数や次数の変更 */
/* kind = -1 : 式の数の変更 */
/* >= : 次数の変更（式番号-1） */
/* sw = 1 : 分子 */
/* 2 : 分母 */
/*****************************************/
function change(kind, sw)
{
// 式の数
if (kind < 0) {
let str = document.getElementById("n_g").value;
// 初期状態
if (str == "") {
for (let i1 = 0; i1 < max_g; i1++) {
str = "div" + i1;
document.getElementById(str).style.display = "none";
}
}
// 式の数が入力されたとき
else {
n_g = parseInt(str);
for (let i1 = 0; i1 < n_g; i1++) {
str = "div" + i1;
document.getElementById(str).style.display = "";
str = "s" + i1;
if (document.getElementById(str).value == "") {
for (let i2 = 0; i2 <= max_order; i2++) {
str = "s" + i1 + i2 + "_t";
document.getElementById(str).style.display = "none";
}
}
else {
let n = parseInt(document.getElementById(str).value);
for (let i2 = 0; i2 < n; i2++) {
str = "s" + i1 + i2 + "_t";
document.getElementById(str).style.display = "";
}
for (let i2 = n; i2 <= max_order; i2++) {
str = "s" + i1 + i2 + "_t";
document.getElementById(str).style.display = "none";
}
}
str = "b" + i1;
if (document.getElementById(str).value == "") {
for (let i2 = 0; i2 <= max_order; i2++) {
str = "b" + i1 + i2 + "_t";
document.getElementById(str).style.display = "none";
}
}
else {
let n = parseInt(document.getElementById(str).value);
for (let i2 = 0; i2 < n; i2++) {
str = "b" + i1 + i2 + "_t";
document.getElementById(str).style.display = "";
}
for (let i2 = n; i2 <= max_order; i2++) {
str = "b" + i1 + i2 + "_t";
document.getElementById(str).style.display = "none";
}
}
}
for (let i1 = n_g; i1 < max_g; i1++) {
str = "div" + i1;
document.getElementById(str).style.display = "none";
}
}
}
// 分子の次数
else if (sw == 1) {
str = "s" + kind;
if (document.getElementById(str).value == "") {
for (let i1 = 0; i1 <= max_order; i1++) {
str = "s" + kind + i1 + "_t";
document.getElementById(str).style.display = "none";
}
}
else {
let n = parseInt(document.getElementById(str).value);
for (let i1 = 0; i1 <= n; i1++) {
str = "s" + kind + i1 + "_t";
document.getElementById(str).style.display = "";
}
for (let i1 = n+1; i1 <= max_order; i1++) {
str = "s" + kind + i1 + "_t";
document.getElementById(str).style.display = "none";
}
}
}
// 分母の次数
else {
str = "b" + kind;
if (document.getElementById(str).value == "") {
for (let i1 = 0; i1 <= max_order; i1++) {
str = "b" + kind + i1 + "_t";
document.getElementById(str).style.display = "none";
}
}
else {
let n = parseInt(document.getElementById(str).value);
for (let i1 = 0; i1 <= n; i1++) {
str = "b" + kind + i1 + "_t";
document.getElementById(str).style.display = "";
}
for (let i1 = n+1; i1 <= max_order; i1++) {
str = "b" + kind + i1 + "_t";
document.getElementById(str).style.display = "none";
}
}
}
}
/**************/
/* 画面の構成 */
/**************/
function s_area()
{
for (let i1 = 0; i1 < max_g; i1++) {
document.write(' <DIV ID="div' + i1 + '" STYLE="font-size: 100%; text-align:center; background-color: #ffffff; width: 800px; height:200px; margin-right: auto; margin-left: auto; border: 1px green solid; overflow: auto">\n');
document.write(' ' + (i1+1) + ' 番目の式の説明：<INPUT ID="exp' + i1 + '" STYLE="font-size: 100%" TYPE="text" SIZE="15" VALUE=""><BR>\n');
document.write(' <DIV ID="div' + i1 + '1" STYLE="font-size: 100%; text-align:center; background-color: #eeffee; width: 395px; height:150px; border: 1px black solid; overflow: auto; float: left;">\n');
document.write(' ' + (i1+1) + '.分子の次数：<INPUT ID="s' + i1 + '" STYLE="font-size: 100%" TYPE="text" SIZE="2" VALUE="" onBlur="change(' + i1 + ',1)"><BR>\n');
document.write(' <SPAN ID="s' + i1 + '0_t">定数項：<INPUT ID="s' + i1 + '0" STYLE="font-size: 100%" TYPE="text" SIZE="3" VALUE=""></SPAN><BR>\n');
for (let i2 = 1; i2 <= max_order; i2++)
document.write(' <SPAN ID="s' + i1 + i2 + '_t">s の ' + i2 + ' 次の項：<INPUT ID="s' + i1 + i2 + '" STYLE="font-size: 100%" TYPE="text" SIZE="3" VALUE=""></SPAN><BR>\n');
document.write(' </DIV>\n');
document.write(' <DIV ID="div' + i1 + '2" STYLE="font-size: 100%; text-align:center; background-color: #eeffee; width: 395px; height:150px; border: 1px black solid; overflow: auto; float: right;">\n');
document.write(' ' + (i1+1) + '.分母の次数：<INPUT ID="b' + i1 + '" STYLE="font-size: 100%" TYPE="text" SIZE="2" VALUE="" onBlur="change(' + i1 + ',2)"><BR>\n');
document.write(' <SPAN ID="b' + i1 + '0_t">定数項：<INPUT ID="b' + i1 + '0" STYLE="font-size: 100%" TYPE="text" SIZE="3" VALUE=""></SPAN><BR>\n');
for (let i2 = 1; i2 <= max_order; i2++)
document.write(' <SPAN ID="b' + i1 + i2 + '_t">s の ' + i2 + ' 次の項：<INPUT ID="b' + i1 + i2 + '" STYLE="font-size: 100%" TYPE="text" SIZE="3" VALUE=""></SPAN><BR>\n');
document.write(' </DIV>\n');
document.write(' </DIV>\n');
}
}
</SCRIPT>
</HEAD>
<BODY STYLE="font-size: 130%; text-align:center; background-color: #eeffee;">
<H2 STYLE="text-align:center"><B>ボード線図</B></H2>
式の数（グラフの数）：<INPUT ID="n_g" STYLE="font-size: 100%" TYPE="text" SIZE="2" VALUE="1" onBlur="change(-1, -1)">
<BUTTON STYLE="font-size: 100%; background-color: pink" onClick="return main()">実行</BUTTON><BR><BR>
周波数下限：<INPUT ID="low" STYLE="font-size: 100%" TYPE="text" SIZE="5" VALUE="0.01">
周波数上限：<INPUT ID="up" STYLE="font-size: 100%" TYPE="text" SIZE="5" VALUE="100">
データ数：<INPUT ID="n_data" STYLE="font-size: 100%" TYPE="text" SIZE="5" VALUE="100"><BR><BR>
<DIV STYLE="font-size: 100%; text-align:center; background-color: #ffffff; width: 820px; height:300px; margin-right: auto; margin-left: auto; border: 2px green solid; padding: 5px; overflow: auto">
<SCRIPT TYPE="text/javascript">
s_area();
change(-1, -1);
</SCRIPT>
</DIV><BR>
<TEXTAREA ID="ans" COLS="70" ROWS="10" STYLE="font-size: 100%;"></TEXTAREA>
</BODY>
</HTML>
<!DOCTYPE HTML>
<HTML>
<HEAD>
<TITLE>グラフの表示</TITLE>
<META HTTP-EQUIV="Content-Type" CONTENT="text/html; charset=utf-8">
<LINK REL="stylesheet" TYPE="text/css" HREF="../master.css">
<SCRIPT TYPE="text/javascript" SRC="graph.js"></SCRIPT>
<SCRIPT TYPE="text/javascript">
function GetParameter()
{
let result = new Array();
if(1 < window.location.search.length) {
// 最初の1文字 (?記号) を除いた文字列を取得する
let str = window.location.search.substring(1);
// 区切り記号 (&) で文字列を配列に分割する
let param = str.split('&');
for(let i1 = 0; i1 < param.length; i1++ ) {
// パラメータ名とパラメータ値に分割する
let element = param[i1].split('=');
let Name = decodeURIComponent(element[0]);
let Value = decodeURIComponent(element[1]);
// パラメータ名をキーとして連想配列に追加する
result[Name] = Value;
}
}
return result;
}
</SCRIPT>
</HEAD>
<BODY CLASS="white" STYLE="text-align: center">
<DIV ID="cl_line" STYLE="text-align: center; display: none">
<FORM>
<DIV ID="c0">
<INPUT ID="rgb0" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r0" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(0)">
緑<INPUT ID="g0" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(0)">
青<INPUT ID="b0" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(0)">
</DIV>
<DIV ID="c1">
<INPUT ID="rgb1" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r1" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(1)">
緑<INPUT ID="g1" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(1)">
青<INPUT ID="b1" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(1)">
</DIV>
<DIV ID="c2">
<INPUT ID="rgb2" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r2" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(2)">
緑<INPUT ID="g2" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(2)">
青<INPUT ID="b2" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(2)">
</DIV>
<DIV ID="c3">
<INPUT ID="rgb3" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r3" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(3)">
緑<INPUT ID="g3" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(3)">
青<INPUT ID="b3" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(3)">
</DIV>
<DIV ID="c4">
<INPUT ID="rgb4" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r4" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(4)">
緑<INPUT ID="g4" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(4)">
青<INPUT ID="b4" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(4)">
</DIV>
<DIV ID="c5">
<INPUT ID="rgb5" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r5" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(5)">
緑<INPUT ID="g5" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(5)">
青<INPUT ID="b5" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(5)">
</DIV>
<DIV ID="c6">
<INPUT ID="rgb6" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r6" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(6)">
緑<INPUT ID="g6" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(6)">
青<INPUT ID="b6" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(6)">
</DIV>
<DIV ID="c7">
<INPUT ID="rgb7" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r7" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(7)">
緑<INPUT ID="g7" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(7)">
青<INPUT ID="b7" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(7)">
</DIV>
<DIV ID="c8">
<INPUT ID="rgb8" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r8" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(8)">
緑<INPUT ID="g8" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(8)">
青<INPUT ID="b8" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(8)">
</DIV>
<DIV ID="c9">
<INPUT ID="rgb9" TYPE="text" SIZE="3" STYLE="font-size: 90%">
赤<INPUT ID="r9" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(9)">
緑<INPUT ID="g9" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(9)">
青<INPUT ID="b9" TYPE="text" SIZE="1" STYLE="font-size: 90%" onblur="c_change(9)">
</DIV>
<DIV ID="line_m">
マーク：<INPUT ID="l_m1" TYPE="radio" NAME="mark" STYLE="font-size: 90%" CHECKED>付ける
<INPUT ID="l_m2" TYPE="radio" NAME="mark" STYLE="font-size: 90%">付けない
</DIV>
<DIV ID="line_w">
線の太さ：<INPUT ID="l_w" TYPE="text" SIZE="3" STYLE="font-size: 90%" VALUE="2">
</DIV>
<DIV>
<SPAN STYLE="background-color: pink; font-size: 100%" onClick="D_Change()">OK</SPAN>
</DIV>
</FORM>
</DIV>
<BR>
<DIV STYLE="text-align: center">
<CANVAS ID="canvas_e" STYLE="background-color: #eeffee;" WIDTH="900" HEIGHT="600" onClick="Click(event)"></CANVAS>
</DIV>
<SCRIPT TYPE="text/javascript">
let result = GetParameter();
graph(result['gp']);
</SCRIPT>
</BODY>
</HTML>
<?php
/********************************/
/* 伝達関数のゲインと位相の計算 */
/* coded by Y.Suganuma */
/********************************/
/****************************/
/* クラスComplex */
/* coded by Y.Suganuma */
/****************************/
class Complex {
public $r; // 実部
public $i; // 虚部
/******************/
/* コンストラクタ */
/* a : 実部 */
/* b : 虚部 */
/******************/
function Complex($a, $b = 0.0)
{
$this->r = $a;
$this->i = $b;
}
}
/******************/
/* 複素数の絶対値 */
/******************/
function abs_c($a)
{
return sqrt($a->r * $a->r + $a->i * $a->i);
}
/****************/
/* 複素数の角度 */
/****************/
function angle($a)
{
$x = 0.0;
if ($a->r != 0.0 || $a->i != 0.0)
$x = atan2($a->i, $a->r);
return $x;
}
/****************/
/* 複素数のｎ乗 */
/****************/
function pow_c($a, $n)
{
$c = new Complex(1.0);
if ($n >= 0) {
for ($i1 = 0; $i1 < $n; $i1++)
$c = mul($c, $a);
}
else {
for ($i1 = 0; $i1 < -$n; $i1++)
$c = mul($c, $a);
$c = div(new Complex(1.0), $c);
}
return $c;
}
/****************/
/* 複素数の加算 */
/****************/
function add($a, $b)
{
$c = new Complex(0.0);
$c->r = $a->r + $b->r;
$c->i = $a->i + $b->i;
return $c;
}
/****************/
/* 複素数の減算 */
/****************/
function sub($a, $b)
{
$c = new Complex(0.0);
$c->r = $a->r - $b->r;
$c->i = $a->i - $b->i;
return $c;
}
/****************/
/* 複素数の乗算 */
/****************/
function mul($a, $b)
{
$c = new Complex(0.0);
$c->r = $a->r * $b->r - $a->i * $b->i;
$c->i = $a->i * $b->r + $a->r * $b->i;
return $c;
}
/****************/
/* 複素数の除算 */
/****************/
function div($a, $b)
{
$c = new Complex(0.0);
$x = 1.0 / ($b->r * $b->r + $b->i * $b->i);
$c->r = ($a->r * $b->r + $a->i * $b->i) * $x;
$c->i = ($a->i * $b->r - $a->r * $b->i) * $x;
return $c;
}
/********************************************************************/
/* 伝達関数のｓにｊωを代入した値の計算 */
/* ff : ω（周波数） */
/* ms : 分子の表現方法 */
/* =0 : 多項式（si[0]+si[1]*s+si[1]*s^2+･･･） */
/* =1 : 因数分解（(si[0]+si[1]*s)*(si[2]+si[3]*s)*･･･） */
/* m : 分子の次数 */
/* si : 分子多項式の係数 */
/* mb : 分母の表現方法 */
/* =0 : 多項式（bo[0]+bo[1]*s+bo[1]*s^2+･･･） */
/* =1 : 因数分解（(bo[0]+bo[1]*s)*(bo[2]+bo[3]*s)*･･･） */
/* n : 分母の次数 */
/* bo : 分母多項式の係数 */
/* return : 結果 */
/********************************************************************/
function G_s($ff, $ms, $m, $si, $mb, $n, $bo)
{
// 周波数を複素数に変換
$f = new Complex(0.0, $ff);
// 分子
$x = value($f, $ms, $m, $si);
// 分母
$y = value($f, $mb, $n, $bo);
return div($x, $y);
}
/****************************************************************/
/* 多項式のｓにｊωを代入した値の計算 */
/* f : ｊω（周波数，複素数） */
/* ms : 多項式の表現方法 */
/* =0 : 多項式（a[0]+a[1]*s+a[1]*s^2+･･･） */
/* =1 : 因数分解（(a[0]+a[1]*s)*(a[2]+a[3]*s)*･･･） */
/* n : 多項式の次数 */
/* a : 多項式の係数 */
/* return : 結果 */
/****************************************************************/
function value($f, $ms, $n, $a)
{
$u0 = new Complex(0.0);
$u1 = new Complex(1.0);
// 因数分解されていない
if ($ms == 0) {
$x = $u0;
for ($i1 = 0; $i1 <= $n; $i1++)
$x = add($x, mul(new Complex($a[$i1]), pow_c($f, $i1)));
}
// 因数分解されている
else {
if ($n == 0)
$x = mul(new Complex($a[0]), $u1);
else {
$x = $u1;
$k1 = 2 * $n;
for ($i1 = 0; $i1 < $k1; $i1 += 2)
$x = mul($x, add(new Complex($a[$i1]), mul(new Complex($a[$i1+1]), $f)));
}
}
return $x;
}
/*******************/
/* main プログラム */
/*******************/
$phase = 0.0;
$uc = 90.0 / asin(1.0);
/*
データの入力
*/
// 出力ファイル名の入力とファイルのオープン
fscanf(STDIN, "%*s %s %*s %s", $o_fg, $o_fp);
$og = fopen($o_fg, "w");
$op = fopen($o_fp, "w");
// 伝達関数データ
fscanf(STDIN, "%*s %lf %lf %*s %d", $fl, $fu, $k);
// 分子
$str = trim(fgets(STDIN));
strtok($str, " ");
$ms = intval(strtok(" "));
strtok(" ");
$m = intval(strtok(" "));
strtok(" ");
$si = array();
// 因数分解されていない
if ($ms == 0) {
$k1 = $m + 1;
for ($i1 = 0; $i1 < $k1; $i1++)
$si[$i1] = floatval(strtok(" "));
}
// 因数分解されている
else {
$k1 = ($m == 0) ? 1 : 2 * $m;
for ($i1 = 0; $i1 < $k1; $i1++)
$si[$i1] = floatval(strtok(" "));
}
// 分母
$str = trim(fgets(STDIN));
strtok($str, " ");
$mb = intval(strtok(" "));
strtok(" ");
$n = intval(strtok(" "));
strtok(" ");
$bo = array();
// 因数分解されていない
if ($mb == 0) {
$k1 = $n + 1;
for ($i1 = 0; $i1 < $k1; $i1++)
$bo[$i1] = floatval(strtok(" "));
}
// 因数分解されている
else {
$k1 = ($n == 0) ? 1 : 2 * $n;
for ($i1 = 0; $i1 < $k1; $i1++)
$bo[$i1] = floatval(strtok(" "));
}
/*
ゲインと位相の計算
*/
$h = (log10($fu) - log10($fl)) / $k;
$ff = log10($fl);
for ($i1 = 0; $i1 <= $k; $i1++) {
// 周波数の対数を元に戻す
$f = pow(10.0, $ff);
// 値の計算
$g = G_s($f, $ms, $m, $si, $mb, $n, $bo);
// ゲインと位相の計算
$gain = 20.0 * log10(abs_c($g));
$x = angle($g) * $uc;
while (abs($phase-$x) > 180.0) {
if ($x-$phase > 180.0)
$x -= 360.0;
else {
if ($x-$phase < -180.0)
$x += 360.0;
}
}
$phase = $x;
// 出力
fwrite($og, $f." ".$gain."\n");
fwrite($op, $f." ".$phase."\n");
// 次の周波数
$ff += $h;
}
/*
------------------------data1.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 0 次数 0 係数（次数が低い順） 1.0
分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0
------------------------data2.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 1 次数 0 係数（次数が低い順） 1.0
分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0
*/
?>
#*******************************/
# 伝達関数のゲインと位相の計算 */
# coded by Y.Suganuma */
#*******************************/
#***************************/
# クラスComp */
# coded by Y.Suganuma */
#***************************/
class Comp
# コンストラクタ
def initialize(_real, _imag)
@_real = _real
@_imag = _imag
end
# 演算子のオーバーロード
def +(obj)
c = Comp.new(@_real+obj._real, @_imag+obj._imag)
return c
end
def -(obj)
c = Comp.new(@_real-obj._real, @_imag-obj._imag)
return c
end
def *(obj)
r = @_real * obj._real - @_imag * obj._imag
i = @_imag * obj._real + @_real * obj._imag
c = Comp.new(r, i)
return c
end
def /(obj)
x = 1.0 / (obj._real * obj._real + obj._imag * obj._imag)
r = (@_real * obj._real + @_imag * obj._imag) * x
i = (@_imag * obj._real - @_real * obj._imag) * x
c = Comp.new(r, i)
return c
end
# 複素数の絶対値
def c_abs()
x = Math.sqrt(@_real * @_real + @_imag * @_imag)
return x
end
# 複素数の角度
def angle()
x = 0.0
if @_real != 0.0 || @_imag != 0.0
x = Math.atan2(@_imag, @_real)
end
return x
end
# 出力
def out(cr = "")
printf "(%f, %f)%s", @_real, @_imag, cr
end
attr_accessor("_real", "_imag")
end
#**************/
# 複素数のｎ乗
#**************/
def pow(a, n)
c = Comp.new(1, 0)
if n >= 0
for i1 in 0 ... n
c = c * a
end
else
for i1 in 0 ... -n
c = c * a
end
c1 = Comp.new(1, 0)
c = c1 / c
end
return c
end
#*******************************************************************/
# 伝達関数のｓにｊωを代入した値の計算 */
# ff : ω（周波数） */
# ms : 分子の表現方法 */
# =0 : 多項式（si[0]+si[1]*s+si[1]*s^2+･･･） */
# =1 : 因数分解（(si[0]+si[1]*s)*(si[2]+si[3]*s)*･･･） */
# m : 分子の次数 */
# si : 分子多項式の係数 */
# mb : 分母の表現方法 */
# =0 : 多項式（bo[0]+bo[1]*s+bo[1]*s^2+･･･） */
# =1 : 因数分解（(bo[0]+bo[1]*s)*(bo[2]+bo[3]*s)*･･･） */
# n : 分母の次数 */
# bo : 分母多項式の係数 */
# return : 結果 */
#*******************************************************************/
def G_s(ff, ms, m, si, mb, n, bo)
# 周波数を複素数に変換
f = Comp.new(0.0, ff)
# 分子
x = value(f, ms, m, si)
# 分母
y = value(f, mb, n, bo)
return x / y
end
#***************************************************************/
# 多項式のｓにｊωを代入した値の計算 */
# f : ｊω（周波数，複素数） */
# ms : 多項式の表現方法 */
# =0 : 多項式（a[0]+a[1]*s+a[1]*s^2+･･･） */
# =1 : 因数分解（(a[0]+a[1]*s)*(a[2]+a[3]*s)*･･･） */
# n : 多項式の次数 */
# a : 多項式の係数 */
# return : 結果 */
#***************************************************************/
def value(f, ms, n, a)
u0 = Comp.new(0.0, 0.0)
u1 = Comp.new(1.0, 0.0)
# 因数分解されていない
if ms == 0
x = u0
for i1 in 0 ... n+1
x = x + Comp.new(a[i1], 0.0) * pow(f, i1)
end
# 因数分解されている
else
if n == 0
x = Comp.new(a[0], 0.0) * u1
else
x = u1
k1 = 2 * n
i1 = 0
while i1 < k1
x = x * (Comp.new(a[i1], 0.0) + Comp.new(a[i1+1], 0.0) * f)
i1 += 2
end
end
end
return x
end
uc = 90.0 / Math.asin(1.0) # 単位変換用
#
# データの入力
#
# 出力ファイル名の入力とファイルのオープン
s = gets().split(" ")
o_fg = s[1]
o_fp = s[3]
og = open(o_fg, "w")
op = open(o_fp, "w")
# 伝達関数データ
s = gets().split(" ")
fl = Float(s[1])
fu = Float(s[2])
k = Integer(s[4])
# 分子
s = gets().split(" ")
ms = Integer(s[1])
m = Integer(s[3])
# 因数分解されていない
if ms == 0
k1 = m + 1
si = Array.new(k1)
for i1 in 0 ... k1
si[i1] = Float(s[i1+5])
end
# 因数分解されている
else
k1 = (m == 0) ? 1 : 2 * m
si = Array.new(k1)
for i1 in 0 ... k1
si[i1] = Float(s[i1+5])
end
end
# 分母
s = gets().split(" ")
mb = Integer(s[1])
n = Integer(s[3])
# 因数分解されていない
if mb == 0
k1 = n + 1
bo = Array.new(k1)
for i1 in 0 ... k1
bo[i1] = Float(s[i1+5])
end
# 因数分解されている
else
k1 = (n == 0) ? 1 : 2 * n
bo = Array.new(k1)
for i1 in 0 ... k1
bo[i1] = Float(s[i1+5])
end
end
#
# ゲインと位相の計算
#
phase = 0.0
h = (Math.log10(fu) - Math.log10(fl)) / k
ff = Math.log10(fl)
for i1 in 0 ... k+1
# 周波数の対数を元に戻す
f = 10.0 ** ff
# 値の計算
g = G_s(f, ms, m, si, mb, n, bo)
# ゲインと位相の計算
gain = 20.0 * Math.log10(g.c_abs())
x = g.angle * uc
while (phase-x).abs() > 180.0
if (x-phase > 180.0)
x -= 360.0
else
if (x-phase < -180.0)
x += 360.0
end
end
end
phase = x
# 出力
og.printf("%f %f\n", f, gain)
op.printf("%f %f\n", f, phase)
# 次の周波数
ff += h
end
=begin
------------------------data1.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 0 次数 0 係数（次数が低い順） 1.0
分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0
------------------------data2.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 1 次数 0 係数（次数が低い順） 1.0
分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0
=end
#*******************************/
# 伝達関数のゲインと位相の計算 */
# coded by Y.Suganuma */
#*******************************/
require 'complex'
#*******************************************************************/
# 伝達関数のｓにｊωを代入した値の計算 */
# ff : ω（周波数） */
# ms : 分子の表現方法 */
# =0 : 多項式（si[0]+si[1]*s+si[1]*s^2+･･･） */
# =1 : 因数分解（(si[0]+si[1]*s)*(si[2]+si[3]*s)*･･･） */
# m : 分子の次数 */
# si : 分子多項式の係数 */
# mb : 分母の表現方法 */
# =0 : 多項式（bo[0]+bo[1]*s+bo[1]*s^2+･･･） */
# =1 : 因数分解（(bo[0]+bo[1]*s)*(bo[2]+bo[3]*s)*･･･） */
# n : 分母の次数 */
# bo : 分母多項式の係数 */
# return : 結果 */
#*******************************************************************/
def G_s(ff, ms, m, si, mb, n, bo)
# 周波数を複素数に変換
f = Complex(0.0, ff)
# 分子
x = value(f, ms, m, si)
# 分母
y = value(f, mb, n, bo)
return x / y
end
#***************************************************************/
# 多項式のｓにｊωを代入した値の計算 */
# f : ｊω（周波数，複素数） */
# ms : 多項式の表現方法 */
# =0 : 多項式（a[0]+a[1]*s+a[1]*s^2+･･･） */
# =1 : 因数分解（(a[0]+a[1]*s)*(a[2]+a[3]*s)*･･･） */
# n : 多項式の次数 */
# a : 多項式の係数 */
# return : 結果 */
#***************************************************************/
def value(f, ms, n, a)
u0 = Complex(0.0, 0.0)
u1 = Complex(1.0, 0.0)
# 因数分解されていない
if ms == 0
x = u0
for i1 in 0 ... n+1
x = x + a[i1] * (f ** i1)
end
# 因数分解されている
else
if n == 0
x = a[0] * u1
else
x = u1
k1 = 2 * n
i1 = 0
while i1 < k1
x = x * (a[i1] + a[i1+1] * f)
i1 += 2
end
end
end
return x
end
uc = 90.0 / Math.asin(1.0) # 単位変換用
#
# データの入力
#
# 出力ファイル名の入力とファイルのオープン
s = gets().split(" ")
o_fg = s[1]
o_fp = s[3]
og = open(o_fg, "w")
op = open(o_fp, "w")
# 伝達関数データ
s = gets().split(" ")
fl = Float(s[1])
fu = Float(s[2])
k = Integer(s[4])
# 分子
s = gets().split(" ")
ms = Integer(s[1])
m = Integer(s[3])
# 因数分解されていない
if ms == 0
k1 = m + 1
si = Array.new(k1)
for i1 in 0 ... k1
si[i1] = Float(s[i1+5])
end
# 因数分解されている
else
k1 = (m == 0) ? 1 : 2 * m
si = Array.new(k1)
for i1 in 0 ... k1
si[i1] = Float(s[i1+5])
end
end
# 分母
s = gets().split(" ")
mb = Integer(s[1])
n = Integer(s[3])
# 因数分解されていない
if mb == 0
k1 = n + 1
bo = Array.new(k1)
for i1 in 0 ... k1
bo[i1] = Float(s[i1+5])
end
# 因数分解されている
else
k1 = (n == 0) ? 1 : 2 * n
bo = Array.new(k1)
for i1 in 0 ... k1
bo[i1] = Float(s[i1+5])
end
end
#
# ゲインと位相の計算
#
phase = 0.0
h = (Math.log10(fu) - Math.log10(fl)) / k
ff = Math.log10(fl)
for i1 in 0 ... k+1
# 周波数の対数を元に戻す
f = 10.0 ** ff
# 値の計算
g = G_s(f, ms, m, si, mb, n, bo)
# ゲインと位相の計算
gain = 20.0 * Math.log10(g.abs())
x = g.arg * uc
while (phase-x).abs() > 180.0
if (x-phase > 180.0)
x -= 360.0
else
if (x-phase < -180.0)
x += 360.0
end
end
end
phase = x
# 出力
og.printf("%f %f\n", f, gain)
op.printf("%f %f\n", f, phase)
# 次の周波数
ff += h
end
=begin
------------------------data1.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 0 次数 0 係数（次数が低い順） 1.0
分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0
------------------------data2.txt----------------------------
ゲイン gain.txt 位相 phase.txt
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 1 次数 0 係数（次数が低い順） 1.0
分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0
=end
# -*- coding: UTF-8 -*- import numpy as np import sys import math import cmath ################################ # 伝達関数のゲインと位相の計算 # coded by Y.Suganuma ################################ ################################################################ # 多項式のｓにｊωを代入した値の計算 # f : ｊω（周波数，複素数） # ms : 多項式の表現方法 # =0 : 多項式（a[0]+a[1]*s+a[1]*s^2+･･･） # =1 : 因数分解（(a[0]+a[1]*s)*(a[2]+a[3]*s)*･･･） # n : 多項式の次数 # a : 多項式の係数 # return : 結果 ################################################################ def value(f, ms, n, a) : u0 = complex(0.0, 0.0) u1 = complex(1.0, 0.0) # 因数分解されていない if ms == 0 : x = u0 for i1 in range(0, n+1) : x = x + a[i1] * (f ** i1) # 因数分解されている else : if n == 0 : x = a[0] * u1 else : x = u1 k1 = 2 * n for i1 in range(0, k1, 2) : x = x * (a[i1] + a[i1+1] * f) return x #################################################################### # 伝達関数のｓにｊωを代入した値の計算 # ff : ω（周波数） # ms : 分子の表現方法 # =0 : 多項式（si[0]+si[1]*s+si[1]*s^2+･･･） # =1 : 因数分解（(si[0]+si[1]*s)*(si[2]+si[3]*s)*･･･） # m : 分子の次数 # si : 分子多項式の係数 # mb : 分母の表現方法 # =0 : 多項式（bo[0]+bo[1]*s+bo[1]*s^2+･･･） # =1 : 因数分解（(bo[0]+bo[1]*s)*(bo[2]+bo[3]*s)*･･･） # n : 分母の次数 # bo : 分母多項式の係数 # return : 結果 #################################################################### def G_s(ff, ms, m, si, mb, n, bo) : # 周波数を複素数に変換 f = complex (0.0, ff) # 分子 x = value(f, ms, m, si) # 分母 y = value(f, mb, n, bo) return x / y phase = 0.0 # データの入力 # 出力ファイル名の入力とファイルのオープン line = sys.stdin.readline() s = line.split() o_fg = s[1] o_fp = s[3] og = open(o_fg, "w") op = open(o_fp, "w") # 伝達関数データ line = sys.stdin.readline() s = line.split() fl = float(s[1]) fu = float(s[2]) k = int(s[4]) # 分子 line = sys.stdin.readline() s = line.split() ms = int(s[1]) m = int(s[3]) # 因数分解されていない if ms == 0 : k1 = m + 1 si = np.empty(k1, np.float) for i1 in range(0, k1) : si[i1] = float(s[i1+5]) # 因数分解されている else : if m == 0 : k1 = 1 else : k1 = 2 * m si = np.empty(k1, np.float) for i1 in range(0, k1) : si[i1] = float(s[i1+5]) # 分母 line = sys.stdin.readline() s = line.split() mb = int(s[1]) n = int(s[3]) # 因数分解されていない if mb == 0 : k1 = n + 1 bo = np.empty(k1, np.float) for i1 in range(0, k1) : bo[i1] = float(s[i1+5]) # 因数分解されている else : if n == 0 : k1 = 1 else : k1 = 2 * n bo = np.empty(k1, np.float) for i1 in range(0, k1) : bo[i1] = float(s[i1+5]) # ゲインと位相の計算 h = (math.log10(fu) - math.log10(fl)) / k ff = math.log10(fl) for i1 in range(0, k+1) : # 周波数の対数を元に戻す f = 10.0 ** ff # 値の計算 g = G_s(f, ms, m, si, mb, n, bo) # ゲインと位相の計算 gain = 20.0 * math.log10(abs(g)) x = math.degrees(cmath.phase(g)) while abs(phase-x) > 180.0 : if x-phase > 180.0 : x -= 360.0 else : if x-phase < -180.0 : x += 360.0 phase = x # 出力 og.write(str(f) + " " + str(gain) + "\n") op.write(str(f) + " " + str(phase) + "\n") # 次の周波数 ff += h og.close() op.close() """ ------------------------data1.txt---------------------------- ゲイン gain.txt 位相 phase.txt 周波数の下限と上限 0.01 100 分割数 100 分子の表現方法 0 次数 0 係数（次数が低い順） 1.0 分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0 ------------------------data2.txt---------------------------- ゲイン gain.txt 位相 phase.txt 周波数の下限と上限 0.01 100 分割数 100 分子の表現方法 1 次数 0 係数（次数が低い順） 1.0 分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0 """
/********************************/
/* 伝達関数のゲインと位相の計算 */
/* coded by Y.Suganuma */
/********************************/
using System;
using System.IO;
class Program
{
/****************/
/* main program */
/****************/
static void Main(String[] args)
{
double uc = 90.0 / Math.Asin(1.0);
char[] charSep = new char[] {' '};
//
// データの入力
//
// 出力ファイル名とファイルのオープン
String[] str = Console.ReadLine().Split(charSep, StringSplitOptions.RemoveEmptyEntries);
StreamWriter g_o = new StreamWriter(str[1]);
StreamWriter p_o = new StreamWriter(str[3]);
// 伝達関数データ
str = Console.ReadLine().Split(charSep, StringSplitOptions.RemoveEmptyEntries);
double fl = double.Parse(str[1]);
double fu = double.Parse(str[2]);
int k = int.Parse(str[4]);
// 分子
str = Console.ReadLine().Split(charSep, StringSplitOptions.RemoveEmptyEntries);
int ms = int.Parse(str[1]);
int m = int.Parse(str[3]);
// 因数分解されていない
double[] si;
if (ms == 0) {
int k1 = m + 1;
si = new double [k1];
for (int i1 = 0; i1 < k1; i1++)
si[i1] = double.Parse(str[i1+5]);
}
// 因数分解されている
else {
int k1 = (m == 0) ? 1 : 2 * m;
si = new double [k1];
for (int i1 = 0; i1 < k1; i1++)
si[i1] = double.Parse(str[i1+5]);
}
// 分母
str = Console.ReadLine().Split(charSep, StringSplitOptions.RemoveEmptyEntries);
int mb = int.Parse(str[1]);
int n = int.Parse(str[3]);
// 因数分解されていない
double[] bo;
if (mb == 0) {
int k1 = n + 1;
bo = new double [k1];
for (int i1 = 0; i1 < k1; i1++)
bo[i1] = double.Parse(str[i1+5]);
}
// 因数分解されている
else {
int k1 = (n == 0) ? 1 : 2 * n;
bo = new double [k1];
for (int i1 = 0; i1 < k1; i1++)
bo[i1] = double.Parse(str[i1+5]);
}
//
// インと位相の計算
//
double h = (Math.Log10(fu) - Math.Log10(fl)) / k;
double ff = Math.Log10(fl);
double phase = 0.0;
for (int i1 = 0; i1 <= k; i1++) {
// 周波数の対数を元に戻す
double f = Math.Pow(10.0, ff);
// 値の計算
Complex g = G_s(f, ms, m, si, mb, n, bo);
// ゲインと位相の計算
double gain = 20.0 * Math.Log10(Complex.abs(g));
double x = Complex.angle(g) * uc;
while (Math.Abs(phase-x) > 180.0) {
if (x-phase > 180.0)
x -= 360.0;
else {
if (x-phase < -180.0)
x += 360.0;
}
}
phase = x;
// 出力
g_o.WriteLine(f + " " + gain);
p_o.WriteLine(f + " " + phase);
// 次の周波数
ff += h;
}
g_o.Close();
p_o.Close();
}
/********************************************************************/
/* 伝達関数のｓにｊωを代入した値の計算 */
/* ff : ω（周波数） */
/* ms : 分子の表現方法 */
/* =0 : 多項式（si[0]+si[1]*s+si[1]*s^2+･･･） */
/* =1 : 因数分解（(si[0]+si[1]*s)*(si[2]+si[3]*s)*･･･） */
/* m : 分子の次数 */
/* si : 分子多項式の係数 */
/* mb : 分母の表現方法 */
/* =0 : 多項式（bo[0]+bo[1]*s+bo[1]*s^2+･･･） */
/* =1 : 因数分解（(bo[0]+bo[1]*s)*(bo[2]+bo[3]*s)*･･･） */
/* n : 分母の次数 */
/* bo : 分母多項式の係数 */
/* return : 結果 */
/********************************************************************/
static Complex G_s(double ff, int ms, int m, double[] si, int mb, int n, double[] bo)
{
// 周波数を複素数に変換
Complex f = new Complex (0.0, ff);
// 分子
Complex x = value(f, ms, m, si);
// 分母
Complex y = value(f, mb, n, bo);
return Complex.dev(x, y);
}
/****************************************************************/
/* 多項式のｓにｊωを代入した値の計算 */
/* f : ｊω（周波数，複素数） */
/* ms : 多項式の表現方法 */
/* =0 : 多項式（a[0]+a[1]*s+a[1]*s^2+･･･） */
/* =1 : 因数分解（(a[0]+a[1]*s)*(a[2]+a[3]*s)*･･･） */
/* n : 多項式の次数 */
/* a : 多項式の係数 */
/* return : 結果 */
/****************************************************************/
static Complex value(Complex f, int ms, int n, double[] a)
{
Complex x;
// 因数分解されていない
if (ms == 0) {
x = new Complex (0.0, 0.0);
for (int i1 = 0; i1 <= n; i1++)
x = Complex.add(x, Complex.mul(new Complex(a[i1]), Complex.pow(f, i1)));
}
// 因数分解されている
else {
if (n == 0)
x = new Complex (a[0]);
else {
x = new Complex (1.0, 0.0);
int k1 = 2 * n;
for (int i1 = 0; i1 < k1; i1 += 2)
x = Complex.mul(x, Complex.add(new Complex(a[i1]), Complex.mul(new Complex(a[i1+1]), f)));
}
}
return x;
}
}
/****************************/
/* クラスComplex */
/* coded by Y.Suganuma */
/****************************/
class Complex {
double r; // 実部
double i; // 虚部
/******************/
/* コンストラクタ */
/* a : 実部 */
/* b : 虚部 */
/******************/
public Complex(double a = 0.0, double b = 0.0)
{
r = a;
i = b;
}
/********/
/* 出力 */
/********/
public void output(StreamWriter OUT, Complex a)
{
OUT.Write("(" + a.r + ", " + a.i + ")");
}
/******************/
/* 複素数の絶対値 */
/******************/
public static double abs(Complex a)
{
return Math.Sqrt(a.r * a.r + a.i * a.i);
}
/****************/
/* 複素数の角度 */
/****************/
public static double angle(Complex a)
{
double x = 0.0;
if (a.r != 0.0 || a.i != 0.0)
x = Math.Atan2(a.i, a.r);
return x;
}
/****************/
/* 複素数のｎ乗 */
/****************/
public static Complex pow(Complex a, int n)
{
Complex c = new Complex (1);
if (n >= 0) {
for (int i1 = 0; i1 < n; i1++)
c = Complex.mul(c, a);
}
else {
for (int i1 = 0; i1 < -n; i1++)
c = Complex.mul(c, a);
c = Complex.dev(new Complex(1.0), c);
}
return c;
}
/****************/
/* 複素数の加算 */
/****************/
public static Complex add(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r + b.r;
c.i = a.i + b.i;
return c;
}
/****************/
/* 複素数の減算 */
/****************/
public static Complex sub(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r - b.r;
c.i = a.i - b.i;
return c;
}
/****************/
/* 複素数の乗算 */
/****************/
public static Complex mul(Complex a, Complex b)
{
Complex c = new Complex();
c.r = a.r * b.r - a.i * b.i;
c.i = a.i * b.r + a.r * b.i;
return c;
}
/****************/
/* 複素数の除算 */
/****************/
public static Complex dev(Complex a, Complex b)
{
Complex c = new Complex();
double x = 1.0 / (b.r * b.r + b.i * b.i);
c.r = (a.r * b.r + a.i * b.i) * x;
c.i = (a.i * b.r - a.r * b.i) * x;
return c;
}
}
/*
------------------------data1----------------------------
ゲイン gain 位相 phase
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 0 次数 0 係数（次数が低い順） 1.0
分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0
------------------------data2----------------------------
ゲイン gain 位相 phase
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 1 次数 0 係数（次数が低い順） 1.0
分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0
*/
'******************************'
' 伝達関数のゲインと位相の計算 '
' coded by Y.Suganuma '
'******************************'
Imports System.IO
Imports System.Text.RegularExpressions
Module Test
Sub Main()
Dim uc As Double = 90.0 / Math.Asin(1.0)
'
' データの入力
'
' 出力ファイル名とファイルのオープン
Dim MS1 As Regex = New Regex("\s+")
Dim str() As String = MS1.Split(Console.ReadLine().Trim())
Dim g_o As StreamWriter = new StreamWriter(str(1))
Dim p_o As StreamWriter = new StreamWriter(str(3))
' 伝達関数データ
str = MS1.Split(Console.ReadLine().Trim())
Dim fl As Double = Double.Parse(str(1))
Dim fu As Double = Double.Parse(str(2))
Dim k As Integer = Integer.Parse(str(4))
' 分子
str = MS1.Split(Console.ReadLine().Trim())
Dim ms As Integer = Integer.Parse(str(1))
Dim m As Integer = Integer.Parse(str(3))
' 因数分解されていない
Dim si() As Double
If ms = 0
Dim k1 As Integer = m + 1
ReDim si(k1)
For i1 As Integer = 0 To k1-1
si(i1) = Double.Parse(str(i1+5))
Next
' 因数分解されている
Else
Dim k1 As Integer
If m = 0
k1 = 1
Else
k1 = 2 * m
End If
ReDim si(k1)
For i1 As Integer = 0 To k1-1
si(i1) = Double.Parse(str(i1+5))
Next
End If
' 分母
str = MS1.Split(Console.ReadLine().Trim())
Dim mb As Integer = Integer.Parse(str(1))
Dim n As Integer = Integer.Parse(str(3))
' 因数分解されていない
Dim bo() As Double
If mb = 0
Dim k1 As Integer = n + 1
ReDim bo(k1)
For i1 As Integer = 0 To k1-1
bo(i1) = Double.Parse(str(i1+5))
Next
' 因数分解されている
Else
Dim k1 As Integer
If n = 0
k1 = 1
Else
k1 = 2 * n
End If
ReDim bo(k1)
For i1 As Integer = 0 To k1-1
bo(i1) = Double.Parse(str(i1+5))
Next
End If
'
' インと位相の計算
'
Dim h As Double = (Math.Log10(fu) - Math.Log10(fl)) / k
Dim ff As Double = Math.Log10(fl)
Dim phase As Double = 0.0
For i1 As Integer = 0 To k
' 周波数の対数を元に戻す
Dim f As Double = Math.Pow(10.0, ff)
' 値の計算
Dim g As Complex = G_s(f, ms, m, si, mb, n, bo)
' ゲインと位相の計算
Dim gain As Double = 20.0 * Math.Log10(g.abs())
Dim x As Double = g.angle() * uc
Do While Math.Abs(phase-x) > 180.0
If x-phase > 180.0
x -= 360.0
Else
If x-phase < -180.0
x += 360.0
End If
End If
Loop
phase = x
' 出力
g_o.WriteLine(f & " " & gain)
p_o.WriteLine(f & " " & phase)
' 次の周波数
ff += h
Next
g_o.Close()
p_o.Close()
End Sub
'******************************************************************'
' 伝達関数のｓにｊωを代入した値の計算 '
' ff : ω（周波数） '
' ms : 分子の表現方法 '
' =0 : 多項式（si(0)+si(1)*s+si(1)*s^2+･･･） '
' =1 : 因数分解（(si(0)+si(1)*s)*(si(2)+si(3)*s)*･･･） '
' m : 分子の次数 '
' si : 分子多項式の係数 '
' mb : 分母の表現方法 '
' =0 : 多項式（bo(0)+bo(1)*s+bo(1)*s^2+･･･） '
' =1 : 因数分解（(bo(0)+bo(1)*s)*(bo(2)+bo(3)*s)*･･･） '
' n : 分母の次数 '
' bo : 分母多項式の係数 '
' return : 結果 '
'******************************************************************'
Function G_s(ff As Double, ms As Integer, m As Integer, si() As Double, mb As Integer, n As Integer, bo() As Double)
' 周波数を複素数に変換
Dim f As Complex = new Complex (0.0, ff)
' 分子
Dim x As Complex = value(f, ms, m, si)
' 分母
Dim y As Complex = value(f, mb, n, bo)
Return c_dev(x, y)
End Function
'**************************************************************'
' 多項式のｓにｊωを代入した値の計算 '
' f : ｊω（周波数，複素数） '
' ms : 多項式の表現方法 '
' =0 : 多項式（a(0)+a(1)*s+a(1)*s^2+･･･） '
' =1 : 因数分解（(a(0)+a(1)*s)*(a(2)+a(3)*s)*･･･） '
' n : 多項式の次数 '
' a : 多項式の係数 '
' return : 結果 '
'**************************************************************'
Function value(f As Complex, ms As Integer, n As Integer, a() As Double)
Dim x As Complex
' 因数分解されていない
If ms = 0
x = new Complex (0.0, 0.0)
For i1 As Integer = 0 To n
x = c_add(x, c_mul(new Complex(a(i1)), c_pow(f, i1)))
Next
' 因数分解されている
Else
If n = 0
x = new Complex (a(0))
Else
x = new Complex (1.0, 0.0)
Dim k1 As Integer = 2 * n
For i1 As Integer = 0 To k1-1 Step 2
x = c_mul(x, c_add(new Complex(a(i1)), c_mul(new Complex(a(i1+1)), f)))
Next
End If
End If
Return x
End Function
'**************'
' 複素数のｎ乗 '
'**************'
Function c_pow(a As Complex, n As Integer)
Dim c As Complex = new Complex (1)
If n >= 0
For i1 As Integer = 0 To n-1
c = c_mul(c, a)
Next
Else
For i1 As Integer = 0 To -n-1
c = c_mul(c, a)
Next
c = c_dev(new Complex(1.0), c)
End If
Return c
End Function
'**************'
' 複素数の加算 '
'**************'
Function c_add(a As Complex, b As Complex)
Dim c As Complex = new Complex()
c.r = a.r + b.r
c.i = a.i + b.i
Return c
End Function
'**************'
' 複素数の減算 '
'**************'
Function c_sub(a As Complex, b As Complex)
Dim c As Complex = new Complex()
c.r = a.r - b.r
c.i = a.i - b.i
Return c
End Function
'**************'
' 複素数の乗算 '
'**************'
Function c_mul(a As Complex, b As Complex)
Dim c As Complex = new Complex()
c.r = a.r * b.r - a.i * b.i
c.i = a.i * b.r + a.r * b.i
Return c
End Function
'**************'
' 複素数の除算 '
'**************'
Function c_dev(a As Complex, b As Complex)
Dim c As Complex = new Complex()
Dim x As Double = 1.0 / (b.r * b.r + b.i * b.i)
c.r = (a.r * b.r + a.i * b.i) * x
c.i = (a.i * b.r - a.r * b.i) * x
Return c
End Function
'**************************'
' クラスComplex '
' coded by Y.Suganuma '
'**************************'
Class Complex
Public r As Double ' 実部
Public i As Double ' 虚部
'****************'
' コンストラクタ '
' a : 実部 '
' b : 虚部 '
'****************'
Public Sub New(Optional a As Double = 0.0, Optional b As Double = 0.0)
r = a
i = b
End Sub
'******'
' 出力 '
'******'
Sub output(OUT As StreamWriter)
OUT.Write("(" & r & ", " & i & ")")
End Sub
'****************'
' 複素数の絶対値 '
'****************'
Public Function abs()
Return Math.Sqrt(r * r + i * i)
End Function
'**************'
' 複素数の角度 '
'**************'
Public Function angle()
Dim x As Double = 0.0
If r <> 0.0 or i <> 0.0
x = Math.Atan2(i, r)
End If
Return x
End Function
End Class
End Module
/*
------------------------data1----------------------------
ゲイン gain 位相 phase
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 0 次数 0 係数（次数が低い順） 1.0
分母の表現方法 0 次数 3 係数（次数が低い順） 1.0 11.1 11.1 1.0
------------------------data2----------------------------
ゲイン gain 位相 phase
周波数の下限と上限 0.01 100 分割数 100
分子の表現方法 1 次数 0 係数（次数が低い順） 1.0
分母の表現方法 1 次数 3 係数（次数が低い順） 1.0 0.1 1.0 1.0 1.0 10.0
*/
| 静岡理工科大学 | 菅沼ホーム | 目次 | 索引 |