Wednesday, August 5, 2009

Making functions

A important resource in any programming language is to make specific functions.
Using Scilab we may make our own functions.

For example, if we need a function that adds two numbers and divides the result by two, we can call that function as add_d2(v1, v2).

The follow script is made in the editor.


function [result] = add_d2(v1, v2)
aux = v1 + v2;
result = aux/2;
endfunction;


Or, we can write just one line of code:


function [result] = add_d2(v1, v2)
result = (v1 + v2)/2;
endfunction;


Obvious, the example is too simple, but I'd like to give a example for present the syntax.

Ok, let's make a better example:

We need a function that finds the peaks in a signal and returns two variable, one with the position and the other with the value of each peak.


function [positions, values] = find_peaks(x)
positions = [];
values = [];

n_peaks = 0;
for i = 2:length(x) - 1, // x should be a vector, like: x = [x_1 x_2 x_3 ... x_n]
if (x(i) > x(i - 1)) & (x(i) > x(i + 1)),
n_peaks = n_peaks + 1;
positions(n_peaks) = i;
values(n_peaks) = x(i);
end;
end;
endfunction;


Test the code and say me what happens.

21 comments:

Unknown said...

Hi, Sheep;

I am new in Scilab and very interested in knowing more about Scilab. And I found your blog where you are intended to teach newbie like me. I hope I can learn more from you.

I have tried your find-peaks function. Somehow it's doesn't really work out. Can you teach me how should I do it?

I have save your find_peaks function and load it to the Scilab console when I use it.

I created a signal as below:

t=0:0.001:10;
x=0.3*sin(t) + sin(1.3*t) + 0.9*sin(4.2*t) + 0.02*rand(1, 10001);
plot2d(x);

Try to find the peaks of the signal

[position, value]=find_peaks(x);

but I get nothing from the function.

Am I doing something wrong? Please guide me.

One more question, how should I modify the code if I would like to detect all the peaks & valleys generated?

Looking forward to hear from you soon.

Alex Carneiro said...

Hi CV, the problem in your example is the rand(.) component in 'x'.

Try eliminate the rand(.) component, like this:

x=0.3*sin(t) + sin(1.3*t) + 0.9*sin(4.2*t);

And make the graphs:

plot(t, x);

plot(t(position), values, '.');


Now, for the valleys:

function [positions, values] = find_valleys(x)
positions = [];
values = [];

n_valleys = 0;
for i = 2:length(x) - 1, // x should be a vector, like: x = [x_1 x_2 x_3 ... x_n]
if (x(i) < x(i - 1)) & (x(i) < x(i + 1)),
n_valleys = n_valleys + 1;
positions(n_valleys) = i;
values(n_valleys) = x(i);
end;
end;
endfunction;


Is it ok? You can ask me for anything.

Unknown said...

Hi Sheep;

Thanks for you soon reply.
Thank you so much.
Is very kind of you in guiding me.
I have made the changes that you have suggested on the rand(.)
I could plot the graph.
But when I use the function to look for the peaks, it was fail.

I try to check the value of both [position] and [values], but I get nothing

This is what I have key in

-->t=0:0.001:10;

-->x=0.3*sin(t) + sin(1.3*t) + 0.9*sin(4.2*t);

-->plot(t,x);

-->find_peaks(x);

-->plot(t(position), values, '.');
!--error 4
Undefined variable: position

-->position
!--error 4
Undefined variable: position

-->values
!--error 4
Undefined variable: values

What have I done wrong?
Something wrong with the way I save the function?
Please guide me.

Looking forward to hear from you soon

Alex Carneiro said...

Hi again CV, you forgot this:

[position, value] = find_peaks(x); //to set the variables 'position' and 'value'

Change and try it.

See you.

Unknown said...

Hi Sheep,

Gracias!! Thank you!!
Thank you for you kindness in guiding a slow learner like me.
Because of you, I love Scilab more!

Now I manage to find the peaks and valleys now!!

Stupid me, you're absolutely right, I forgot to assign the position and value into the array.

The peaks&valley topic is solved. Thanks Sheep!

May I ask other question? Is it possible to find the full width at half maximum(FWHM) of the waveform?

Looking forward to hear from you soon

Alex Carneiro said...

Hello CV, I think it's possible to find the FWHD of the waveform.

Look, you have 't', where f(t) is a peak, so you need only to find 'Dt' that f(t+Dt) is a half of f(t).

Try something like this:
// 'x' is f(t)
// 'tp' is the position of a peak
t_aux = 1;
while x(tp) > 2*x(tp + t_aux),
t_aux = t_aux + 1;
end;

plot(x);
plot((tp - t_aux):(tp + t_aux), x((tp - t_aux):(tp + t_aux)), 'r-.');


Say me if it works.
See you soon.

Alex Carneiro said...

I made a mistake, change the signal '>' to '<' in the line that starts with while.

Unknown said...

Hi Sheep,

Thanks for your kind reply.I try to make a fwhm function by referring to the method you have suggested.This is what I have wrote:

function [fwhm_positions, fwhm_values] = fwhm(positions)
fwhm_positions[];
fwhm_values[];

tp=posionts;//assign the peak position into tp

t_aux = 1;
while x(tp) < 2*x(tp + t_aux),
t_aux = t_aux + 1;
end;
endfunction;

When I try to call it in the console, this is the output that I got:

[fwhm_positions, fwhm_values] = fwhm(positions);
!--error 4
Undefined variable: fwhm


fwhm_positions[];
!--error 276
Missing operator, comma, or semicolon.

at line 2 of function fwhm called by :
endfunction;
at line 11 of exec file called by :
xec); disp(msprintf(gettext("E
while executing a callback

Can you guide me?

Looking forward to hear from you soon.

Unknown said...

Hi Sheep,

There was a typo in my previous syntax.

tp=posionts;
change to
tp=positions;

some of the errors were cleared but the below still remain:

fwhm_positions[];
!--error 276
Missing operator, comma, or semicolon.

at line 2 of function fwhm called by :
endfunction;
at line 11 of exec file called by :
xec); disp(msprintf(gettext("E
while executing a callback

what other errors that I have made?

Looking forward to hear from you soon

Alex Carneiro said...

Dear CV, in Scilab you don't declare variables like:

vector[];

in Scilab you only use the variable like this:

vector = [1 3 5 7 10];

Now, let's see your code:

function [fwhm_positions, fwhm_values] = fwhm(positions)
fwhm_positions[]; // not right
fwhm_values[]; // not right

tp=posionts;//assign the peak position into tp // it's changed

t_aux = 1;
while x(tp) < 2*x(tp + t_aux),
t_aux = t_aux + 1;
end;
//where do you set 'fwhm_positions' and 'fwhm_values'?
// try it: fwhm_positions = (tp - t_aux):(tp + t_aux); fwhm_values = x((tp - t_aux):(tp + t_aux))
endfunction;


And call the fwhm function like this:

[fwhm_positions fwhm_values] = fwhm(x, positions);

Don't use global variables.

That's all. See you.

Unknown said...

Dear Sheep;

Thanks for your kindness. Thank you for pointing out my mistake. It helps for a beginner like me! This round, I can load it to console without trouble.

I have rewrite the function as below:

function [fwhm_positions, fwhm_values] = fwhm(x,positions)


tp=positions;//assign the peak position into tp

t_aux = 1;
while x(tp) < 2*x(tp + t_aux),
t_aux = t_aux + 1;
end;

fwhm_positions = (tp - t_aux):(tp + t_aux);
fwhm_values = x((tp - t_aux):(tp + t_aux));
endfunction;

when I call the function as below:

[fwhm_positions, fwhm_values] = fwhm(x, positions);

I get this output:

!--error 204
':': Wrong type for argument 1: Scalar expected.

at line 12 of function fwhm called by :
s] = fwhm(x, positions)

what are the mistakes that I have made?

Please guide me.

Looking forward to hear from you soon

Alex Carneiro said...

Dear CV, is everything fine? I think so.

Say me, how are your variables 'x' and 'positions'?

'x' must be a vector, e. g. x = exp(-0.1*([1:100] - 15).^2); //'x' can't have negative values

and 'positions' must be a scalar, e. g. positions = 15;

Now, try it:

[fwhm_positions, fwhm_values] = fwhm(x, positions);

plot(x);
plot(fwhm_positions, fwhm_values, 'r.-');

If you've studied the Gaussian function, you can make changes.
Else, whatever (no problem).


Bye my friend.

Unknown said...

Dear Sheep;

How was your weekend?

Now, I could find the fwhm of the peaks!
Thanks for you advices.
I have learned so much from you. Thanks for always patience with such a slow learner like me.

Please keep up the good work in the blog. I will always visit your blog.

Omi said...

Missing operator, comma, or semicolon.
i m getting this error on using the given code what is wrong could u plzz help me?
x=4;y=2;(x>y)?4:2

Alex Carneiro said...

Omi, Scilab does not have ternary operator.

I suggest you try it:

x = 4;
y = 2;
max([x y])

Best regards

Unknown said...

actually i have been asked this question n dre r oders too....could u plzz help me with their answers...actually i m new to scilab n i dont know any of it...but i need to find out dis answers ASAP
its MCQs...

1. Predict the output of: x=4;y=2;(x>y)?4;2 in scilab
a. 4 b. 2 c. 6 d. 5

2.which of the following is true about functions in scilab
a. they are all user defined
b. none of these
c. user defined functions are a subset of all functions
d. functions do not

3. identify the file format generated from scilab screen dumps to verify the commands used during a sesion
a) .exe b) .scc c) .sce d).sss

4. scilab comment begins with
a) % b) // c)both a & b d) none of these

5. in scilab vectors are a particular case of matrices, where the no of rows(or the no of columns) is equal to
a) 2 b)1 c)3 d)4

Unknown said...

6. the correct way to define complex no in scilab is
a) no.=2+3i b)no.=2+3*i c)no.=2+3%i d)no.=complex(2,3)

my answers to the respective questions and the (response i got in paranthesis)
(1) 2(wrong)
(2) they are all user defined(wrong)
(3) .sss(wrong)
(4) //(wrong)
(5) 2(wrong)
(6) no.=complex(2,3)(wrong)

its omi here

thank you for replying me...but could u plzzz help me wid dis its reallyy imp for me...i need answers to it asap

Unknown said...

Plzzz rep me

Unknown said...

Hi Alex Carneiro,
please help me wid my doubts....

anyone...

Unknown said...

plzz reply me

hannabeprakash said...

I didnt get the graph for the 2d plot with values in it. I see only the marker'.'. Then I wanted to know how to get the minimum and maximum value of a 3d plot :

//Peaks
function values = find_peaks(z)
posx = [];
posy = [];
values = [];
[posx,posy]=size(z);
n_peaks = 0;
for i=2:posx - 1,
for j = 2:posy - 1, // x should be a vector, like: x = [x_1 x_2 x_3 ... x_n]
if (z(i,j) > z(i-1,j)) & (z(i,j) > z(i+1,j)) then
n_peaks = n_peaks + 1;
posx(n_peaks) = i;
posy(n_peak) = j;
//disp('max',x(i));
values(n_peaks) = z(i,j);//this is the value to be seen on the plot
disp('max', z(i,j));
end;

end;
end;
endfunction;
//Valleys
function values = find_valleys(z)
posx = [];
posy = [];
values = [];
[posx,posy]=size(z);
n_valleys = 0;
for i=2:posx - 1,
for j = 2:posy - 1, // x should be a vector, like: x = [x_1 x_2 x_3 ... x_n]
if (z(i,j) < z(i-1,j)) & (z(i,j) < z(i+1,j)) then
n_valleys = n_valleys + 1;
posx(n_valleys) = i;
posy(n_valleys) = j;
//disp('min',x(i));
values(n_valleys) = z(i,j);//this is the value to be seen on the plot
disp('min',z(i,j));
end;

end;
end
endfunction;
xd=linspace(-1,1);
yd=linspace(0,2*%pi);
[x,y]=meshgrid(xd,yd);
z=x.*cos(y);
value1=find_peaks(z);
value2=find_valleys(z);
disp('peaks', value1);
disp('valleys', value2);
clf
surf(x,y,z)