texts:
I am analyzing some program responsible for detecting the correct signal in the channel on the receiver side. It seems to me that such an implementation allows for the detection of the correct packets (because I mean packet transmission, from transmitter to receiver) from the received samples. It is known that under normal conditions the energy of the signal is greater than the energy of the noise, so the simplest algorithm to find the beginning of an incoming packet is to measure the energy of the received signal. When no packet is received, the received signal consists of only noise. When a packet begins to be received, the energy of the received signal is increased by the signal component, so the packet can be detected as a change in the energy level of the received signal.
So, the program I am analyzing looks like this (it's MathScript in LabView, it's almost identical to Matlab):

So, the example only divide the samples into groups with size “Bin Size”, and then take the sum of each
group, so that we can get the average power. The logic function can be interpreted as this:

%Keep the powerin the lower fractions of bins
low_lim = max(1,floor(0.1*bin_num));
%z = sort(power);
%noise_power = mean(z(1:low_lim));
noise_power = min(power);
signal_power = max(power);
%signal_power = mean(z(bin_size-low_lim:end));

%energy_threshold = max(.01, min(0.009, noise_power*1000));
energy_threshold = max([noise_power*100, (noise_power+signal_power)/2, 1e-6]);

noise_estimation_error = 0;
if 20*log(signal_power/noise_power) < 3
   noise_estimation_error = 1;
end

Now the one issue I don't fully understand is "Energy Threshold". It seems to me that the value 1E-6 was chosen based on brute force analysis (but correct me if this is some value that is often used to determine noise or something like that, because I have no idea), as some sort of reference point. Noise_power is just that. well, the noise power (although I don't know why multiplied by 100, it is possible that a higher order of magnitude allows a more fair detection - I don't know, sadly)
But the thing that bothers me the most is the value of (signal_power+noise_power)/2 - as far as I can see, it's a simple arithmetic average of signal power and noise power. Is this how energy threshold is calculated? What is it then responsible for? Is it just information about where the packet begins to be received and where the noise ends? The processed signal is later "chopped" of unnecessary noise in the receiver, so that only the proper signal containing the packet comes to the demodulation blocks (well, usually it is so, but we know how it works in practice). I've heard that such an implementation allows for greater efficiency, because it relieves the algorithms used in synchronization - they don't have to process portions of the signal, which e.g. are not the actual packets but just part of the noise. This speeds up processing of packets and simplifies their reception.
So is this the thing? Does the "Energy Threshold" carry information about the end of the noise and the beginning of the packet due to the calculated arithmetic mean? Or am I not fully understanding something here?


texts: amplitude compensation after filtering and mixing audio

texts:
I am working on a mixer which lets me mix dry and wet signal in a user defined ratio. Wet signal is basically just the all pass filtered signal. The mixing takes place in frequency domain and so does the filtering. So the requirement here is, after mixing I have to compensate for amplitude loss, i.e. the output signal should have same peak value over a frame as input signal.

So basically

output= x*wet +(1-x)*dry;

where dry is the input signal and wet is the output of filter, and all the signals are in frequency domain. x is a scalar ratio(0 to 1).

Is there any way to apply this compensation. I tried to multiply by ratio of RMS values of input and output
i.e. output = RMS(input)/RMS(output) *output (again everything in frequency domain)

but even this doesn't produce the desired result. (I am seeing meters in audacity after processing in aforementioned manner).


texts: Processing an ECG signal with a median filter

texts:
I have read in a couple of papers that the noise from an ECG signal can be removed via median filter. One such example I found was on stackoverflow, where multiple methods were suggested and one of them being the median filter. The following image is taken from the post on stackoverflow.

The main idea of the median filter is to run through the signal entry by entry, replacing
each entry with the median of neighboring entries. The pattern of neighbors is called the
"window", which slides, entry by entry, over the entire signal.

What I do not understand is, why doesn't the QRS complex get removed as well? When I tried it using matlab, the bigger I set the window, the better the signal got filtered, should't it be the opposite?


texts:
This is the shape of the bipolar amplitude function for a FIR filter with odd m:

I have 2 questions :

If i want to use this FIR as a low pass filter between 0 and $\frac{F_{c}}{2}$ , how can i set the value of $G(f)$ to 0 between $\frac{F_{c}}{2}$ and $F_{c}$ ?
My prof said that i can use the this FIR as a lowpass between 0 and $\frac{F_{c}}{2}$, but if i want to use the filter between $\frac{F_{c}}{2}$ and $F_{c}$ does it still behave as a lowpass ? (what does it mean that in that range $G(f)$ takes negative values ? )


texts: What is the benefit of adding zeros subcarriers in frequency-domain OFDM signal?

texts:
I have noticed in some OFDM systems, some null subcarriers are added in frequency domain. What's the benefit of those zero subcarriers, and is that mandatory ?
An example I have read about that adding of zero subcarriers is here This paper


texts: Ideal sampling of audio signals-Periodic spectrum

texts: When applying ideal sampling on audio singals, why is the spectrum periodic ?


texts: FFT analysis for series of different length

texts:
i have some data series in time domain which are in different length. i'm using FFT to transform them to frequency domain and now i want to calculate their average.
Since their lengths are not equal i'm doing:
(1) normalized the series before applying FFT
(2) appending zeros in the end of sequences to complete the max series length

however, once i'm applying inverse FFT after averaging the series, i'm getting discontinuous results.  the result up to the point of equal length is exactly the average, but from the point of the the difference it's completely out of bounds.

could you please advice what i'm doing wrong?


texts: How to segment multiple overlapping coins (ellipses)

texts:
I am doing a program that counts coins (and tells their value by their size) with OpenCV (C++).

Let's say we have an image of some coins.

If there is an overlap of two coins at the same time, it is relatively easy to segment the coins doing the distance transform, thresholding so that coins are sepparated, doing a skeleton of the background, and then drawing that skeleton into the original image.

But what if there is an overlap of 3 coins like in the image below:

What is the best way to detect 3 coins in that case?  Hough circles?


texts:
Give the equation
$$y[n]-y[n-1] +\frac14y[n-2]= x[n]$$
With initial conditions $$y[0]=0\\ y[1]= 0$$
Find the step response.

I have been able to get the solution $$y[n] =C_1(\frac12)^n +C_2(\frac12)^n + 1,,$$ I am have problem finding $C_1$ and $C_2$. How do I find $C_1$ and $C_2$?


texts: Fourier transform's sine and cosine count as N grows

texts:
I remember reading that in (discrete) Fourier transform for signals with even numbered N for length, the sine and cosine count is equal. Is this correct?

A bit of analysis:

N=1, there is only DC offset, which is a cosine wave of unlimited length.

N=2, Now, in addition to the DC offset term there is a other wave. But it's also a cosine wave, since trying to represent a sine with 2 points becomes impossible as only the 0 values of it can be sampled.

N=3, Sine can be represented as well as cos so there are 1 cos, 1 sine of same frequency and the DC offset cosine.

N=4, 3 cos, 1 sine, 3rd bin added with only cos wave possible.
...

So to me it seems that there are always going to be more cosine waves than sines. Did I make a mistake?

The interpretation that DC offset is a cosine could perhaps be debated though, has there ever been debate over this term? The fact that it surpasses the boundaries of N is a bit mind boggling.


texts: Is this an error in Oppenheim and Schafer's Discrete-Time Signal Processing?

texts:
In Discrete-Time Signal Processing by Alan V. Oppenheim and Ronald W. Schafer (3rd Ed.), in Figure 4.47 the input of D/A converter is $\hat{y}[n]$ but later in Figure 4.64 the input of D/A converter is $\hat{x}[n]$. Is this a mistake?
Normally, based on Figure 4.47 $\hat{y}[n]$ is the output of the discrete-time system with input $\hat{x}[n]$.


texts:
I am a new learner of OFDM. One source I refer explained OFDM in this way:

  Consider 4 symbols $X_0$, $X_1$, $X_2$, $X_3$ $\in \mathbb C$. They
  are modulated using Multi-Carrier Modulation. The time domain equation
  is given as

  $$c(t)=\sum_{n=0}^{3}X_ne^{j2\pi f_nt}$$

for the duration $0\leq t
> \leq 4T$, and $f_n=\dfrac{n}{4T}$.

Now, he samples the signal at a sampling rate, $F_s=\dfrac{1}{T}$ so that

$$c(mT)=\sum_{n=0}^{3}X_ne^{j2\pi \dfrac{n}{4T}mT}=\sum_{n=0}^{3}X_ne^{j \dfrac{2\pi nm}{4}}$$

This means that the sampled version of $c(t)$ is the IDFT of the symbols $X_i$, $i=0,1,2,3$.

What I am thinking is, since the spectrum of $c(t)$ is the shifted version of sinc pulses, considering only the main lobe of each sinc pulse, the spectrum is band limited to $\dfrac{-1}{4T}\leq f \leq \dfrac{4}{4T}$. So the sampling frequency should be twice the maximum which is $2*\dfrac{4}{4T}=\dfrac{2}{T}$. This is not taking me to the IDFT equation. Can you please correct where I am doing wrong?


texts: Complex Demodulation (CDM) and Fixed Frequency Complex Demodulation (FFCDM)

texts:
I do not quite understand the difference between the CDM and FFCDM.
How can I implement CDM in a program?


texts: should I apply low-pass filter when calculating central derivative?

texts:
Let us say that we have a discrete signal $I_n$, $n=0, 1, 2, ...$. According to Nyquist theorem the maximum frequency for such discretization is $f_{max} = 0.5$.

Now imagine that I want to calculate the derivative $D_n$ for this signal. The simplest approximation is the right hand side derivative

$$D_n = I_{n+1} - I_n$$

But what if I want to use a central derivative instead?

$$D_n = \frac{1}{2}(I_{n+1} - I_{n-1})$$

Just looking at this relation, it seems like I am increasing the sampling time from $1$ to $2$. Does this mean that before calculating central derivative I should apply low-pass filter to kill all the frequencies $f > 0.25$? Or am I being wrong and the simple fact that I am using the spacing = 2 for calculating derivative does not mean that I increase the sampling distance, because I can still calculate this derivative for each value of $n$?

No matter what the correct answer is (yes or no), can you please explain in more detail why? If the low-pass filter should be applied, can you explain if there is a method, how given discrete time relation I can infer to what degree I should smoothen the signal before using such relation?


texts: Proper rounding when converting between digitized values and floating point

texts:
I have some very large files that are unsigned 8 bit (0 to 255), real valued.  I need to do some simple tuning and decimation.  In the past, when I've dealt with low bit rate samples, I've run into issues with artifacts.  I want to make sure I get the rounding correct.

The data will be processed in Matlab using all native double precision.  The first question is should I subtract 127.5 or 128 when converting, or does it depend?

To keep things simple, I'm not doing anything fancy with the filter gain or multiple passes to make sure I don't clip.  Instead I'll notify the user if there was a problem.  The output also needs to be the same format, so I can't use any additional effective bits gained by decimation.  The next thing I'm unsure about is do I round 255.3 down to 255 or do I consider this an overflow condition?

EDIT

I thought I could just process the data with the offset in there doing something like this at the output.

y = y * scale + (1-scale)/2*255

The trouble is the tuning operation may filter out the offset before I correct for it.  So instead, I decided the simplest thing to do is subtract 127.5 and add it back in at the end.  Then round as I described above.

I'm not sure how the different methods of adjusting for the offset and scaling effect the resultant artifacts, but as arnfinn said there is probably nothing I can do about in any case (unless I'm willing to use dither).

EDIT 2

This is actually more along the lines of what I was thinking would be an issue.  In the distant past I have experienced this exact problem with requantization.  Certain combinations of quantization parameters can really mess up quality (worse artifacts).  My situation is a little different because I'm processing the data before requantization.  But I believe it can still be an issue.

A requantization algorithm for the transcoding of JPEG images

Signal Processing: Image Communication.
Volume 21, Issue 1, January 2006, Pages 13–21.
Jae Won Moona, Jong Seok Leeb, Nam Ik Cho.


texts:
I'm currently working on this problem from Oppenheimer's book:

Given $x_{1}(t)$ and $y_{1}(t)$ I should figure out $y_{2}(t)$ and $y_{3}(t)$ given that the system is LTI.
My progress so far:

I'm specifically asking about the graph for $y_{3}(t)$ because it looks like my solution doesn't agree with what's in the solution manual as for $0<t<1$ it's $y_{3}(t)=2$ instead of $y_{3}(t)=4t$ so am I missing anything?


texts: How to compare many bandpass filters (for their signal propogation delay)?

texts:
Let's say we have an audio signal sampled at 96 kHz, and we want to compare several bandpass filters to find the one with the lowest rectangular-envelope-propagation-delay (see graph below).

The band is interest is $[f_0-30\,\text{Hz}, f_0+30\,\text{Hz}]$, with $f_0 = 1000 \,\text{Hz}$, and we aim for an attenuation of, say $-50\,\text{dB}$ at $1300\,\text{Hz}$.

Is there a rule of thumb or a ready-to-use calculator (either as online website, or in Python) to find the best filter with the lowest propagation delay, for example in the case of a signal modulated by a rectangular envelope?

Of course I can run tests in Python, with many filters, and many variations of these parameters:

FIR filters (how many taps?)
Butterworth filters (which order, ...)
Elliptic, etc.

but it's a bit long, and a bit random.

Example: here with a Butterworth of order 2, we have a delay of about ~ 25ms for the propagation of the envelope (the x-axis is in ms indeed).

PS: is this more or less close to group delay?


texts: Amplitude not attenuated after bandpass filtering?

texts:
I want to make a bandpass filter for an IQ signal recorded from an Software Defined Radio with center frequency "Fc" & sampling frequency "Fs". For futhur analysis I wanted to focus upon signal received between a fixed frequency range & filter out rest frequencies.

I generated "sinc" funtion using Python Scipy firwin then I convoluted it with my signal using Python fftconvolve. But as you see in the plots amplitude of the frequency bins outside the pass band of my are still not attenuated.

Approach

Read the data,convert it from uint8 & strore it as complex notation.
Generate a "sinc" function of width half of difference of cutoff frequencies eg. for 20MHz to 40MHz width of filter will be 10MHz. Since "sinc" is symmetric it will filter frequencies from -10MHz to +10MHz.
Move the centre of filter to mean of cutoff frequencies using frequency shift property.
Convolve the filter with the data.

Relevant Code

data_raw=signal.flatten()
data_raw=data_raw-127.5

data = np.empty(data_raw.shape[0]//2, dtype=np.complex128)
data.real = data_raw[::2]
data.imag = data_raw[1::2]

nyq = 0.5 * fs

# Normalized frequency to be used for frequency shift
centre=(highcut+lowcut)/2-fc

width=np.ceil((highcut-lowcut)/2)
h = firwin(513,width,nyq=nyq,window='hanning')
h = np.append(h,np.zeros(1024-len(h)))

t=np.arange(len(h))
h=h*(np.exp(1j*2*np.pi*t*(centre/fs)))
filtered_data = fftconvolve(data, h)

Details:

Fcentre = 137.65 MHz
Fsample = 2 MHz
Lowcut = 137.65 MHz
Hightcut = 138.15 MHz

Plots:

Firwin filter linear plot

Firwin filter semilog plot

Output after applying filter


texts: How to warp a pixel patch to another camera

texts:
I have two cameras with rotation and translation starting from the world origin and a patch in the first camera's reference coordinate system. I want to warp this patch into the second camera's reference coordinate system. My idea is to transform every pixel like this:

$$T=T_2-T_1$$

$$R=R_1^{T}R_2$$

$$x' = \left[(x  K_{inv} R) + T\right]\;K$$

Where $R$ and $T$ are the transformation matrices between the cameras and $K$ is the affine transformation matrix that describes the camera properties.

Is this correct? I know that when you go from 3D to 2D you need to project, but that matrix is not invertible, so you cannot really do it here? Am I missing something really big?


texts: Detecting fundamental zero crossing point

texts:
I am looking for ideas from the more experienced members on algorithms for determining the fundamental zero crossing point of a 50 or 60 Hz signal. The fundamental zero crossing point is then used as the commencement point to calculate the RMS value of the original (voltage) signal AND other (current) signals.

Typically the frequency content of the Voltage is below 5% THD. At the moment I can't give u an exact figure on the SINAD but the typical Voltage signal will almost always be at 75% full scale (~250mV), or close to, of a delta sigma ADC (SINAD 110db) channel, with a sampling rate set to 6.4kHz. The target is an embedded DSP ARM M7 ~200Mhz with floating point hardware.

Ideally i would like to have 3-4 algorithms to prototype in Matlab and find the most appropriate. FIR vs. IIR vs. Bandpass vs. Zero-phase? Smarter approach? Would like to understand issues relating to sensitivity to noise, harmonics, group delay, edge effect, etc...

The IEC 61000-4-30:2015 standard on Power Quality states:
"r.m.s. voltage, refreshed each half-cycle
value of the r.m.s. voltage measured over 1 cycle, commencing at a fundamental zero crossing, and refreshed each half-cycle.

r.m.s. current, refreshed each half-cycle
value of the r.m.s. current measured over 1 cycle, commencing at a fundamental zero crossing on an associated voltage channel, and refreshed each half-cycle."

Thanks
Alex


texts:
Follow up to question and comments posted here

The SNR calculation method posted as answers were implemented in matlab

Implementation of Deve's answer:

% the starting and ending point of speech in samples obtain through Audacity;
% select a portion of the waveform where speech is seen to be present
begin_speech = 25400; end_speech = 26800;

snr_before = mean( data_s(begin_speech:end_speech) .^2) / mean( data_n(begin_speech:end_speech) .^2);
db_snr_before = 10*log10( snr_before );     % same thing, but in dB

% calculation on data after noise reduction follows:
residual_noise = data_s(begin_speech:end_speech) - output(begin_speech:end_speech);
snr_after = mean( data_s(begin_speech:end_speech) .^ 2)/mean( residual_noise .^ 2);
db_snr_after = 10*log10( snr_after );

The result is

db_snr_before = 15.6853
db_snr_after = -22.0388

Implementation of dspGuru's answer:

begin_speech = 25400; end_speech = 26800;

% before noise reduction
DB1 = 10*log10( var(data_n(begin_speech:end_speech)));
DB2 = 10*log10( var(data_s(begin_speech:end_speech)));
db_SNR_before = DB2 - DB1;

% after noise reduction
[data_stnr, fs_stnr, nbits_stnr] = wavread('AdMStry_thru_nr.wav');
Z = output(begin_speech:end_speech) - data_stnr(begin_speech:end_speech);
db_residual_noise = 10*log10( var(Z) );
db_speech = 10*log10( var(data_stnr(begin_speech:end_speech)));
db_SNR_after = db_speech - db_residual_noise;

The result is :

db_SNR_before = 15.6848
db_SNR_after = 0.9252

The begin_speech and end_speech represent those parts of the wavfile, in number of samples, where there is presence of speech for certain (obtained through audacity by observing the waveform directly)

Why is the SNR less after passing through noise reduction?

In the second method (dspGuru), to calculate SNR after noise reduction, to get the Signal value on the numerator, I have passed just the pure speech through the noise reduction algorithm as well. If I do the same thing in the first (Deve) method, the results of the two methods become almost equal

Why is this happening? Why is SNR after less than SNR_before ?

Here is a screenshot of the waveforms:


texts: Cross-correlation of time compressed 1D signal

texts:
If I ave 2 otherwise identical 1D signals that are phase shifted from each other then cross-correlation is a perfect way to identify the time lag between the 2 signals.

Now supposing I double the speed of one of those signals (eg sampled at every other sample).  The cross correlation will fail miserably at this point.

However, in image we can account for this kind of scaling (or even rotation) by converting the images to log polar coordinates and then performing the [phase] correlation.  This has the added advantage of being able to identify the angle of rotation and amount of scaling, I believe.

Is there a similar, simple transform for 1D signals to allow for identification of the time lag as well as identifying the amount of time compression of the signal?

I guess you could use log polar coordinates directly on the 1D signal but it would, instantly, square the memory requirements.  Any other thoughts?


texts: Carrier frequency offset and carrier leakage

texts:
Recently I study frequency Offset estimation due to mismatch in carrier frequency from the oscillator at the receiver. I understand how this cause problem can cause a phase shift. But I also saw that this cause problem called CARRIER LEAKAGE. So I look online and I saw the following:

Certain transmitter implementations may cause leakage of
the center frequency component. Such carrier Leakage can
occur in some transmitter because of DC offset. This issue
is seen on the receiver side as energy in the transmit center
frequency

Does carrier leakage mean an increase in interference power?? I thought that carrier offset happens at receiver and not transmitter. I also saw that carrier leakage can cause increase in the Error Vector Magnitude (EVM), I can't see the relationship between carrier leakage and EVM.

Looking forward for your discussion.


texts:
I am working in a project to record audio and process in real-time in C.For starting, I wanted to check how it works.
I downloaded the pa_stable_v190600_20161030.tgz from here.
I ran ,/configure and then ran make using WSL(Ubuntu) and then to test I ran ./pa_devs from the bin folder, but it returns 0 devices

I dont know what the problem is, for the moment I just want to get this working, is there any other method which I can try?


texts: Analyzing a particular discrete-time LTI system for input signal $x[n]=(1/3)^n$ for *all* $n$

texts:
I'm considering the following problem from some course notes.

  Suppose the following is known about a discrete-time LTI system:

  Given the input $x[n]=(1/3)^n$ for all $n$, the system produces output $y[n]=(2/3)(1/3)^n$ for all $n$.
  Given the input $x[n]=(1/2)^nu[n]$, where $u[n]$ is the unit step, the output is $y[n]=\delta[n]+a(1/4)^nu[n]$.

  The problem asks for the value of $a$ and for the response $y[n]$
  given the input $x[n]=(1/6)^nu[n]$.

I thought from (1) I could determine the system function by taking a Fourier or Z transform, and then use that with (2) to deduce the value of $a$. For example, the Z transform of (2) says $Y(z)=H(z)X(z)$, where $Y$ and $X$ are easily computed (or obtainable from a table). If I could deduce $H$ from (1), I could get $a$.

But, as discussed here, the signals $(1/3)^n$ and $(2/3)(1/3)^n$ for all $n$ do not have Fourier or Z transforms, so I'm stuck. I could imagine there being a typo: perhaps $(1/3)^n$ should be $(1/3)^{|n|}$. But if that's not the case, I'm not sure how to proceed.


texts: Where should triangle, sawtooth and square waves start?

texts:
I am reading a text on DSP, which shows generation of sin, cos, triangle, sawtooth and square waves.  Generated cos waves start at 1 and sin waves start at 0, which makes sense (cos(0)=1, sin(0)=0).  Generated triangle waves start at 1, square waves at -1, sawtooth at -1.  Wikipedia gives a number of constructions for square wave, which would lead to starts at 1, -1 or 0.

Is there a convention in the DPS/electrical engineering/oscillosope world on where triangle, sawtooth and square waves should start?


texts:
I have posted a question on Stackoverflow and was suggested to post here.

Please refer this question for the images: Similar Question:

I have used SURF Keypoints and a FLANN Based matcher and plotted the matching Keypoints as circles and then drawn the rectangle grouping all the keypoints.
As you can see, there are multiple vehicles matching the Template. I want to group the keypoints in a way such that I can get the count of the matching Vehicles.
How can I go about doing this?
I know that I can't just use a clustering method because, I need to group them more than just spatially.


texts: Inverse discrete time Fourier transform with differentiation

texts:
Consider a signal x[n] and its DTFT $X(e^{jω})$ . Assume $X(e^{jω})$ is differentiable. Compute the inverse DTFT of
$j\frac{dX(e^{jω})}{d\omega}$
You should write your answer in terms of $x[n]$ and elementary functions and constants, for example $π/2x[n]$ would be written :
$\pi / 2x[n]$


texts:
I'm having trouble on computing the autocorrelation of the sinc function
I want to compute
$$R_{hh}(\tau)=\int_{-\infty}^{\infty}\operatorname{sinc}(t) \ \operatorname{sinc}(t-\tau) \ \mathrm{d}t$$
where
$$\operatorname{sinc}(t) \triangleq \begin{cases}
\dfrac{\sin(\pi t)}{\pi t} \qquad & t \ne 0\\
\\
\quad \ \ 1 & t = 0
\end{cases}$$

Is the result another sinc function?


texts: explanation of correlation of stationary stochastic processes

texts:
I have some doubts about correlation in stationary stochastic processes.

I know that the expectation of a random variable is
$$E(x)=\int_{-\infty}^{+\infty} a f_x(a)da$$
and the mutual correlation of two variables is $E(x_1 x_2)$.

In two stationary stochastic processes the mutual correlation should be $E(x_1(t+\tau)x_2(t))$ now when and why I can translate this formula in:
$$
\int_{-\infty}^{+\infty} x_1(t+\tau)x_2(\tau)d\tau
$$
I have not understand where this last formula comes from.


texts:
I'm currently working on a project that requires me to classify sleep stages (Awake W, N1, N2, N3 and REM) based on only an EEG. Various algorithms and classifying standards (such as Rechtschaffen & Kales or AASM) actually require EMG and EOG in order to distinguish from REM sleep and wakefulness, which I don't have. Additionally, most algorithms use techniques such as Wavelet analysis, ICA, Neural Networks and such, but I have only learned up to Fourier Transform. Which means that I have to detect and classify these stages and their components (detect sleep spindles in order to classify N2, for example) based on only frequency analysis.

Power spectra analysis is poor since low frequency components tend to have high powers so the program ends up over-classifying 30s epochs as deep sleep stages.

I have a .txt file telling me exactly what's going on in every 30s epoch, so I have a priori knowledge regarding when is the subject awake or asleep (and which sleep stages the subject is in).

I appreciate any insight into this problem.


texts:
I tried to apply opencv OTSU binarization on the following image, but it generate a black image.

Here is the following python code:

import cv2

image = cv2.imread('test_otsu.png', 0)
ret, thresh = cv2.threshold(image, 0,250,cv2.THRESH_BINARY+cv2.THRESH_OTSU)
cv2.imwrite('temp.png', thresh)

Did I miss anything in my code? Please help me how to get correct binarization of the given image.


texts: Differential Equation of a Band Pass Filter

texts:
I'm struggling with a signal analysis electrical engineering college assignment. We need to a derive differential equation for a low-pass filter, high-pass filter and a band pass filter (made by connecting the output of a low pass filter to the input of a high pass filter). Here are diagrams for reference:

Low pass filter

High pass filter

$$$$
Band pass filter (buffer can be removed)

The equations I derived are:

low pass filter:
      $$V_o=V_i-RC \frac{dV_i}{dt}$$
high pass filter:
      $$V_o=RC \frac{dV_i}{dt}$$
band pass filter:
      $$V_o = R_2C_2 ( (1-C_1-R_1)- R_1C_1  \frac{d^2V_i}{dt^2} )$$

The equation for the band pass filter I found by making the input of the high pass filter the output of the low pass filter. So I derived the equation of the low pass filter with respect to time and got:
$$
\frac{dV_o}{dt} = \frac{dV_i}{dt} - R_1\frac{dVi}{dt} - C_1\frac{dV_i}{dt} - R_1C_1\frac{d^2V_i}{dt^2}
$$

I substituted in this rate of change of voltage into the equation of the high pass filter to get the equation I derived for the band pass filter.

Could anyone please tell me if what I have done is valid and if not, please explain to me what the correct strategy is. Thank you kindly.


texts: What algorithm can I use to "match" or "align" (translate and scale) two vectors of different size in C#?

texts:
Let's say I have two vectors A and B, of lengths, m, and n respectively. They are not "signals" but a collection of ordered time marks or indexes from another two "original" longer vectors (let's call them X and Y).
I try to compare A and B to "estimate" the similarity between X and Y. If they are aligned I can use a simple Jaccard index to do it. However, X and Y (and so A and B) can be delayed and/or steched with each other. I've assumed that if X and Y are similar, then B should be a linear function of A (B = c*A+d, c and d scalars). Besides, some elements could not match exactly or be missing in the other vector because A and B are of different lengths.
For instance:
A = [1 6 12 25 36 37 40], B = [3 13 25 73 82]
Should match the following pairs:
C = [[1,3],[6,13],[12,25],[36,73],[40,82]] -> B = 2*A+1
What algorithm or method can I use to find the gain c and offset d to get the best alignment between A and B? In other words, the parameters that maximize the amount of common elements between both these vectors.
Typically A,B are around 150 elements, and 0.8 < c < 1.2.

A little context: I'm trying to make a song dedupe application in C#. In fact, X and Y are audio firngerprints (chromaprints) made with AcousticID, which are basically binary matrixes (every row is a 32-bit integer). This questions could be directly solved to align X and Y but I think it could be harder due to the "noisy" nature of the prints (sample) and the larger size of vectors.


texts: Preamble detection using cross-correlation

texts:
I am generating a preamble with a sequence of [1 -1 1 -1 1 -1....] and spread with a barker-11 sequence. Now, the transmitted signal becomes [B(Barker), -B(Inverse of Barker), B, -B,....]. In channel, there is a random phase added to signal. In order to estimate the phase, First, I am doing the cross-correlation received signal with known barker sequence. When the peak value, I know it is barker sequence, I removed the symbol phases according to the barker sequence, now only random phase will present on the symbols.

However, I am also getting the peak value when it is correlated with '-B', then it inherently introducing 180 degree phase offset. How to avoid this problem??

Cross correlation I am doing by following:

if barker sequence [1 0 1 1 0 1 1 1 ....], and received symbols are complex,[r1 r2 r3 r4 r5 ....].

Corr_Real = real(r1) - real(r2) + real(r3) + real(r4) - real(r5) +......

Corr_Imag = imag(r1) - imag(r2) + imag(r3) + imag(r4) - imag(r5) +......

Correlation Index = Corr_Real^2 + Corr_Imag^2.

Please tell me if any other way to do.


texts: Probability of similarity for two random white noises

texts:
I'm looking for the method to calculate theoretical similarity of two random noises $x$ and $y$ with mean $0$ and standard deviation $1$, I've got also $12500$ samples of signal, sampled $2.5$ Ga/s.

Due to my idea I wanted calculate the probability that the two noises have coherence
$$C_{xy}(f)=\frac{|G(f)_{xy}|^2}{G(f)_{xx}G(f)_{yy}}$$
higher than $0.5$ in frequency range $2$-$20\mathrm{MHz}$? I wanted to compare it with experimental result. I found probably the answer to my question: Statistical significance of coherence values. How to do this for my signal? How $\alpha$ coefficient in this link was calculated?


texts: Find filter coefficients to model a device using its measurement

texts:
I am trying to realize a digital filter that has the same freq. response of an existing speaker. I have fed an audio sine sweep to the speaker and measured the speaker output, both at 48kHz. Then I perform FFTs of the input(X) and output(Y), divide the absolute values point-wise (absY./absX) to get the transfer function/freq. response (H).

Now I would like to determine the filter coefficients B and A for an FIR/IIR filter, so that I can model the speaker's response digitally. I understand my transfer function is in the domain of w, and filter coefficients are the coefficients of the difference equation of the filter.

This is where all my text-book theory seems to fall apart, and I was hoping someone could clarify:

If I understand correctly, since my data is discrete (sampled at 48k), my data is in terms of 'n'. If I performed an ifft of H (which is complex), the resulting vector is basically my filter coefficients?
Matlab suggests using invfreqz to find filter coefficients, given a complex frequency response. Does this function map from w to z domain? If so, am I right in understanding that the FFT needs to be converted to Z domain to derive its filter coefficients?
What parameters does a auto-regressive filter fit to? for eg: Yule_Walker, or Levinson, or if I want to run a gradient descent to fit a filter? What is the error calculated between?


texts: How to reconstruct an image from s-parameters in Matlab?

texts:
I have a set of data which is the s-parameter (s11), consisting of 2 columns ($1001 \times 2$). The first column represents the frequency and the second column is the magnitude in dB. I need to use these data to reconstruct an image using Matlab but I have no idea how to do it.

I know that I have to use ifft function first to transform from frequency domain to time domain but the result I had consists of real and imaginary part and I do not know how to deal with this.


texts: How to use matlab toolbox for 2qam dfe equalizer?

texts:
I made Channel equalization with matlab,
You can find the code below.

clc; clear all; N = 1e3; MQAM = 2;

d = randi([0 1],N,1);
x = qammod(d,MQAM,'gray','InputType','bit');
h = [0.1 0.5, 0.1 0.2];
r = awgn(complex(filter(h, 1, x), zeros(size(x))), 10, 'measured');

eq = comm.DecisionFeedbackEqualizer;
eq.Constellation = qammod(0:1, 2);
eq.StepSize = maxstep(eq,r);
eq.ReferenceTap = 2;

[y,e,w] = eq(r, x);

dec = qamdemod(y,MQAM,'gray','OutputType','bit');

errors = sum( dec~= d)

For some reasons it is not working and I don't know why the error is still the half of the transmitted symbols number?

I think it is the delay of the equalizer that I am not understanding ... somehow the outputs are delayed and I need to find the correct delay


texts: trying to understand fast fourier transform in matlab

texts:
Y = fft(X);
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
P1(2:end-1) = 2*P1(2:end-1);

1) I dont understand what lines 2-4 are acheiving, why must lines 2-4 be performed?

2) I realised sometimes input signal X must be windowed using a hamming window, what is the purpose of it?

currently an undergrad student, pls explain in a more basic form


texts:
How to find MOD-2 Circular convolution for the two sequences $h =[-1,3,-2,1]$ and $x = [1,-1,-2,1,3,2,1,2]$.

I know the answer is $7$ $0$ from matlab but I don't know how to find it graphicly or mathematically


texts: How to efficiently apply spectrum sensing algorithm

texts:
Lets assume 1 signal (lets call it PPG) that is composed "mostly" of 3
different main frequencies components ($f_1, f_2, f_3$) and their harmonics. Now let's assume we have 16 copies of this signal (let's call them $S_1,...S_n...S_{16}$), which are essentially the same (PPG) signal, but affected by some random delays (a fraction of smallest period present in  ), and "A LOT" of noise. These 16 signals, might themselves be also an aggregation of noisy-delayed copies of the original signal (PPG) itself.
Now, the version that we have now does the following:

It get's the addition/aggregation of this 16 noisy signal copies ($S_1,...S_n...S_{16} = S_{agg}$ )

We apply MUSIC (to $S_{agg}$) to generate the pseudospectrum, by creating an autocorrelation matrix made by artificially delayed signals.

Factorizing/diagonalizing it and assuming that only "main" 3 frequencies components are in the signal, while all the other frequencies are noise, we can extract and reconstruct the original signal (PPG).

From the Pseudospectrum, we get the main frequencies and recompose the signal based on them.

Now, this is very much similar to a typical communication scenario where the receiver can detect a "signal" which is essentially composed of and aggreagation of noisy randomly delayed copies of the original signal (PPG) due to e.g. multipath propagation.
However, in our problem, we DO have access to the 16 noisy delayed signals without the aggregation, so it does not make sense for us to solve this problem starting from the mixed signal.
We could apply MUSIC to each one of this 16 signals ($S_n$) and then find a way of aggregating the pseudospectrums (adding them or multiplying them, e.g.). This approach kind of works ok with the FFT (get FFT of each one of the signals, discard phase information, multiply or add the freq information and get a frequency representation of the real-signal).
I don't see why this would not work also with MUSIC. This makes our algorithm n times slower (n=16)
However, I am assuming that this is not the optimal approach, since I am only calculating the correlation of the individual versions of the signal with themselves, ignoring that I KNOW that they are essentially
the same signal and they should correlate with each other. To acknowledge for this, in practice, we could apply the same as single-signal-based MUSIC, but keeping Sn fixed and delaying Sn+1.
But....this would make the computation n-factorial (n!, with n=16, this is 2.1xE13) times slower, and we would also need to find something to do with the multiple versions of the frequencies extracted,
since they might not be "exactly" the same (addition or multiplication of the pseudospectrum might work here too)
I believe there should be some type of MUSIC-type algorithm that actually accounts for multiple received signals and it is able to extract the fundamental frequencies by computing e.g. a matrix of cross
correlations instead of autocorrelations. How should this matrix be composed? What would I need to delay? How to factorize/diagonalize the resulting matrix? and specially.... how to implement this in a simple
manner starting from a simple MUSIC implementation.


texts: If a note has a fundamental frequency of 100Hz, then its second harmonic is located at 200Hz and its tenth harmonic is located at 1000Hz. But for acoustic instruments, the harmonics will not occur at exact multiples of the fundamental. My question is, how much error is typical and what is the relationship between the degree of error and the harmonic number. So for example, is it common for the error in the tenth harmonic to be as much as, say, 50Hz. 1050Hz in other words. That's a 5% error. But would a 5% error be just as likely for the second harmonic?


texts: Amount of exactly approximated samples with Padé approximation

texts:
I've already had some good discussion with Fat32 in this question yesterday.

Today I'm confused again. The author of Statistical digital signal processing and modeling M.h. Hayes and Fat32 both stated that Padé can approximate p+q+1 samples. Then I found this exercise

and the solution

Do I count wrong or are there 5 samples used for the approximation? q+p+1 should be $2+3+1=6$ and therefor 6 perfectly approximated samples. Okay I've calculated it. Due to the periodicity the sixth sample is correct too. But what if the vector would be e.g. 10 at $x[5]$ (sixth sample)?

To solve for the a's I'd then have 3 equations for the vector $[1, a_1, a_2]^T$ which is impossible to solve. Where do I make the mistake?


texts:
I want to find the line between color blue and color gray in the following picture:

More exactly, the end points of that line.

Is the Hough transform a choice? If yes how should I set up the picture before using Hough?

Any way of solving this would be helpful.


texts: How to scale up the theoretical graph BER vs. Eb/N0

texts:
I would like to draw a graph presenting a correlation between BER (bit error rate) vs. distance for different bit rates for 16 QAM with bandwidth of 2.4 GHz.
$$
k=\log_2 M, \quad\text{where}\quad M=16
$$
I have $E_b/N_0$ theoretical graph based on the formula:
$$
\frac 1k \frac 32 \cdot \operatorname{erfc}\left(\sqrt{k \cdot 0.1 \cdot10^{\dfrac{{E_b/N_0}_{dB}}{10}}}\right)
$$

I know also that $E_b/N_0$ for $BER=10^{-4}$ equals $12.2\ \rm dB$.
$E_b$ is in ratio with bit rate and distance in $\sim \dfrac{1}{B_r d^2}$, where $B_r$ is the bit rate and $d$ the distance in meters.
I don't know how to scale up the theoretical graph changing only distance and bit rate (not knowing the theoretical distance and bit rate I have for $BER=10^{-4}$ and $E_b/N_0= 12.2\ \rm dB$).
What should be my approach?


texts: MATLAB: How to copy column vectors from a matrix into a cell array?

texts:
I have some sensor data stored in a 1000 x 5 matrix. I'd like to copy each column into a 5 x 1 cell array such that each cell contains a 1000 x 1 vector. I tried the num2cell(mymatrixname,1) function, but this created a cell array of 1 x 5, not 5 x 1. Using MATLAB, how can I do this in a compact way - without using a for loop?
Thanks,
kilosaw


texts: How to calculate biquad coefficients for the partial (1%, 10%) applying of pre-emphasis effect

texts:
I'm using sox to emulate RIAA pre-emphasis and de-emphasis effects on audio files. For pre-emphasis I call it like this:

sox --plot octave in.wav out_biquad.wav biquad 1 -1.700724 0.7029382 0.2380952 -0.1718791 -0.0442981

The biquad coefficients are taken from http://jiiteepee.tripod.com/de-pre-emphasis.html. I wonder if there is any easy way to calculate the biquad coefficients if I want to apply pre-emphasis only partially - in the range from 0% (no pre-emphasis) to 100% (full pre-emphasis).

The case for 0% seems to be easy, coefficients 1 0 0 1 0 0 would work just fine. But I'm not sure how to approach calculation of biquad coefficients for the combination of zero_filter + float_fraction * full_filter.


texts: Information Theory - units of channel capacity

texts:
In a first course in Information Theory, when the operational interpretation of channel capacity is introduced, it is said to be the highest data rate (in bits/channel-use) of reliable communication. While reading a few papers, I came across channel capacity being expressed in units of bits/s/Hz. So I was thinking about the connection between the two units and came up with the following explanation. Please let me know if this is wrong.

For a bandlimited channel (bandwidth = $W$ Hz), you can transmit at $2W$ symbols/sec by the Nyquist sampling theorem. So the rate "per bandwidth" (spectral efficiency) can be written as 2 symbols/sec/Hz. If each symbol is 1 bit, then you are transmitting 1 bit in each of the samples. So is 1 bit/channel-use equivalent to 2 bits/sec/Hz?

What is one "channel-use"?


texts:
I am a graduate student in applied math, specializing mostly in harmonic analysis. My undergraduate studies were mostly in pure math but I realized I prefer the applied/computational side. I was hoping to undertake a reasonable project which will lead to a short, but publishable paper in the area of signal processing. I have toyed around with a few ideas but I just cannot come up with something concrete as of yet.

Topics of interest:

Data compression
Distortion/noise
Optimization problems
Image processing


texts:
What is the difference between the weight vector (w) and the steering vector a(θ)?

I initially assumed they were the same thing, but the literature seems to be to be a bit fuzzy with this (or maybe its just me).


texts: Transfer function of resonant filter with 2 poles, peak at $f_0 = 500\text{ Hz}$, and $\Delta f = 32\text{ Hz}$

texts:
This a contest question. I'd like some help because I can't find any materials related to this topic.

https://www.qconcursos.com/questoes-de-concursos/questao/ecd6c966-51

My english translation:

Find the transfer function of a resonant filter with 2 poles, peak on $f_0 = 500\text{ Hz}$, $\Delta f = 32\text{ Hz}$, and sampling frequency $f_s = 10 \text{ kHz}$.

-> a) $H(z) = \frac{0.062}{1-1.8831z^{-1}+0.09801z^{-2}}$

b) $H(z) = \frac{0.081}{3-2.641z^{-3}+0.09801z^{-2}}$

c) $H(z) = \frac{0.082}{1-1.8831z^{-1}+0.09801z^{-2}}$

d) $H(z) = \frac{0.762}{1-2.8831z^{-1}+0.09801z^{-2}}$

e) $H(z) = \frac{0.262}{2-1.8831z^{-1}+0.09801z^{-2}}$

Any hints?


texts: What is the need for prediction filter in PCM and DPCM?

texts:
The DPCM works by the difference between actual samples and predicted samples.

If we already have the actual samples, we can quantize it and encode it later.

But why we use prediction filter?
And by the definition "predicts the present samples with previous samples" how does the filter do that?


texts: Data representation of largest DCT coefficients

texts:
In many examples of data compression, we take the S largest coefficients of a signal in a sparse basis (eg. 1% of largest DCT coefficients) and set the rest to zero.

How do we practically represent the compressed signal? In the above example, once the smallest coefficients have been set to zero, the signal is directly converted back into the time domain and compared to the original. What if we wanted to store the data?

In other words, how do we concisely store the indexes of the zeroed out coefficients?


texts: Finding set of orthogonal basis functions for composite signals

texts:
I've been given a list of 5 composite signals, where each is composed of 10 sinusoids of different frequencies. For instance, the first composite signal $S_1$ is given by
$$
S_1 = \sum_{i=1}^{10} A_i \sin (2 \pi f_i t- \phi_i)
$$
where $A_i, f_i$ and $\phi_i$ are known values. The other 4 composite signals are similar to $S_1$. In total, I have 50 different values for $A_i, f_i$ and $\phi_i$ (10 for each $S_i$).
My goal is to find the orthogonal basis functions for $\{S_1, S_2, \cdots, S_5\}$. I've decomposed them into $\sin$ and $\cos$:
$$
\sum_{i=1}^{10} A_i \sin (2 \pi f_i t- \phi_i) = \sum_{i=1}^{10} A_i [\sin(2\pi f_i t)\cos(\phi_i) - \cos(2\pi f_i t)\sin(\phi_i)]
$$
and from here I can see that $A_i, \cos(\phi_i), \sin(\phi_i)$ are constants, and a spanning set would be
$$
\bigcup_{k=1}^{50} \{\sin(2\pi f_k t), \cos(2\pi f_k t) \}
$$
so I end up with a 100-element set. I presume these are my orthogonal basis functions, since each $f_k$ is unique. Can I check if my train of thought is right?


texts: What is the $\mathcal Z$-transform of Bessel function $J_0(\alpha n)$ sequence

texts:
What is the $\mathcal Z$-transform of the sequence $J_0(\alpha n)$ for $n \in \mathbb{Z}$?

The Fourier transform of zero$^{\rm th}$ order Bessel function $J_0(\alpha x)$ is known to be $\frac{2}{\sqrt{\alpha^2 - \omega^2}}$ for $|\omega| < \alpha$.  This has a pole at $\omega = \alpha$. Does this imply that the $\mathcal Z$-transform will also have a pole on the unit circle?

EDIT:

The problem I'm looking at involves discrete samples of Bessel function i.e. $J_0(n)$. How should I proceed to determine its $\mathcal Z$-transform?


texts: Locally Adaptive Thresholding - How to Detect if Image Block is Uniform

texts:
I need to implement a locally adaptive threshold algorithm on a grayscale image in matlab. I use a sliding window over the image and calculate the threshold for the block using Otsu's method, which works well for blocks with non-uniform intensities.

However, when I apply this method to blocks with uniform intensities, I get very busy (junk) blocks. I want to set these blocks to either black or white in the result. I can decide whether to use black and white by comparing the global threshold value to the mean of the block (is this a good approach?).

My main question is: how to decide whether a block is uniform or not, or rather, how do I know which blocks to threshold and which to completely set to black or white?


texts:
I have series of $ N =133 $ data sampled symmetrically around the origin in an actual space at $ x(n \tau) = x_n $ with $ n \in \{ -(N-1)/2, \, ..., \, 0, \, ... , \, (N-1)/2 \} $. Now, I'm a bit unsure of how a typical FFT implementation performed on my data will be correct as a result of $ N $ being odd. My data will be zero padded to a length of $ (2N-1)_2 = 512 $, which is even. Ordering the sample from $ 0 $ to $ N $ as is done in an implementation,
$$
X_f = \frac{1}{N} \sum_{n=0}^{N-1} x_n e^{-2 \pi i f n / N}
$$
I know that $ x_0 = 0 $, and $ X_0 \stackrel{!}{=} 0 $ needs to be my DC term, as $ f $ is sampled at intervals of $ m / N \tau, \, m \in \{0, \, ..., \, N-1 \}$. But now my data has one more element on one side of my zero-datum than the other since all the spatial frequencies above $ N/2 = 256 $ are negative frequencies due to the periodicity $ X_{-m} = X_{N-m} $: there are $ 256 $ points to the right of zero, and $ 255 $ points to the left (as $ X_{257} \rightarrow -255 $ cycles per length, ..., $ X_{511} \rightarrow -1 $ cycles per length).
So my real data is sample-symmetric about the origin, but a zero-padded data set isn't. It's hard to explain exactly what I'm having difficulty with, but it seems like my transform is going to have some kind of "fencepost" error in it due to some kind of offset by one in my data about the center. There will be one extra zero-pad data point on the right side than the left. In fact, if I had taken an even number of real samples from the start, there is no datum actually associated with the origin at all, so shouldn't an $ n = 0 $ term zero out some poor piece of data that doesn't deserve it? I know my fear is unjustified, but I can't seem to work through whatever imaginary problem it is I'm having.


texts:
This is my first attempt at using GNU Radio.  I am trying to do this example, which is a sine wave generator.
The example works, so I am trying to connect the throttle block to the audio sink block.  I am getting an error that doesn't seem to match up with any of the options in the blocks.  The error is
   Source IO size "8" does not match sink IO size "4"

I googled around, but I didn't find anything related to the blocks I am using.  I also don't know what this error means, so I am hoping someone can match up the error to the parameters in the blocks.


texts: How to estimate the modulation transfer function of images?

texts:
I'd like some way to estimate the modulation transfer function in images.  This is the observed drop in amplitude vs increase in frequency of spatial details - nice illustration here:

(from https://www.imatest.com/docs/sharpness/)
In my case, I have a large collection of digital images of the same kind of thing; not literally the same object, but a fairly uniform set of objects with consistent overall properties (eg amount of spatial detail).  The imaging system is also exactly the same, but (because this is a fairly challening setup) the image quality varies considerably from image to image.  One of the ways in which it varies is that some parts of some images are seriously out of focus.  These are large images, around 10MP, and there is a lot of them (tens of thousands) so the statistical properties over the whole set can be determined with some confidence.  I think the typical image can resolve spatial detail on the scale of a couple of pixels (let's say 2-3 times the pixel pitch) but in some examples it looks like what it would look like if that area had a gaussian blur with radius ~10 pixels.  The images contain some sharp edges (features which are expected to be a fairly sharp edge in reality - they don't always look very sharp in the images) but no particular area is guaranteed to have an edge.
I'd like to be able to calculate an approximation to the real modulation transfer function just from the real images (or parts of images), without having something like a test chart in the image.  I think the consistency of the image statistics should make that possible in principle, but I'm not sure how to approach it.
How can this be done?


texts: How to obtain round peaks instead of triangular peaks when smoothing an input?

texts:
I have a signal (in blue) that I would like to smooth out (in red, moving average).

I would like to get red peaks to look more round than triangular and at the same time their width would become closer to the original ones.

Do you know / can you suggest a filter that would yield such result ?

Background:

The input is the zero-crossing rate from a song, for building a colored waveform as I've tried in this question : Coloring a waveform with spectral centroid or by other means

Here's the result of the new approach:

(unprocessed)

(ideal result, cheated somehow here using a post, expensive gaussian blur)

Update

Here's the result, using @Laurent Duval answer:

Also, I still need to try all of your suggestions again as my input/initial take was buggy, surprisingly now the ZCR yields better output than using FFTW out of the box (pics on the right):

Update 2

Simple moving average (green: 1-pass, gold: 2-passes)


texts: Ideal Discrete Time Hilbert Transformer with constant group delay

texts:
Let's assume we have an ideal discrete time Hilbert Trasnformation system (90-degree phase shifter) with a frequency response over one period:

$$H(e^{jw}) =
\left\{
 \begin{array}{ll}
  -j  & \mbox{if}\;\;\;\;\;\; 0 < w < \pi \\
  \;\;j & \mbox{if}\;\; -\pi < w < 0
 \end{array}
\right.$$

I want to define an ideal frequency response $H_d(e^{jw})$ of an ideal discrete time Hilbert Trasnformation system that has a non-zero constant group delay.

I consider that $$H_d(e^{jw})= [1-2u(w)]e^{j(\frac{\pi}{2}-τ\cdot w)}$$ since $$\left | H_d(e^{jw}) \right | = 1 \; \; \; \forall w$$ and
$$\measuredangle H_d(e^{jw})  =
\left\{
 \begin{array}{ll}
  \frac{\pi}{2}-τ\cdot w  & \mbox{if} \;-\pi < w < 0 \\
  -\frac{\pi}{2}-τ\cdot w & \mbox{if}\;\;\;\;\;\;  0 < w < \pi
 \end{array}
\right.$$

So, it seems an acceptable answer. However, I didn't make any calculations to find it. Just guessed it. So my question is: which is the approach to such a problem is order to define $H_d(e^{jw})$?
Thanks in advance!


texts: General form of dilations and translations of wavelet function

texts:
In some papers, the dilations and translations of a wavelet function is written as follows:

$$\psi_{j,k}=\frac {1}{\sqrt {2^j}} \psi\left({2^{-j}t-k}\right),\quad\text{where $j$ and $k$ are integers}$$

However, in some literature, the above formula is written in a different way:

$$\psi_{j,k}=\ {\sqrt {2^j}} \psi\left({2^{j}x-k}\right), \quad\text{where $j$ and $k$ are integers}$$

Are these two formula equivalent?


texts: When listening in to an AM signals of various frequencies, how do we exactly tune in?

texts:
Alright, so I read on AM and Double Sideband AM, but I don't the more fundamental idea - the frequencies.

Frequency is just: $ \frac{1}{T} $, where T is period

So difference between 5 Hz and 30 Hz is that period of 1st signal is $\frac{1}{5}$ = 0.2 seconds, and of 2nd signal = $\frac{1}{30} = 0.03 seconds$

and I know in AM you need to demodulate by multiplying by cos(wt) to bring message to lower frequency and then low pass filter to only have message signal etc

but for simplicity sake, let's assume two signals are transmitted, one is at 5Hz the other one is at 30Hz

How do you "listen in"? Because the frequency spectrum, where frequency is x-axis is a bit misleading!

It shows that you can magically shift to the right and just listen in to the needed signal! That's not the case! Frequency is just how fast a signal repeats.
Amplitude is changing according to sine function.

So if you have two signals, one at 5Hz and second one at 30Hz, then they both "are in the air", you can't just magically listen in every $\frac{1}{30Hz} = 0.03 seconds$ and only get signal at 30Hz, you'll also be getting the signal at 5Hz, they'll be overlapping because frequency is just how fast signal is changing.

So how do you listen in? Do you just listen in for everything?

Then you have a signal that's changing every "0.2 seconds (5Hz)" and another one every "0.03 seconds", they'll overlap of course. The only way they don't is if there's only one signal is in the air. I don't see how different frequencies allow to fine tune to one particular signal and ignore others.

Besides, even if there's just one signal in the air, say it's the 5Hz signal, how do you "listen in"? Do you just "listen in" continuously, or a receiver just turns on and shuts off every $\frac{1}{5Hz} = 0.2 sec$????

I watched videos, and they just say "receiver" receives, but I don't get how it receives...


texts: How to reconstruct the signal using overlapping frames in MATLAB

texts:
I understand that we can use the buffer function in matlab to split a long signal into overlapping frames. But is there any functions that we can use to put the frames back to reconstruct the signal?
Here's my current code in Matlab
n = 0.1:0.1:100;
y = sin(n*pi*0.5)';
Nw = 100;  % window length
Nh = 25;   % hop size
Y = buffer(y,Nw,Nw-Nh, 'nodelay');  % Y is the matrix containing the overlapping frames

What I want to do next is to use Y to reconstruct the original signal y. Something like this
yNew = unbuffer(Y, Nw, Nw-Nh)

Of course there is no such function as unbuffer, not that I've known of. So my question is, is there any function available to do the job of unbuffer ?
The situation is similar to this question here: https://www.reddit.com/r/matlab/comments/361y9w/reverse_buffer_function/, which was 5 years ago and still received no answers.


texts:
I'd like to do a fade out of a sound in a few seconds, but in a "spectral way", with a special curve :

(here x-axis is time, y-axis is frequency from a few hz to 20khz).

All the part "to be removed" has to be removed !

Do you know how to do this ? (an EQ / simple IIR/FIR filter  won't work for this, etc.)


texts: In Kalman filters why is it necessary to transform the systems dynamics matrix to the state transition matrix?

texts:
Previously when I have implemented Kalman filters I have used the transformation

$$
\mathbf{A(t)} = \mathcal{L}^{-1} \left( s \mathbf{I} - \mathbf{F} \right) ^{-1}
$$

to calculate the state transition matrix $\bf{A}$ from the system dynamics matrix $\bf{F}$.

Where the equation describing the state of the system is

$$
\vec{x_t} = \mathbf{F_t} \vec{x_{t-1}} + \mathbf{B_t} \vec{u_t} + \vec{w_t}
$$

and my prediction step calculation is

$$
\vec{x_{k|k-1}} = \mathbf{A} \vec{x_{k-1|k-1}} + B\vec{u_k}
$$

For example for simple sinusoidal Kalman filter started from the following equation of motion:

$ \ddot{x} = -\omega^2x $

So my matrix equation of motion was:

$$
\begin{bmatrix}
  \dot{x} \\
  \ddot{x} \\
\end{bmatrix}
=
\begin{bmatrix}
  0 & 1 \\
  -\omega^2 & 0 \\
\end{bmatrix}
.
\begin{bmatrix}
  x \\
  \dot{x} \\
\end{bmatrix}
$$

and my system dynamics matrix was

$$
F = \begin{bmatrix}
  0 & 1 \\
  -\omega^2 & 0 \\
\end{bmatrix}
$$

I then performed the transformation $
\mathbf{A(t)} = \mathcal{L}^{-1} \left( s \mathbf{I} - \mathbf{F} \right) ^{-1}
$ to get the state transition matrix like so:

$$
A(t) = \begin{bmatrix}
  \cos(\omega t) & \dfrac{\sin(\omega t)}{\omega} \\
  -\omega \sin(\omega t) & \cos(\omega t) \\
\end{bmatrix}
$$

Where I then replaced $t$ with my sample time $T_s$ where $T_s = \dfrac{1}{F_s}$ where $F_s$ was my sample frequency.

I then used this to extract my sinusoidal signal.

My question is did I need to perform this transformation? And if I did, why did I have to do so?

I was reading the paper linked here where they attempt the explain the Kalman filter in simple way and they did no such transformation.


texts:
I am trying to find the way to reduce the standard expression to compute the power of a generic sequence $x(n)$:

$$P_{\text{x}}= \lim\limits_{N \to \infty}\frac{1}{2N + 1}\sum\limits_{n=-N}^{N}|x(n)|^2\text,$$

when $x(n) = x(n + N_0)$ ($N_0$-periodic sequence), into the simplified formula

$$P_{\text{x}}= \frac{1}{N_0}\sum\limits_{n=0}^{N_0-1}|x(n)|^2\text.$$

How can I demonstrate that second one could be obtained from the first one, using the hypothesis of periodicity? I've tried something like this, but I get stuck in after $\text{(iii)}$:

$$\begin{align}
N &= N_0\tag{i}\\
P_{\text{x}}&= \lim\limits_{N_0 \to \infty}\frac{1}{2N_0 + 1}\sum\limits_{n=-N_0}^{N_0}|x(n)|^2\\
&= \lim\limits_{N_0 \to \infty}\frac{1}{2N_0 + 1}\sum\limits_{n=-N_0}^{N_0}|x(n + N_0)|^2\tag{ii}\\
&m = n + N_0\text{, so }\\
P_{\text{x}}&= \lim\limits_{N_0 \to \infty}\frac{1}{2N_0 + 1}\sum_{m=0}^{2N_0}|x(m)|^2\tag{iii}
\end{align}
$$


texts: Can a unit delay in discrete time be represented by exponential functions?

texts:
If we have a signal $y[n]$ and its unit delayed version $y[n-1]$, can we write $y[n-1]$ in terms of $y[n]$ times some exponential?

The reason I want to do this is to then take $y[n]$ common and simplify my analysis that I am doing for my research. The above is a generic form.


texts: How does a shift in time domain result in phase shift in frequency spectrum?

texts:
I am aware of the fact that a time shift of say $t_0$, results in a phase shift in the frequency spectrum. What confuses me is how this scales the rotational part of the transform by $t_0$ and doesn’t add a factor of $t_0$.

For example, $FT(x(t-t_0))=X(j\omega)e^{-j\omega t_0}$

I understand the magnitude isn’t changed, but how is this a phase shift? This results in $cos(\omega t_0)+jsin(\omega t_0)$, but to my understanding this is not a shift, but a scaling of the frequency.  A shift would be $cos(\omega+t_0)+jsin(\omega+t_0)$.

I probably have a misunderstanding with complex exponential, but the lack of info on this is keeping me up at night.


texts:
I would like to know what is the basic procedure to estimate the noise level of a CCD/CMOS sensor from an image/sequence of images acquired. I am going to use MATLAB to perform this operation.

What I did is I acquired a certain number of images with the shutter on, so the images are basically black images. The cameras I tested so far are all colour cameras, so there's a Bayer filter in front of the sensors. I have three channels R,G,B in the image output.

I was thinking that I can compute the average pixel value and the two-dimensional RMS to have an estimation of the noise. But I don't know if this is the recommended approach.

In addition, I guess that the noise estimation with the shutter on does not relate much with the noise level with the shutter on.

What do you think, do you think that the lens aperture also plays a significan role?


texts: DFT and perfect reconstruction of a square wave on a digital computer

texts:
I know that in theory, when reconstructing a square wave from its Fourier coefficients, unless we have an infinite amount of them, the resulting reconstruction will have Gibbs ringing artifacts due to lack of enough harmonics.
On a computer, we can take the Fourier transform X = fft(x) of a square wave x, and reconstruct it without artifact with x_rec = ifft(X), maybe with some rounding error of the order of 1e-17 or something but no visible ringing.
I don't have a satisfying answer for that? I guess there has to be something to do with the fact the "the square wave" x is a digitized version of a continuous wave, and my Fourier basis vectors (complex exponentials are also discretized of course since we are in a computer...) but still... how would you justify the absence of Gibbs ringing artefacts from Fourier reconstruction of the Fourier transform of a digital square wave ?
%%%%%%%%%%%%%%%%%%%% Tought experiment proposed by  Dan Szabo
fs=10;%sampling frequency
t=0:(1/fs):1-(1/fs);

s = [1 1 1 1 1 0 0 0 0 0];
sTr = imtranslate(s,[0.5 0])

sTr =

    0.5000    1.0000    1.0000    1.0000    1.0000    0.5000         0         0         0         0


texts: Identification of a given barcode as either 1 D or 2 D

texts:
I am working on an algorithm (in MATLAB) to identify/classify a given bar code (image) as 1 D or 2 D. I have considered using the contrast of THE IMAGE matrix as a parameter or using the FFT of the image.

I would be grateful for any ideas or suggestions on what you consider the best approach for this problem


texts: A system having impulse response $ h(t)=u(t) $ stable or not?

texts:
I know that for a system to be BIBO stable its impulse response must be absolutely integrable and the impulse response $ h(t)= u(t)$ integrates to approach infinity (i guess)
I proceeded as$$ \int_{-\infty}^{\infty}|u(t)|dt $$
$$ =\int_{0}^{\infty}dt $$
$$ =t\Big|_0^{\infty} $$
$$ =\infty$$
So the system must be unstable right? I m unsure  Can anyone explain in terms of poles of the $H(s)$?


texts: Jacobian of affine warp in Lucas-Kanade image alignment

texts:
I hate this kind of questions, but I'm really stuck. I'm trying to implement Lucas-Kanade algorithm as described in this paper (see pages 4 and 5). Unlike most explanations I've seen they don't assume optical flow, but instead use affine transformation for most examples and in general refer to warp function just as $W(x, p)$, where $x$ is the point and $p = (p_1, p_2, p_3, p_4, p_5, p_6)^T$ is parameter vector of affine transformation.

At the beginning of page 4 authors outline their version of Lucas-Kanade algorithm. I'm stuck at steps (4) and (5), namely, evaluating the Jacobian $\frac{\partial W}{\partial p}$ and calculating the steepest descent images $\nabla I\frac{\partial W}{\partial p}$.

(As far as I can understand, $\nabla I$ is just a pair of 2 matrices $(\frac{\partial I}{\partial x}, \frac{\partial I}{\partial y})$ in this context. Please, correct me if I'm wrong).

Authors carefully describe affine warp and even provide formula for its Jacobian (Equation 8):

$$\frac{\partial W}{\partial p} = \pmatrix{x & 0 & y & 0 & 1 & 0 \\ 0 & x & 0 & y & 0 & 1}$$

However, Jacobian of the warp is (not surprisingly) defined only for a single pixel and not the entire image. But in step (5) we calculate multiplication of image gradient and this Jacobian - $\nabla I \frac{\partial W}{\partial p}$, and as far as I can see from the context (see, for example, Figure 2 on page 5), it is done for the whole image and not per pixel.

My question is, how should I interpret this multiplication and what are the real sizes/formats of $\nabla I$ and $\frac{\partial W}{\partial p}$?

I understand that this question may require reading a lot from that paper, so I will be glad to explain any point that is not clear enough. Also feel free to refine question title or contents if you have an idea of a better wording.


texts: Derivative of equation containing Gaussian filtered image

texts:
Let $w$ be our image. For example, consider the following with the vectorized $w$:
$$ E(w) = \frac 12 \|Aw+b\|_2^2 $$
I know the optimal condition of the above equation:
$$ \nabla_wE = A^TA\hat w + b = 0 $$

Now, we consider the filtered image by Gaussian such that
$$ E(w) = \frac 12 \|A (G * w)+b\|_2^2 $$

In this case, how to get the optimal condition? (i.e., how to derivative w.r.t $w$?)


texts: How to find the difference equation directly from Direct Form II signal flow graph

texts:
I am trying to solve for the difference equation of the following signal flow graph:

I am aware that Direct Form II can be converted to Direct Form I, which finding the difference equation directly is much easier.
What I have come up with is defining an intermediate signal $v[n]$ as the top middle node. I get the following result,

$$\begin{align}
v[n] &= x[n] + \tfrac{3}{5}v[n-1] - \tfrac{38}{75}v[n-2] - \tfrac{2}{15}v[n-3] \\
y[n] & =v[n]-\tfrac{3}{10}v[n-1]+\tfrac{1}{3}v[n-2] \\
\end{align}$$

As you can see, $v[n]$ is recursively dependent on itself. I know I could transform the two difference equations into their system functions and multiply together, find the total system function, and then transform back to their total difference equation. But for practice, I would like to know if their is another way to do this. I would like to somehow find the total difference equation without using their system functions.

Any advice would be appreciated! Thank you!


texts: Compute phase coherence between two eeg signals

texts:
I have two EEG signals/channels and I would like to compute phase coherence in MATLAB. I have this paper: https://www.nature.com/npp/journal/v39/n5/full/npp2013330a.html but it does not say what Sij is so I don't know the formula.

Can anyone advise me on how to do this or does anyone have MATLAB code?

Many thanks.


texts: Welch Method FFT Python - Scaling factor?

texts:
I've been implementing a Welch method FFT and I am trying to work out the correct scaling factor that should be applied to the output of the function so the PSD is accurate because at the moment it's too low.
I am using scipy, signal.welch.
I know what the Power per resolution from a measurement of my signal and ive also calculated manually. I also know my manual calculation is correct because a straight forward FFT produces power around what I calculate but it fluctuates wildly across frequency (+ - 4 dB along what visually appears to be the average in the signal bandwidth).
So I needed a Welch method for a better spectral estimation. However this Welch method output it needs to be shifted vertically in dB. Meaning there must be a scaling factor applied.
Would anyone be able to say what it should be, perhaps the scaling factor varies with fftsize, percent of overlap, length of each segment etc?
Many thanks for your help


texts: Standard deviation of spectral flatness -- so what am I measuring, conceptually?

texts:
In my never-ending quest to identify snores, I've found that "spectral flatness" seems to be a fair measure of signal "quality".

I'm computing spectral flatness as the geometric mean of the FFT power $(R*2 + I*2)$ data points divided by the arithmetic mean of the same points.

I then (a little twist here) am computing the running (over 50 frames) arithmetic mean and standard deviation of the spectral flatness and computing a "normalized" standard deviation as the running standard deviation divided by the running mean.

For my samples I find that this metric is greater than about $0.2$ (ranging up to $0.5$ or so) when audio is "good" (ie, I have reliable tracking of the breathing/snoring sounds of a sleeping subject) and it generally slips down below $0.2$ when audio is "in the mud".  (I can improve on this discrimination somewhat by using a threshold that moves with other factors, but that's presumably a different topic.)  I also observe that the measure goes over $1.0$ when there's substantial background noise (eg, someone enters the room and rustles about).

So, my basic question is:  Is there a name (beyond "normalized standard deviation of spectral flatness") for what I'm measuring, and can anyone offer a conceptual explanation of what the metric "means"?

(I've tried a dozen other metrics for signal "quality", and this one seems to be the best to date.)

Added:  I probably should admit that I don't have a particularly good conceptual handle on what simple spectral flatness is measuring (just the Wikipedia article), so any further explanation of that would be appreciated.


texts: Streaming windowed sinc interpolation/resampling: trying to understand a Rust implementation

texts:
I'm working on a fork of the Rust dasp library, which is intended to be a DSP toolkit that abstracts over samples/frames/signals, and contains a number of functions and tools that work with these. One such bit of functionality is a signal interpolation feature, this supposedly allows for on-the-fly lazy resampling of a signal using sinc resampling. The Rust code that does the heavy sinc work is found here: https://github.com/RustAudio/dasp/blob/master/dasp_interpolate/src/sinc/mod.rs
However, I'm a bit confused. Looking at previous asked questions on this site and other resources on sinc interpolation and resampling, I'm seeing some yellow flags:

I can't seem to find any other "streaming" implementations of sinc resampling: they usually require a fixed-size buffer of samples as input, produce a fixed-size buffer of interpolated samples as output, and don't yield interpolated values one at a time.
Other sinc implementations are significantly more complex in terms of their required parameters. For example, the r8brain resampler requires a cutoff frequency and an attenuation level in the constructor (link), both of which are nowhere to be seen in dasp's resampler/interpolator.
I'm not seeing how the resampled result is lowpass filtered.

I'm admittedly new to the world of DSP, so I'm not sure if this code is correct, and if not, I'd like to understand why and how to fix it as part of my fork.


texts:
I'm trying to implement a normalized cross-correlation algorithm but I don't get what in fact is this measure. What confuses is the wikipedia definition:

$\frac{1}{n} \sum \frac{(f(x,y)- \overline{f})(t(x,y)- \overline{t}) }{\sigma_{f}\sigma_{t}}$

Which result is an scalar (AFAIK)

But then adds other way to measure it:

$\left \langle \frac{F}{\left \| F \right \|},\frac{T}{\left \| T \right \|} \right \rangle $

Where $F$ and $T$ are normalized vectors and $\left \langle . , . \right \rangle$ is the inner product. But the output will be a vector, isn't? Isn't supposed to give me a scalar as well? Am I getting something wrong?

The idea is to implement this formula and use it with matrices with same dimensions.

Thank you


texts:
I am from ECE background and trying to solve channel allocation problem. Let's assume I have three users and three available channels. I would like to allocate channel among them in such a way that average BER would be minimum.
For example, consider BER matrix for each user for each frequency allocation
\begin{bmatrix} User/frequency & User1 & User2 &  User3 \\ f1 &0.01 & 0.05 & 0.3\\f2 & 0.12 & 0.045 & 0.02 \\ f3 &0.02 & 0.01 & 0.01 \end{bmatrix}

Now, I would like to allocate channels f1--> User1, f2--> User3 and f3-->User2, which gives least average BER (Sum of each users BER/no of users). (Note: Same frequency can't allocated to different users.)

Brute force is obliviously one solution, but I am trying to find optimized solution.


texts: What is the physical interpretation of Lowpass/Highpass filtering?

texts: In Image/Signal processing context we have seen/we know that if there are say 4 samples and if we take an average of those 4 samples, then we say that the result sample is a low pass filtered output sample. This is very relevant in Image processing context. Now then, what is the operation which implies a high pass filtering, is it finding some max of the 4 samples and putting that as the output sample, or what other mathematical operations implies high pass filtering as averaging implies low pass filtering.


texts: I have two different audio short samples of an instrument playing. If one of the samples has a reverb effect applied, how easy is it to detect which sample has reverb? Is this possible?


texts: Ambient noise levels in desert and recording studio

texts: Are there mathematical models which suggest which of the two - desert or recording studio - is likely to be quieter? I would think a recording studio would be quieter than a desert. Is there a mathematical model which suggests this?


texts:
Are there any known methods to use existing keypoint detectors (which generally work on grayscale images) on color images? of course, one could always add the results of the detector on all three channels to get $3N$ keypoints, but I was wondering if there exist other established methods.

Specifically, I am looking to use a FAST like detector.

Note, that I am not looking for a color based descriptor, but rather for a color based detector.


texts: find impulse response from step response

texts:
if we have LTI system and we know unit step response of this system(we haven't original signal)
how we can calculate impulse response?


texts: ROC of Z transform Doesnt include a pole on the boundary?

texts:
I cannot figure out what is going on here. I have an example problem in my book that says the ROC of a certain function is $$ 0.5 < |z| < \infty$$
The function's denominator is $$ 1 - z^{-1} + 0.75z^{-2} -0.25z^{-3} + 0.0625z^{-4} $$
By properties of the ROC, there should be a pole on the boundary of the ROC. Since the example given was a right sided signal, the ROC will be the area outside of a circle boundary. So, the signal should have a pole on the boundary of that circle, which has a radius of 0.5 correct?
But when I plug in 0.5 to the denominator equation, I get 1, not zero, which makes 0.5 not a pole in this case. I have tried graphs, and placing different parenthesis, and calculated it multiple times, but the denominator equation above never reaches zero. According to ROC properties, it has to though!
Where am I thinking wrong here?
Thanks,
-Dom


texts: EEG signal processing with wavelet or fft?

texts:
I have confusion about signal processing related with EEG signal. I have done some of my research and that made me more confused about processing and filtering the signal.
Let me jump into the problem -

I have found that people are using FFT for EEG signal processing which I don't understand. Why would you do that? FFT is mainly made for stationary waves and we know that EEG signals are non-stationary waves hence FFT is not so useful for EEG signal processing. Well that can be solved with the SFFT but then again, this can't be done with real-time signal EEG. (I might be wrong, please help)

Wavelet transform is another way to process EEG signals because, first, it preserves time and freq whereas FFT looses time resolution. Also, wavelet can be applied on non-stationary signals. (I might be wrong here too, so please go ahead help me with this too)

As in the figure below - Neurosky shows real-time EEG power bands fluctuations, which I believe is FFT on raw signal and some kind of mathematical operation (maybe frequency averaging) applied on a range of frequencies gives those bands (alpha,beta,gamma,delta,theta), Am I right? if yes, then how can you apply FFT on real time signal? if no, then what is the best way to get those bands?

I have muse, emotiv, neurosky and openbci hardware which I borrowed from my community friends, I have been playing around with those hardware so I would do something with the raw signal so I started learning about these things but as deep as I go in that rabbit hole, more I get confuse. I was gonna apply ML eventually but I seriously got stuck on my first step.


texts: The figure is the FDMA,and i know that FDMA is that dividing the bandwidth into many blocks(frequency band),i want to ask Is there just a carrier in each frequency band in  FDMA?i mean one carrier for User1 ,one carrier for User2,and one carrier for User3,and so on.Or for the user N,there are lots of carrier in it,if i want to send a signal,then i have to use a carrier?


texts: JPEG compression steps after quantization

texts:
I have a 3-channel (for colours) a png image that I opened and

I splitted the image into 8x8 blocks
I applied all of the blocks discrete cosine transform
And then applied quantization
I stored the values in an array by zigzag traverse
I do not know what to do to reduce the size in this step

Now I do not really understand what i am supposed to after this steps. How did I compress this image? If I applied inverse of these steps and put the values in a new image, I believe the size will not change at all. I mean I am supposed to save the values in a txt and decode them with my application so that I actually made a JPEG compression?


texts: QPSK, from Square Root Raised Cosine to Rectangular Impulse Shaping

texts:
I have implemented a QPSK modulator and demodulator using classic text-book examples: Square Root Raised Cosine Impulse shaping and Matched Filtering.

The problem is that I've been told that I could only used square waves for the modulation instead of sine waves on the emitter side (only +/-Vmod) and this really degrades the performance of the communication.

After further investigation it seems that the root cause of the lack of performance is not really the use of square wave but the change of impulse shaping/matched filtering implied by it.

Basically my question is the following: Using only +/-Vmod (and 0V if necessary) on the emitter side, am I wrong in assuming that the best case scenario is using a Rectangular Impulse Shaping with the according Matched Filter or is there some other combination that could provide better results?

If you have any ideas or links I could learn from don't hesitate! I'm struggling because I don't know if my reasoning is wrong or if its my implementation.


texts: How to evaluate fixed-point implementation of LMS filter is correct?

texts:
I am having an LMS block with 6 filter coefficients.
The value of filter coefficients are
0.0001
0.00045
0.2535
0.546536
0.0000243
0.3423

I have tried to implement the LMS algorithm in floating-point with $\mu=0.01$.
LMS algorithm implemented in floating point is based on the below equations:(from figure) :
v(n)=wT(n)u(n),
e(n)=d(n)-v(n),
w(n+1)=w(n)+mu e(n)*u(n) ,which is the weight updation equation.
Here d(n) is the desired signal,u(n) is the input to lms block,v(n) is the lms filter output ,w(n) is the adaptive filter weights and e(n) is the error.
I have given a set of desired samples and input samples to the LMS filter and the LMS filter output data samples were analyzed. Now, to implement the same LMS filter in fixed-point representation, I will have to do scaling and rounding. Is scaling of filter coefficients by some factor $2^{s_1}$ enough and $\mu$ value be represented like $2^{-6}$ or should we consider any other factors while doing scaling operations. Finally, how to check the fixed-point implementation of LMS filter worked correctly?


texts: Discrete Fourier transform of a finite length signal which saturates at non-zero value

texts:
I am performing spectral analysis of a finite length signal that saturates to a non-zero value. The signal ($s(t)$) can, practically, be write as
$s(t) = f(t) \big(1-H(t-t_0)\big) $, where $t_0$ is the time at which the signal ends and $H(t)$ is the Heaviside step function. The sudden step from a non-zero value to zero at $t = t_0$ induces artifacts in the analysis.
I am looking for advice/literature on how to deal with signals of this type.
Question: How does one deal with artifacts cause by signals terminating at a non-zero value when
performing spectral analysis.


texts: What is Causing a Cyclic 1kHz Noise Signal in Audio Recording?

texts:
I dug my old Behringer UCA202 out of my shed and plugged it into a computer. I did some tests and it appears to be working fine... apart from some noise signals on the line input channels. I'm not sure where this type of noise could be coming from; weather it's the computer, the powered USB lead, or the device itself. What components could cause this type of interference?

###Spectrum Analysis and Statistics

###8 Second Spectrogram 750Hz - 5kHz (showing ~3.3Hz Period)

###1 Second 0-10Hz Low-Freq Noise Signal(related? i.e above ~3.3Hz cycle)

I have unplugged everything else from the computer (except power lead and VGA monitor).
I have tried plugging the device into a powered hub. The DC input on the hub states 5v 2.5A, but I'm using a 5v 2.0A power source. I don't know if this matters.
This is the recording:
https://soundcloud.com/n00dlz/uca202-testdevices-unplugged-monitor-off/s-H9zQsQzqArY
Only the monitor, mains lead and UCA202 are plugged into computer, and I unplugged the monitor from the mains while recording. There is nothing plugged into the UCA202.


texts: Does an audio equalizer consist of a perfect-reconstruction filterbank?

texts: Within the theory of perfect reconstruction, I discovered that it is not obvious that the combination of all filters in a filterbank gives you the original signal. Therefore, I was wondering whether an audio equalizer (either analog or digital) distorts the spectrum of a signal, even when the equalizer has 0 dB gain in each frequency band. Is this ever considered a problem?


texts: How to periodically estimate states of a LTI if the output is measured irregularly?

texts:
How can I periodically estimate the states of a discrete linear time-invariant system in the form $$\dot{\vec{x}}=\textbf{A}\vec{x}+\textbf{B}\vec{u}$$
$$\vec{y}=\textbf{C}\vec{x}+\textbf{D}\vec{u} $$if the measurements of its output $y$ are performed in irregular intervals? (suppose the input can always be measured).

My initial approach was to design a Luenberger observer using estimates $\hat{\textbf{A}}$, $\hat{\textbf{B}}$, $\hat{\textbf{C}}$ and $\hat{\textbf{D}}$ of the abovementioned matrices, and then update it periodically every $T_s$ seconds according the following rule:

  If there has been a measurement of $y$ since the last update: $$\dot{\hat{x}}=\hat{\textbf{A}}\hat{x}+\hat{\textbf{B}}\hat{u}+\textbf{L}(y_{measured}-\hat{\textbf{C}}\hat{x})$$
  If not:
  $$\dot{x}=\hat{\textbf{A}}\hat{x}+\hat{\textbf{B}}\hat{u}$$

(I have omitted the superscript arrows for clarity)

I believe that there may be a better way to do this, since I'm updating the observer using an outdated measurement of $y$ (which is outdated by $T_s$ seconds in the worst case).

Thank's in advance.


texts: Best way to reduce data to its symmetric and antisymmetric components

texts:
I have some noisy data that looks like this:

I want to find the components that are symmetric and antisymmetric with respect to magnetic field. The first thing I tried was just using:

Symmetric part $= \displaystyle \frac{1}{2}\left[y(x) + y(-x)\right]$

Antisymmetric part $= \displaystyle \frac{1}{2}\left[y(x) - y(-x)\right]$

Which gives this:

This isn't totally horrible since I was expecting the antisymmetric part to be linear with approximately this slope.

My question: Is there a better way to do this, possibly by using an FFT to break this into sine and cosine components?

(Edit: Fixed math)


texts: Instantaneous Frequency Estimation by Hilbert Transform - Theoretical Justification and Proof

texts:
I would like to better understand why the instantaneous frequency estimation by Hilbert transformation works (and especially why it doesn't work / lead to precise results in many cases).

The motivation is to estimate signal $x(t)$ by decomposing it into an amplitude envelope $m(t)$ and phase of cosine $\omega_c (t)$ (or carrier waveform):

$$x(t) = m(t)  \cos\left(\omega_c (t)\right)$$

Now, assume that $x(t)$ indeed is a result from such a process.

Questions:

1) There are two "parameters" to be estimated for any $t$, as such some constraints are needed. What are the constraints regarding $m(t)$ and $\omega_c (t)$ that are selected when applying the Hilbert transform decomposition?

2) Is there a proof available somewhere that given the constraints, the estimation indeed finds the correct amplitude envelope and carrier (for continuous and also discrete case)?


texts: Detecting edges oriented at certain angle

texts: How would I detect an edge oriented at 22.5° using a spatial filter?


texts: FFT of square wave: maths and actual numbers do not coincide

texts:
this fact is driving me insane... Theory states that the Fourier transform of a square wave is a sinc function, but if I compute the fft of a synthetic (perfect) square wave of amplitude 1 and frequency 1 and I check the numbers, I see them drift from a sinc function. Phases seem to change linearly thru the bins, instead of being all constant as expected, and real values (the odd ones of course) are all strangely equal to 2/framesize instead of being zero as I would expect. In fact, if I synthesize the same square wave in frequency domain by simply using a sinc function and I do the IFFT and I check the output, I notice that the square wave so obtained is almost perfect if it were not for some tiny ringings at the edges (Jibbs ?)
How is it ? What is then the "correct" function representing it in frequency domain if not a complex sinc as it should be ??? Or perhaps what I noticed is simply due to the fact that, since we are dealing with discrete quantities (fft =dft) things work slightly different than in the abstract case (fourier transform in pure math sense), and therefore some little modifications (scaling ?)are required ?  Thanks in advance


texts:
I'm looking at a data augmentation method for training up a neural network of speech data.  Currently I have two version of augmentation.
The first method works by taking an audio file and mixing in a set of noise files.
The second method works by convolving an audio file with another audio file to augment things like room responses and the like.
Once I have these combinations I then run them through my MFCC extractor and pass the MFCCs to the neural network training.
However this is very disk intensive with repeated loads of audio file.
What would be less disk intensive would be able to perform the combinations in memory.
Is there a way I can take 2 sets of MFCCs and combine them in such a way that I would get the same result as mixing or convolving the 2 base audio signals?
To put it another, way is it possible to convolve or mix signals in the MFCC space?


texts: If I have the following power spectral density of the quantization error can I tell that quantization noise is approximately white? Can you please provide some justification?


texts: Using Spectrogram to check whether aliasing will occur

texts: I was wondering if the spectrogram can be used as a MATLAB tool to check whether or not aliasing will occur when the signal is downsampled from let's say, 44.1kHz to 8kHz? And if yes, how do I go about it?


texts: How to smooth the echo duration change to get rid of audio clicks?

texts:
I implemented an echo effect using the following diagrams:

To summarize, an echo effect is basically a positive feedback effect that is implemented with a delay line zD which delays the input signal by D samples. My effects works fine. I get the echo effect.
My issue is I gave the user to adjust D sample amount during runtime just like in Pioneer equiments you can set echo as 4beats, 2beats, 1beats etc. When a user request a beat duration change I immediatly change the D parameter in my audio plugin. The problem is sometimes this causes click due to sudden audio level changes. I am trying to find a way to counter this but so far I can not think of anyways to do it. How may I ensure that the audio level is continous(doesnt spike and result in audio clicks) when user changes the parameter D?


texts:
I was reading pdf by caltech and in one of its section, Fourier transform of Unit step signal is calculated but I am confused, how this can be possible if region of convergence for Laplace transform ($1/s$) of unit step signal does not contain imaginary axis?

And if above case is possible then if it given that impulse response of a system is Unit step then frequency response of it should also exist and equal to          $H(ω)= πδ(ω) + 1/jω$ and then can we
calculate Fourier transform of output by computing $H(ω)X(ω)$ where $X(ω)$ is Fourier transform of input?


texts:
I need to design a stable IIR filter in MATLAB.  Specifically, I am trying to use a butterworth filter to correct the phase and gain of a 2 pole RC filter:

d = fdesign.arbmagnphase('n,f,h', n, f, h);
myFilter = design(d,'iir','butter');

The butterworth part is not really important, there are several options to choose from.  Is it possible to constrain MATLAB to produce an IIR filter that is stable?


texts: Characterising a system from its input and output

texts:
I have 2 data sets sampled at regular intervals which represent the input and the output of a system. A quick analysis of the plots shows the following:

There is a proportional factor between the input and the output.
There is a delay before a change in the input is observable in the output. This delay depends on the scale of the change: larger changes in the input are observable in the output after a longer delay than smaller changes in the input.

I would like to approximate the system with a function, but I am novice in signal processing and I have no idea how to get started. I don't even know if this is a solvable problem. Any idea? I am comfortable using Matlab/Octave, and computer programming in general.

I can get more data from the system if required, but I cannot control the input. This is because the system is only observable when a PID controller keeps it stable.


texts:
I have the next problem.

$H_{c1}(j\omega )$ is the ideal antialising filter and $H_{c2}(j\omega )$ is a real one. I'm asked to design $H(e^{j\Omega })$ so that $y[n]$ in the second diagram (the one on the right) is exactly the same as $x[n]$ in the first diagram (the one on the left).

I did something but I'm pretty sure that it isn't right. I found the equation for $H_{c2}(j\omega )$ as a function of $\omega$. Then I normalized the $\omega$-axis due to the A/D converter and divided the magnitude by $T$. If I didn't make any mistakes while doing the algebra, we get that

$$H_{c2}(e^{j\Omega })=\frac{9}{10}\cdot \frac{1}{\pi - \omega _{p}T}\cdot \Omega + \frac{9}{10}\cdot \frac{1}{\pi - \omega _{p}T}\cdot \pi + \frac{1}{10}\hspace{0.5cm} for\hspace{0.5cm} -\pi < \Omega < -\omega _{p}T$$

After that, I just thought of finding the inverse function of $H_{c2}(e^{j\Omega })$ for $-\pi < \Omega < -\omega _{p}T$... And that was it. Let $B(\Omega)=\frac{1}{\frac{9}{10}\cdot \frac{1}{\pi - \omega _{p}T}\cdot \Omega + \frac{9}{10}\cdot \frac{1}{\pi - \omega _{p}T}\cdot \pi + \frac{1}{10}}$. Then,

$$H(e^{j\Omega }) =
\left\{
 \begin{array}{l1}
  B(\Omega )  & \mbox{if } -\pi <\Omega <-\omega _{p}T \\
  1 & \mbox{if } -\omega _{p}T < \Omega<\omega _{p}T \\
B(-\Omega) & \mbox{if } \omega _{p}T < \Omega<\pi
 \end{array}
\right.$$

The problem with this is that I'm pretty sure that the filter I "designed" has no antitransform. I think that there is a smarter way of approaching this exercise but I just can't see it.


texts: Can we combine the application of SIFT with Kalman filter for 3D reconstruction of structural buildings?

texts: This is just a novice question from somebody who just got into this topic. Can the application of Scale-Invariant Feature Transform be combined with extended Kalman Filter?


texts: Building a Signal Model for Samples from a Sensor

texts:
I have a signal I'm getting from a sensor (shown below in the photo of an excel plot) and I need to process it to get it to look something like the red line I've overdrawn on it.  Sort of a moving average type process that doesn't fall too quickly unless the overall signal is actually falling.

Anyone have any ideas as to what the best method(algo) would be to approach this problem?

Any suggestions on directions I could look towards would be helpful as well, it's been years since I took signal processing in school and I don't remember what methods are even out there that could be tried.

Thanks!

Jonathan


texts: Frequency Analysis (DFT / FFT) of a Signal Without a Constant Sampling Frequency (Non Uniform Sampling in Time Domain)

texts:
I'm a stack exchange user for some time and now I'm registering to ask a simple question (I think!).

I have a vibration signal with an amplitude and time (sampling frequency not constant) in  a $10000\times 2$ double variable.

The data is available at: https://1drv.ms/x/s!AoCOij4si31tzgY89bhr6XH-_gQq

How can I do some sort of frequency analysis (FFT, DFT) or similar. How can I do it in Matlab?

Sorry if the question is duplicated but I couldn't find any answer for my problem.


texts: Power spectral density vs. FFT bin magnitude

texts: What's the difference between these? Both are measurements of some form of signal power, but surely there's some difference between the power they are measuring?


texts:
I am studying filters, in particular the Infinite Impulse Response ones.
I digged into theory and practice and I gain some decent knowledge on the topic. However, I can't really figure out or find good resurces about the so called "Initial conditions" or "State of the IIR".
So far I kind of understood that it is related to the fact that we are using a feedback on previous outputs (y(n)), but I couldn't find a clear and exhaustive explanation on the topic.
In the Digital Signal Processing, Proakis at page 91 of the book (107 if considering the PDF index) it says:

The term y(n_0 -1) is called the initial condition for the systems....

In Digital Signal Processing with Scipy at page 5, the author shows an example of a low pass IIR filter and shows how the initialization of the filter state changes the result of the output.
In the latter paper, the author uses 4 values as initial state, and I can't figure out why.
I'm having some problems in understanding the meaning of these initial conditions /state, where do they come from and how should I compute them. I know that Scipy and Matlab already provide functions to compute those states, but I would like to have a more deep knowledge and understanding of what's going on behind the "scenes".


texts: Confused about applying Routh Hurwitz to $s^2 +s + k$

texts:
I have a closed loop transfer function
$$G_{\rm cl}(s) = \frac{k}{s^{2} + s + k}$$
I am trying to find the critical gain of the system. From using the Routh Hurwitz criterion I get a $k = 0$.


texts:
Reportedly, coherence-based metrics are primarily influenced by linearity in phase. However, Klein et al report that coherence is also:

  ... very  sensitive  to fluctuations  of  linearity  in phase, relatively less so to nonlinear fluctuations of amplitude and completely insensitive to linear fluctuations in amplitude. How much both types of linearity contribute to the final measure, however, remains somewhat unclear and, therefore, some authors have argued that amplitude should not contribute to the final measure at all, using a "phase locking value"[9] instead,which is essentially the square root of (4) with all the moduli set to 1

without citing papers that prove the sensitivity of coherence with non-linear amplitude fluctuations.

Does this stand and if yes, are there any papers supporting it?