How to integrate between 2 points of a given array.

I want to integrate between two points on the array, but when I do it I get 0.

x2 = -0.4678
x1 = -0.4587

% find indices of these points
idx1 = find(field5(:,1) >= x1, 1);
idx2 = find(field5(:,1) >= x2, 1);
a1 = trapz(field_5(idx1:idx2,1),field_5(idx1:idx2,2))

what I want to do is to integrate between calculated points of the graph.

[John D'Errico]

Please do not answer your own question, JUST to make a comment to add information. Worse, you then accepted that answer to your question. And finally, your "answer" did not add any information, since you had that same information in your question.

Answer by @Vahram Voskerchyan moved to a comment:

Problem with the indices should be

x2 = -0.4678
x1 = -0.4587
% find indices of these points
idx1 = find(field_5(:,1) >= x2, 1);
idx2 = find(field_5(:,1) >= x1, 1);
a1 = trapz(field_5(idx1:idx2,1),field_5(idx1:idx2,2)

[Vahram Voskerchyan]

thank you, noted!

Answer 1

[John D'Errico]

idx1
idx1 =
2629

idx2
idx2 =
2622

So then what did you do? In order to understand (and then fix) the problem, you need to recognize what you are doing, and then to understand how MATLAB interptrets your actions.

a1 = trapz(field_5(idx1:idx2,1),field_5(idx1:idx2,2)

What is the vector

idx1:idx2
ans =
  1×0 empty double row vector

That is because idx1 is GREATER than idx2. So the colon operator recognizes this must be an empty vector.

So you used trapz to integrate an empty array. What should it return? Logically, trapz might decide the result must be zero.

Should MATLAB know that perhaps you meant something different than the code you wrote? Should it be able to read your mind, to then decide that it should do differently, effectively rewriting your code? I might bet that if it made the wrong choice here when MATLAB decided to rewrite your code with no warning, your thoughts might validly include an expletive.

So how did idx1 and idx2 arise? To understand that, we need first to consider that field5(:,1) is NOT a sorted vector.

  min(diff(field5(:,1)))
  ans =
  -3.7000
 
  max(diff(field5(:,1)))
  ans =
  0.0100

So merely finding the FIRST element of field5(:,1) that satisfies some property is not a good idea anyway.

Since I do not know what your real intention is here, I cannot give you a complete answer to your problem. Do you expect to integrate from x2 to x1, thus in a positive direction with increasing limits of integration? Are you hoping to integrate the elements from that variable from x2 to x1?

Finally, I might add that simply using trapz to perform that integration is asking for a terribly poor approximation, because if x1 and x2 are not in fact EXACTLY data points in that array, then there is part of the interval you are not using for integration. Far better MIGHT be to integrate a spline interpolant over that interval. I might do this:

  idx = find((field5(:,1) >= x2) & (field5(:,1) <= x1))
  idx =
  2622
  2623
  2624
  2625
  2626
  2627
  2628

Recognizing that this interval does not actually include x1 and x2, I might choose idx different. Thus, find the indicies from the largest point that does not exceed x2, then the smallest point that is not less than x1. And, since field5 is not a sorted array on the first column, better is first to sort it.

  [xs,tags] = sort(field5(:,1));
  ys = field5(tags,2);

This still suffers, because there are replicate elements in xs. While in the vicinity of x2 and x3, there are no reps, but even so, there are a LOT of replicates.

  numel(xs)
  ans =
  6026
  
  numel(unique(xs))
  ans =
  2163

Best might be to average the reps down. Thus, replace each replicate x with a single element, with an average value for y at those points. I've posted the consolidator function on the file exchange for this purpose.

  [xsc,ysc] = consolidator(xs,ys,@mean);
  >> numel(xsc)
  ans =
  2163  

Finally, now we can build a spline interpolant of that curve. Since the curve appears to be a highly spiky thing, best is to use pchip.

  spl = pchip(xsc,ysc);
  fnplt(spl,[x2,x1])
  hold on
  idx = find((xsc >= x2) & (xsc <= x1)); % plotting only points found in that interval
  plot(xsc(idx),ysc(idx),'ro')
  xline(x2,'g')
  xline(x1,'g')
  grid on

You should see the problem now. If we integrate the function over the interval given only the indexes you found, using trapz, your integral will miss a fairly significant portion of the interval. So the true integral over that interval seems to be:

  integral(@(x) ppval(spl,x),x2,x1)
  ans =
    0.0146

Yet, if we use ONLY the data from those points inside x2 and x1, we would get:

  trapz(xsc(idx),ysc(idx))
  ans =
    0.0117

So trapz was in error, underestimating the area by roughly 26%. Integral essentially agrees with trapz, if we restrict the interval. So the difference lies not in the use of a trapezoidal integration versus a spline, but in the use of the wrong interval.

  integral(@(x) ppval(spl,x),xsc(idx(1)),xsc(idx(end)))
  ans =
    0.0116

Find consolidator here:

https://www.mathworks.com/matlabcentral/fileexchange/8354-consolidator?s_tid=srchtitle