# Plot the Bode, Black, Nichols and Nyquist diagrams of a transfer function from an XCOS diagram

**URL:** <https://scilab.discourse.group/t/plot-the-bode-black-nichols-and-nyquist-diagrams-of-a-transfer-function-from-an-xcos-diagram/541>\
**Category:** Xcos\
**Created:** [June 5, 2024, 7:34pm UTC](https://scilab.discourse.group/t/plot-the-bode-black-nichols-and-nyquist-diagrams-of-a-transfer-function-from-an-xcos-diagram/541 "2024-06-05T19:34:19Z")\
**Posts on this page:** 5\
**Page:** 1

<div class="post-metadata">

**Author:** ![Pascal77](https://avatars.discourse-cdn.com/v4/letter/p/a9a28c/32.png) [@Pascal77](https://scilab.discourse.group/u/Pascal77)\
**Post date:** [June 5, 2024, 7:34pm UTC](https://scilab.discourse.group/t/plot-the-bode-black-nichols-and-nyquist-diagrams-of-a-transfer-function-from-an-xcos-diagram/541/1 "2024-06-05T19:34:19Z")

</div>

This code can be used to plot the Bode, Black, Nichols and Nyquist diagrams of a transfer function from an XCOS diagram. It is not necessary to open the XCOS file before running the script.  
The post processing “post\_xcos\_simulate” functionality is used for this operation.

```auto
function post_xcos_simulate(%cpr, scs_m, needcompile)
// Do Nothing
endfunction

//We load the XCOS file to be simulated
importXcosDiagram("C:\\...path...\Exemple diagramme de Bode depuis XCOS/Essai_BODE_plot.zcos")

// Clear any pre_xcos_simulation;
clear pre_xcos_simulate;

function post_xcos_simulate(%cpr, scs_m, needcompile)
    // Retrieve all objects
    objs = scs_m.objs;

    clrBlock = [];
    //Looking for CLR block
    for i=1:size(objs)
        if objs(i).gui == "CLR" then
            clrBlock = objs(i);
            break;
        end
    end

    // Retrieve exprs
    s = poly(0,'s');
    expr = clrBlock.graphics.exprs;
    num = evstr(expr(1));
    den = evstr(expr(2));
    h = syslin('c', num/den);
    
    // Open new figure then plot Bode, Black, Nyquist
    clf();//Clear all previous curves
    f = gcf();// Create a figure
    f.figure_size=[1200,1200];
    // f window background color
    f.background=color('lightsteelblue');
    f.figure_name='BODE BLACK NICHOLS NYQUIST DIAGRAMS FROM XCOS FILE'

    subplot(221);
    bode(h, 0.01, 100,"Mon filtre");
   
    subplot(222);
    black(h,0.01,10,"Mon filtre");
    nicholschart(colors=color('gray')*[1 1]);
  
    subplot(223);
    nyquist(h,0.01,10,"Mon filtre");
    hallchart(colors=color('blue')*[1 1])
endfunction

//We run the simulation of the Bode and other diagrams 
xcos_simulate(scs_m, 4);

```

XCOS diagram is the following one:  
 ![2024-06-05 21_31_23-Essai_BODE_plot (C__Users_Pascal_Documents_SciLab files_SCILAB from 2024 version](https://global.discourse-cdn.com/free1/uploads/utc/original/1X/efe00f1216f7408ad6f7b9df081042f6cf211778.png)

The plot is depicted below:

 ![2024-06-05 21_29_19-BODE BLACK NICHOLS NYQUIST DIAGRAMS FROM XCOS FILE](https://global.discourse-cdn.com/free1/uploads/utc/original/1X/91edbabc7bddc037614f3e723d2e471c90e42965.png)

Of course, it is possible to split the code in order to keep only one kind of diagram.

Example:

```auto
importXcosDiagram("C:\\...path...\Essai_BODE_plot.zcos")

// Clear any pre_xcos_simulation;
clear pre_xcos_simulate;

function post_xcos_simulate(%cpr, scs_m, needcompile)
    // Retrieve all objects
    objs = scs_m.objs;

    clrBlock = [];
    //Looking for CLR block
    for i=1:size(objs)
        if objs(i).gui == "CLR" then
            clrBlock = objs(i);
            break;
        end
    end

    // Retrieve exprs
    s = poly(0,'s');
    expr = clrBlock.graphics.exprs;
    num = evstr(expr(1));
    den = evstr(expr(2));
    h = syslin('c', num/den);
    // Open new figure then plot Bode
    scf(max(winsid())+1);
    bode(h, 0.01, 100,"Mon filtre");
endfunction

//Bode simulation running
xcos_simulate(scs_m, 4);

```

![Bode plot with xcos_simulate](https://global.discourse-cdn.com/free1/uploads/utc/original/1X/4a97b055bfe2f3fe8e861c6565adf1b6650951b6.png)

---

<div class="post-metadata">

**Author:** ![pascal](https://yyz2.discourse-cdn.com/free1/user_avatar/scilab.discourse.group/pascal/32/40_2.png) [@pascal](https://scilab.discourse.group/u/pascal)\
**Post date:** [June 6, 2024, 6:09am UTC](https://scilab.discourse.group/t/plot-the-bode-black-nichols-and-nyquist-diagrams-of-a-transfer-function-from-an-xcos-diagram/541/2 "2024-06-06T06:09:14Z")

</div>

actually one plot with Magnitude and Phase as second Y-Axis would be nice, so you have all the information at once. 😁

---

<div class="post-metadata">

**Author:** ![Pascal77](https://avatars.discourse-cdn.com/v4/letter/p/a9a28c/32.png) [@Pascal77](https://scilab.discourse.group/u/Pascal77)\
**Post date:** [June 6, 2024, 10:43am UTC](https://scilab.discourse.group/t/plot-the-bode-black-nichols-and-nyquist-diagrams-of-a-transfer-function-from-an-xcos-diagram/541/3 "2024-06-06T10:43:45Z")

</div>

Hello, I have prepared a first try (not perfect) with a common single axis for both plots : Magnitude and Phase.

I used `gainplot` and `phaseplot` instructions because it is very convenient but it is also possible to do this with the classical formula 20.log (|h|) and arg(h), but I did not explore this way to do.

Here is the code (I am started from the previous file) :

```auto
// Loading the XCOS diagram:
importXcosDiagram("C:\\..path...\Exemple diagramme de Bode depuis XCOS/Essai_BODE_plot.zcos")

// Clear any pre_xcos_simulation;
clear pre_xcos_simulate;

function post_xcos_simulate(%cpr, scs_m, needcompile)
    // Retrieve all objects
    objs = scs_m.objs;

    clrBlock = [];
    //Looking for CLR block
    for i=1:size(objs)
        if objs(i).gui == "CLR" then
            clrBlock = objs(i);
            break;
        end
    end

    // Retrieve exprs
    s = poly(0,'s');
    expr = clrBlock.graphics.exprs;
    num = evstr(expr(1));
    den = evstr(expr(2));
    h = syslin('c', num/den);
    // Open new figure then plot Bode
    //scf(max(winsid())+1);
    //bode(h, 0.01, 100,"Mon filtre");
    //----------------------------------------------------------------------
    //Plot magnitude and phase on 1 plot with 2 y axis
    clf;
    f=gcf();
    f.figure_size=[1000,800];
        a=gca();
    a.background=color('dimgrey');// or any other predefined color name (in color_list)
    a.foreground=color('blue');
    
    gainplot(h,0.01,100);
    phaseplot(h,0.01,100);
    
    //first curve
    poly1= a.children.children(1); //store polyline handle into poly1
    poly1.foreground = 2; // change the color...
    poly1.thickness = 2; // ...and the thickness of a curve.

    ylabel(". Magnitude (dB) & Phase (deg) . ");ylabel edgecolor b background grey95;  
   
    //second curve
    poly2= a.children.children(2); //store polyline handle into poly1
    poly2.foreground = 21; // change the color...
    poly2.thickness = 2; // ...and the thickness of a curve.    
    
    legend("Magnitude (dB)","Phase (deg)");
    t = gce(); // get the handle of the newly created object
    t.font_foreground=32; // change font properties
    t.font_size=2;
    t.font_style=9;  

    //title('Bode plot of h(s) = s^2+2*0.9*10*s+100)/(s^2+2*0.3*10.1*s+102.01');//identifying the curve 1 on the plot
    title(["$\frac{1000}{1+20*s+s^2}$"]);
    title font_foreground gold
    title edgecolor red
    title fontsize 4
    title backgroundcolor dimgrey
    title font_style 8
    
    //----------------------------------------------------------------------
endfunction

//Bode simulation running
xcos_simulate(scs_m, 4);

```

 ![2024-06-06 12_38_58-Bode Diagram on the same plot with 2 y axis (post_xcos_simulate)png](https://global.discourse-cdn.com/free1/uploads/utc/original/1X/63ddccdfc32c93a69cd649a4a8fe1abbc703800b.png)

---

<div class="post-metadata">

**Author:** ![pascal](https://yyz2.discourse-cdn.com/free1/user_avatar/scilab.discourse.group/pascal/32/40_2.png) [@pascal](https://scilab.discourse.group/u/pascal)\
**Post date:** [June 6, 2024, 12:06pm UTC](https://scilab.discourse.group/t/plot-the-bode-black-nichols-and-nyquist-diagrams-of-a-transfer-function-from-an-xcos-diagram/541/4 "2024-06-06T12:06:46Z")

</div>

is it possible to do the second Y (Phase on the right side of the graph?)

---

<div class="post-metadata">

**Author:** ![Pascal77](https://avatars.discourse-cdn.com/v4/letter/p/a9a28c/32.png) [@Pascal77](https://scilab.discourse.group/u/Pascal77)\
**Post date:** [June 6, 2024, 9:50pm UTC](https://scilab.discourse.group/t/plot-the-bode-black-nichols-and-nyquist-diagrams-of-a-transfer-function-from-an-xcos-diagram/541/5 "2024-06-06T21:50:02Z")

</div>

Here is the code to plot magnitude and phase on the same plot with mag axis at letf and phase axis at right.

```auto
// plot the magnitude and phase versus frequency
//cf. repfreq help 

clf;
f1=gcf();//
f1.figure_size=[1000,1000];//figure size
f1.figure_position=[100,100];//figure position
f1.background=color('gray90');//color surrounding the plots

//Magnitude and phase computation

/*
f=linspace(fmin,fmax,n); // frequency vector (Hz) (fmin,fmax,number of points)
w=f*2*%pi;//w=pulsation in rad/sec
Pmag=zeros(1,n);Pph=zeros(1,n);//we define the two vectors magnitude and phase (filled with zeros at the start)
for k=1:n, // compute P(jw) at each frequency point
  P=1/(-10*w(k)^2+0.1*%i*w(k)+1); //transfer function in jω (complex expression)
  [Pmag(k),Pph(k)]=dbphi(P);//get magnitude in db and phase in degrees from P (we fill in the two vectors Pmag and Php defined above))
end
*/

s=poly(0,'s');//definition of Laplace variable
h=syslin('c',1/(10*s^2+0.1*s+1));//definition of transfer function with a polynome (rational fraction)

fmin=0.001;//Hz
fmax=10;//Hz
f=fmin:0.001:fmax;//Hz (fmin:step:fmax)
w=f*2*%pi;// conversion in degrees
[fr1,rep] =repfreq(h,f);//definition of frequency response (rep)
[db,phi]=dbphi(rep);//magnitude and phase extraction from rep 

//Plotting curves
//1st curve
ymin = -180;
ymax = 40;
plot2d(f,db,logflag = "ln",rect=[fmin,ymin,fmax,ymax]); // magnitude plot (same scale than phase plot) rect=[xmin,ymin,xmax,ymax]
xlabel("Frequency (Hz)");
ylabel("Magnitude (dB)");
a = gca();//The left axis is defined (magnitude in dB)
aa = a.children(1).children(1);//this represents the magnitude track itself
aa.foreground=color('red');//we decided to paint it red
aa.thickness=2;//with a thickness of 2
xgrid();//plot the grid
legend("Magnitude (dB)",1);

title("$\frac{1}{10s^2+0.1s+1}$","FontSize",4');//We express the formula for the transfer function in s and jω (in Latex style)
title color blue
title edgecolor blue

//2nd curve
b=newaxes();//We define the second axis (phase in degrees) which will be on the right
plot2d(f,phi,logflag = "ln",rect=[fmin,ymin,fmax,ymax]); // phase plot (same scale than magnitude plot)
b.filled="off";
b.axes_visible(1)="off";//make y2 axis transparent to avoid masking y1 axis
b.y_location="right";//put y2 axis on the right
ylabel("Phase (degree)");
bb = b.children(1).children(1);//this represents the phase track itself
bb.foreground=color('blue');//we decided to paint it blue
bb.thickness=2;//with a thickness of 2
legend("Phase (°)",4);

```

 ![2024-06-06 23_44_43-Mag & phase curves on the same plots with 2 y axis](https://global.discourse-cdn.com/free1/uploads/utc/original/1X/eac2f3cf29467e83af34d4d8d4cac092122c54e9.png)

Below the code with `bode` instruction. Very simple compared to the previous one.

```auto
s = poly(0, 's');
h1 = syslin('c', 1/(10*s^2+0.1*s+1));

clf(); 
bode(h1, 0.01, 10);

```

And the corresponding plots:

 ![2024-06-06 23_48_28-bode plot](https://global.discourse-cdn.com/free1/uploads/utc/original/1X/c2a7fe2b12f953474ddfc6d3772752e282e6531f.png)

```auto

```