3
$\begingroup$

I have a signal containing a variation of digital ASK type on-off keyed sine wave among which some stochastic irrelevant harmonic waves occur, like you can see on this figure: Original signal

In its current form each pulse duration is the same, but it's not always the case.

I'd like to calculate each pulse duration. To do that, I was thinking of converting this signal to a square wave signal like this:

Expected signal

I'm not sure of the different steps needed for this conversion but I have some leads:

  1. Remove the noise => Low pass filter I presume?
  2. To have the signal between 0 and the amplitude, I multiply the signal by itself
  3. I apply a magic filter to convert the signal to square signal.

So my questions are:

  • Do I need a low pass filter to remove the noise?
  • What kind of filter do I need to convert the signal to square signal? Is it possible?

For information, I'm using Python and Numpy.

If there is any other method to calculate the duration of each pulse in the sequence, I'm open to all ideas.

Thanks!

Edit: Here's my signal:

   freq = 1000000
   fs = int(2*freq) # sample rate
   x = np.linspace(0, 10, fs)
   noise1 = 0.01*np.sin(2*np.pi * 1*x) + 0.05*np.sin(2*np.pi * 1.5*x) +     0.04*np.sin(2*np.pi * 0.2*x)
   noise2 = 0.1*np.sin(2*np.pi * 0.1*x) + 0.2*np.sin(2*np.pi * 1.8*x) + 0.1*np.sin(2*np.pi * 0.4*x)

  realSignal = 2*np.sin(2*np.pi * freq*x)
  sig =  noise1

  sig = np.append(sig, [realSignal, noise2, realSignal, noise1*noise2])
$\endgroup$
  • $\begingroup$ Can you provide your example signal ? $\endgroup$ – MaximGi Feb 26 '16 at 10:59
  • $\begingroup$ Ok, thanks. I made an answer, I'm also upvotting the question since it is interesting and shows good redaction effort from you. $\endgroup$ – MaximGi Feb 26 '16 at 14:46
1
$\begingroup$

What you can do is process this signal with a windowed RMS function (I'm writing in matlab code, but it can easily be translated in python) :

windowlength = 5; %you can change that, but keep it odd
hw = (windowlength-1)/2; %hw stands for "half window"
%pre-assigning the output RMS signal
sig_RMS = zeros(1,length(sig)); 
%I'm adding zeros the the signal's left and right for boundary calculation
sig_zero_padded = [zeros(hw,1), sig, zeros(hw,1)];

sig_zero_padded = sig_zero_padded.^2; %squaring the samples
for i=(hw+1):(length(sig_zero_padded)-hw)
   % calculating square-mean inside the window and assigning the result
   %to the appropriate sample in sig_RMS
  sig_RMS(i-hw) = sqrt(mean(sig_zero_padded(i-hw:i+hw)));
end

Then you set a minimum value for which superior values of the RMS signal will be assigned to 1, and the other to zero :

trigger = 0.5; %I chose arbitrarily
square_sig = (sig_RMS > trigger);

And you should get your square signal.

What we are doing here, in simple language, is taking little bits of the signal and calculating the RMS value of each bit. If this calculated value is significant enough to be considered as "not noise", we assign '1' to it.

EDIT : I managed to compute it and it worked. Note that I had to take Fs = 10^5 because my computer does not have enough capacity for Fs = 10^6.

Here is the signal (blue), and the calculated RMS signal (red) for windowlength = 8191 :

enter image description here

Then I set trigger = 0.4 and obtained two crystal clear zones identifying the useful signal :

enter image description here

So in conclusion :

  1. Set windowlength to get the RMS signal mean value far enough of the noise while keeping accuracy. In my case it was approximately 1 % of the total number of samples.
  2. Set trigger between the RMS crests and the noise

However, I'm having that feeling that I'm reinventing something so feel free to look up for some literature about this method, I'm sure someone did the math to obtain optimal values of windowslength and trigger. ;)

$\endgroup$
  • $\begingroup$ Thanks so much for your answer. I try to apply it to my problem but I still have some problem (I have many square waves instead of two "big"). I think I'm going to play with the windowlength value and the trigger a little more. $\endgroup$ – Athius Feb 26 '16 at 16:06
  • $\begingroup$ You're welcome. If the answer suits you, don't forget to accept it. If you need a few more information, you can also ask here. $\endgroup$ – MaximGi Feb 26 '16 at 16:08
  • $\begingroup$ Sorry, i didn't see your edited comment. If it detects too many square wave, its probably because the windows is too small. You can also do a low-pass filtering (smoothing) operation after the RMS calculation and keep and small window $\endgroup$ – MaximGi Feb 27 '16 at 17:28
  • $\begingroup$ @Athius see edit, I managed to make it work. $\endgroup$ – MaximGi Feb 27 '16 at 19:24
  • $\begingroup$ It works great! Thanks so much!! I just put the max value instead of True if we have the signal or not. I'm going to do more tests with different signals. I don't know if you reinventing something, but it's working pretty well and pretty fast too! $\endgroup$ – Athius Feb 29 '16 at 9:45
1
$\begingroup$

Create the envelope of the signal (analytic signal) using the Hilbert transform and then threshold.

import scipy.signal.signaltools as sigtool

env = numpy.abs(sigtool.hilbert(sig))
threshold = 1
square_sig = (env > threshold)

(I haven't got python installed ATM so I can't test this)

$\endgroup$
  • $\begingroup$ Thank you for your answer. But unfortunately I don't have scipy and I can't build it with pip (some errors and I don't have time (and take time) to fix those compilation errors) $\endgroup$ – Athius Feb 29 '16 at 9:47
  • $\begingroup$ You should install anaconda or some other python distro for scicomp. In particular, numpy and scipy will already be there. $\endgroup$ – dohmatob Feb 29 '16 at 18:03
  • $\begingroup$ Do you have FFT? Take the FFT, set half the spectrum to 0, then IFFT back. It is the same as hilbert above. $\endgroup$ – geometrikal Mar 1 '16 at 10:55
  • $\begingroup$ I just succeed to install scipy and hilbert method is just awesome! I also implement the FFT and IFFT method and it's worked great too. Thanks a lot guys! $\endgroup$ – Athius Mar 11 '16 at 17:25
  • $\begingroup$ @Athius no problem, glad it worked. $\endgroup$ – geometrikal Mar 12 '16 at 1:40

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.