1

I implemented a function that gives the envelope curve of discrete values. I think there might be an error as when I tested it over a date I will provide you with at the bottom of the post, I get dissimilarities between the real data points ant the envelope curve as sown in this figure example of the code applied to a given data

from scipy.interpolate import interp1d
import numpy as np
import matplotlib.pyplot as plt

def enveloppe(s):
    u_x = [0,]        
    u_y = [s[0],]
    q_u = np.zeros(s.shape)
    for k in xrange(1,len(s)-1):
        if (np.sign(s[k]-s[k-1])==1) and (np.sign(s[k]-s[k+1])==1):
            u_x.append(k)
            u_y.append(s[k])
    u_x.append(len(s)-1)
    u_y.append(s[-1]) 
    u_p = interp1d(u_x,u_y, kind = 'cubic',bounds_error = False, fill_value=0.0)
    #Evaluate each model over the domain of (s)
    for k in xrange(0,len(s)):
        q_u[k] = u_p(k)
    return q_u   

fig, ax = plt.subplots()
ax.plot(S, '-o', label = 'magnitude')
ax.plot(envelope(S), '-o', label = 'enveloppe magnitude')
ax.legend()

Data S : array([  9.12348621e-11,   6.69568658e-10,   6.55973768e-09,
         1.26822485e-06,   4.50553316e-09,   5.06526113e-07,
         2.96728433e-09,   2.36088205e-07,   1.90802318e-09,
         1.15867354e-07,   1.18504790e-09,   5.72888034e-08,
         6.98672478e-10,   2.75361324e-08,   3.82391643e-10,
         1.25393143e-08,   1.96697343e-10,   5.96979943e-09,
         1.27009013e-10,   4.46365555e-09,   1.31769958e-10,
         4.42024233e-09,   1.42514400e-10,   4.17757107e-09,
         1.41640360e-10,   3.65170558e-09,   1.29784598e-10,
         2.99790514e-09,   1.11732461e-10])
skr
  • 2,146
  • 19
  • 22
Tassou
  • 439
  • 1
  • 5
  • 16
  • 1
    Why do you interpolate? Cubic spline interpolation is causing the effect that you observe. – Dima Chubarov Aug 08 '17 at 13:29
  • would you please suggest how to do it without interpolation ? – Tassou Aug 08 '17 at 13:37
  • What about linear interpolation, which does not overshoot? Or if this is in a signal processing background check out [analytic signal](https://en.wikipedia.org/wiki/Analytic_signal#Envelope_and_instantaneous_phase) and [hilbert transform](https://en.wikipedia.org/wiki/Hilbert_transform). – MB-F Aug 08 '17 at 13:44
  • could you provide me with what I can change in my code so there is no overshooting in the interpolation? – Tassou Aug 08 '17 at 14:04

1 Answers1

1

I would make two modifications to your enveloppe function to get a more monotone output

enter image description here

The idea is to avoid implicit addition of the left and right ends to the list of peaks that is used to construct the envelope

def enveloppe(s):
    u_x = [] # do not add 0
    u_y = []
    q_u = np.zeros(s.shape)
    for k in range(1,len(s)-1):
        if (np.sign(s[k]-s[k-1])==1) and (np.sign(s[k]-s[k+1])==1):
            u_x.append(k)
            u_y.append(s[k])
    print(u_x)
    u_p = interp1d(u_x,u_y, kind = 'cubic',
              bounds_error = False, 
              fill_value="extrapolate") # use fill_value="extrapolate"
    for k in range(0,len(s)):
        q_u[k] = u_p(k)
    return q_u   
Dima Chubarov
  • 16,199
  • 6
  • 40
  • 76