matlab calculation area circular

I have a signal and want to calculate the area. The problem is that i want to make a loop with decimal step, because the length is 1189,58

Any ideas how to proceed with the loop?

[Rik]

Comment posted as answer (and erroneously flagged by the spam filter) by Sirius8:

Noone? Come on guys a little help!

[Rik]

Why did you delete all your comments? Now most of my comments have a greatly reduced value for future readers. You were told the same thing 2 years ago here.

Answer 1

[Rik]

You can use trapz for a numerical approximation of an integration. By selecting your bounds you can use this to calculate the area under a curve.

Edit:

Using the code you posted and adding the loop I'm getting something similar to what you describe. If you don't want the legend entries, you should use explicit handles in your call to legend.

S=load('test1.mat');
[x1,y1,y11]=deal(S.x1,S.y1,S.y11);
figure(1),clf(1)

subplot(2,1,1)
plot(x1,y1)
grid on
hold all 
plot(x1,y11,'m')
grid on
hold all 
legend ('normalized fft','smoothed')
subplot(2,1,2)
plot(x1,y11,'m')
grid on
hold all 
legend ('smoothed')

cmap=colormap('prism');n_col=0;
freq_range=5;
for x_start=0:freq_range:max(x1)
    x_end=x_start+freq_range;
    L=x1>=x_start & x1<=x_end;
    if sum(L)<=1,continue,end%skip if there is only 1 point
    X=[x_start x1(L)  x_end];
    Y=[0       y11(L) 0];
    A=trapz(x1(L),y11(L));
    n_col=n_col+1;C=cmap(n_col,:);
    patch(X,Y,C,'DisplayName',sprintf('A(%d-%dHz)=%.2f',x_start,x_end,A))
end

[Rik]

Decimal values in Matlab follow the American notation: 3/10 is written as 0.3 not as 0,3.

You should be aware that you can't use decimal values as indices, if either of those are arrays.

It also looks like you're inputting scalars in trapz, instead of a range of values.

[Rik]

What are you trying to do? Why did you write your code this way? I have the feeling it wouldn't help if I make the edits to this code to make it run without errors. And how did you create these two variables or functions?

[Rik]

What code did you use to plot your lower line chart? If you used variables, can you attach a mat file with the amplitude and frequency?

[Rik]

Thank you for posting the data. I will get back to you on this later.

[Rik]

The wedge is there because of how the x-values are selected in the loop. If you want to avoid them, you have to make sure there is a y-value for those x-values.

If you want to store the calculated area: nothing is stopping you. If you read the code you can already see how it is calculated. You can simply store it in a vector.

[Rik]

Yes, but there isn't an x-value of 10 with the corresponding y-value. That means that my code will leave this gap. If you can't resample the curve, you could also round the boundaries to the nearest limit value. Note that this will actually increase the difference between the true value and the calculated value.

%untested coded
for n=1:(max(x1)-freq_range)/freq_range
    x_start=(n-1)*freq_range
    x_end=n*freq_range;
    
    [~,ind]=min(abs(x_start-x1));
    x_start=x1(ind);
    [~,ind]=min(abs(x_end-x1));
    x_end=x1(ind);
    
    %rest of the loop
end

[Rik]

That question is not about Matlab. The unit would be amplitude/s (i.e. the product of the two axes).

[John D'Errico]

I would add that a highly osccilatory curve like that is usually poorly integrated using a spline.

This is because the spline is trying to fit noise. Worse, if the spline actually predicts a negative result in some places, you have negative area being summed in there, when it is surely not appropriate.

In the end, a spline, even a smoothing spline, just adds variance to the result.

These are the curves where you want to use trapz - a big reason for the presence of MATLAB. trapz is NOT just another pretty toy, put there for no good reason except for some students to use. In fact, trapz is a minimal variance estimator for an integral on noisy data. So a very useful tool.

If you insist on using the smoothing spline though, there is no reason why you could not have just used fnint, which does come with the curve fitting toolbox and should work on a smoothing spline.

[Rik]

That looks like you still need to divide by the interval width. If that is the case, you can simply square your data and use tools like trapz to calculate the sum divided by width.

[Rik]

Exactly like what I posted in my answer, but replacing this

A=trapz(x1(L),y11(L));

with this

A=trapz(x1(L),y11(L).^2);

I'm assuming y11 only contains real numbers, otherwise you would still need abs to get the absolute value.

[Rik]

The division by width of the frequency interval (NB: not the total energy) is happening inside trapz, as that is the whole point of using it. You aren't dividing by any total in the formula you posted either, so I assumed that wasn't required. I'm not familiar enough with the mathematics involved in your specific application, so you need to check the math on your own (or with peers). I can help you with implementing it in Matlab.

[Rik]

You really shouldn't comment that line, since it is meaningles to calculate the area of a line.

I'm not sure I understand what your problem is. If you want to divide by the total area, what is stopping you? You're storing the partial areas is a vector, so it trivial to divide that vector after the loop by its sum.

[Rik]

First calculate E(1), E(2), etc, the divide E by sum(E).

[Rik]

What happens when you set freq_range=1;?

[Rik]

I would have assumed that the area this code calculates also changes. That is what happens in my copy of Matlab. Just store those in an array and devide that array by its sum after filling all elements.

[Rik]

No, because that will overwrite A in every loop iteration. In a previous comment you already showed you understand how to index a variable.

After the loop you can indeed use A/sum(A).

[Rik]

Since you don't index A anywhere in that code: no.

for n=1:45
    x_start=(n-1)*freq_range
    x_end=n*freq_range
        [~,ind]=min(abs(x_start-x1));
    x_start=x1(ind);
    [~,ind]=min(abs(x_end-x1));
    x_end=x1(ind);
    L=x1>=x_start & x1<=x_end;
    if sum(L)<=1,continue,end%skip if there is only 1 point
    X=[x_start x1(L)  x_end];
    Y=[0       y11(L) 0];
    A=trapz(x1(L),y11(L));
    %^^^^^^
    %store this result in a vector and don't square your y here, because then it is not longer actually the area
    A_new(n)=trapz(x1(L),y11(L).^2);
    n_col=n_col+1;C=cmap(n_col,:);
    patch(X,Y,C,'DisplayName',sprintf('A(%d-%dHz)=%.2f',x_start,x_end,A))
end
A_new=A_new/sum(A_new);

[Rik]

You're welcome.

I doubt your legend is correct now, but that shouldn't be important for the rest of the calculation.

[Rik]

Because A is a vector once you put it in sprintf to create the legend entry. You should also change the frequency part from '%d' to something like '%.1f'.