0

I am trying to calculate the location of maximum curvature on many curves, and I notice that many times the curvature calculation is not maximum at the location I expect, whereas the rest of the time it seems to work.

For an example case in which this does not work as expected, I have included my curvature estimate calculation below (using diff to estimate first derivative and second derivative from given values describing a curve), plus the plots of each of the steps to determine max curvature. As you can see, the max curvature seems to occur way low on the curve, in contrast to expectation. Can someone please explain why this is and how to fix it?

first_deriv=(np.gradient(y))
second_deriv=(np.gradient(np.gradient(y)))
curvature=second_deriv/((1+first_deriv**2)**1.5)

fig, axs = plt.subplots(4,1,figsize=(6,10))
axs[0].plot(x,y)
axs[0].axvline(1000,linestyle='--',color='red')
axs[0].set_ylabel('Fit at \n Current')

axs[1].plot(x,first_deriv)
axs[1].axvline(1000,linestyle='--',color='red')
axs[1].set_ylabel('First \n Derivative')

axs[2].plot(x,second_deriv)
axs[2].axvline(1000,linestyle='--',color='red')
axs[2].set_ylabel('Second \n Derivative')

axs[3].plot(x,curvature)
axs[3].axvline(1000,linestyle='--',color='red')
axs[3].set_ylabel('Curvature')

Original curve vs calculations leading to curvature, with red vertical line placed where I would have identified maximum curvature by eye:

plots

lrm
  • 1
  • 1
  • if the steps in x are not uniform then you are not calculating the gradient properly, and should divide by the step size at each step. – Ahmed AEK Oct 22 '22 at 21:21
  • Steps in the x are uniform – lrm Oct 22 '22 at 21:31
  • Are you asking why the red line is slightly to the left of the peak in the third plot, or why the red line is far to the right of the peak in the fourth plot? – mkrieger1 Oct 22 '22 at 21:35
  • The magnitude of the first derivative (and then possibly also of the second derivative) seems way off. The original curve changes by less than 500 from x=1000 to x=1500, so the 1st derivative should be <1, not ~25. – mkrieger1 Oct 22 '22 at 21:49
  • @mkrieger1 I'm asking why there seems to be a peak and max magnitude of curvature around x=625 when there is no curvature to speak of in the original plot at 625. I would expect the max curvature to be around 1000 which is demarcated by the red line in all plots for reference. The shape of the first and second derivatives seem correct at least – lrm Oct 22 '22 at 22:41
  • It's because you are not calculating the derivatives properly. See my earlier comment (and also Ahmed's first comment). You need to take the x step sizes into account. – mkrieger1 Oct 22 '22 at 23:21
  • 1
    Does this answer your question? [Weird behavior of np.gradient with increased resolution](https://stackoverflow.com/questions/67971076/weird-behavior-of-np-gradient-with-increased-resolution) – mkrieger1 Oct 22 '22 at 23:24
  • @mkrieger1 Thanks! I doubted it was the x step size at first, since I was under the impression that it shouldn't matter given that it was evenly spaced and shared between first and second derivative derivations, and since I thought it was the shape that mattered rather than the amplitude for finding peaks. I have to admit I'm still confused as to why it matters, but I'll figure it out now I see that it works when I adjust for x spacing magnitude. I appreciate it! – lrm Oct 23 '22 at 04:49

0 Answers0