texts: ORB Algorithm: What Is Intensity Centroid?

texts:
So, I am reading a paper about ORB: An Efficient Alternative to SIFT or SURF and it happens to calculate the direction of  corner from it center . According to Rosin the author of the paper;moments of a patch is given by

and the centroid is of the patch from the center is calculated using:

The direction of the patch is computed from

My question is what are these value of p and q and why is it out there?


texts: Can the numerical variance be lower than Cramer-Rao Lower Bound?

texts:
For a sample mean unbiased estimator, which achieves the CRLB, I plotted the variance and compared it to the theoretical CRLB. The plot below shows what I got. I have small number of samples to calculate this variance.

The black line is the theoretical bound and the red line is the variance calculated using 10 samples.
Why is the calculation lower than CRLB? If CRLB gives the lowest possible value, should it not mean that the red line should always be on or above the black one, even considering low sample size and numerical errors? An insightful answer will be welcome.


texts:
I have two questions:

If the FM radio at the receiver mixes with a another frequency to bring down the carrier frequency, wouldn't that affect the message signal?
In today's digital world, why the FM radio still uses frequency modulation (FM) and not FSK?


texts:
I am currently trying to evaluate a high resolution ADC for which I need to window my output data for the FFT. Up until now I have been using the default 4 term Blackman harris window in MATLAB. However, I came across a few app notes that recommend using a seven term Blackman Harris for the SNR ranges that I am interested in.
I have tried looking quite a bit but I cannot seem to find an equation for implementing the seven term BH window, neither is it available in MATLAB. Can someone point me to a source where I can find this?
As far as I am aware Harris' classic paper on windows doesn't provide the equation either.


texts: Fourier Transform of a signal using direct integration and properties

texts:
Am trying to compute the Fourier Transform of a function using the properties of the Fourier Transform once and checking my answer using direct integration. My problem is that am not getting the same for both ways.

Sorry for not being able to type the equations as I am not familiar on how to do that, It would be great if somebody can tell me how. Your help is appreciated. Thanks.


texts:
This is a question I found on internet under DSP section so that's why I am posting it here. Help me understand it please.

  A computer adds 1000 random numbers that have each been rounded off to the nearest 10$^{\rm th}$. Find the probability that the total round-off error for the sum is $\ge 1$.

How probability of any system can be greater than 1?


texts: Effect of cyclic prefix and zero-padding in OFDM

texts:
What happens if zeros are added instead of cyclic prefix in OFDM system?
How does it affect performance of system?


texts: Windowing after filtering a discrete signal to get it back to its original length

texts: What happens if you filter a discrete bandpass signal of length $N$ with a bandpass filter (Hann truncated sinc) of length $L$ and apply a window of length $N$ to the filtered signal of length $N+L-1$ to get back the signal to its original length $N$?


texts: How to validate that an audio algorithm result is independent of microphone device?

texts:
I have a machine learning algorithm that takes speech sample audio recordings collected from mechanical Turk. During processing it was shown that some audio from certain OS/microphone devices have been preprocessed by noise suppression algorithm by default which affects the results I am getting.
As part of my writeup I need to ensure that OS and Microphone preprocessing is accounted for. Ideally I would like my algorithm to be independent of the recording source, however I am not even sure how to go about detecting the kinds of preprocessing that is being applied and give a lower and upper error bounds for how it affect my algorithm.
My question is how does other people deal with preprocessing on audio files that are outside your control? Are there some common preprocessing that I should know about (echo cancellation, noise suppression, etc.)?
Something I can do for example is to add noise suppression myself on some reference files and add +/- for the max and minimum deviation from the results from the reference file.  However it's a bit unsatisfactory as I do not know how different my preprocessing is to other devices.


texts:
I have a matrix of 10 column, each column has measurements at different times. I have a second matrix also of 10 columns, each consisting of measurements at different times.
I want to carry out cross correlations of every combination within this matrix. Example, xcorr row 1 matrix 1 with row 1 of matrix 2, then row 1 of matrix 1 of row 2 of matrix 2 and so on. I would like correlations and time lag for each and every combination of columns.
I currently do this, but it means I have to do each row separately which is very time consuming.

[xc,lags] = xcorr(matrix1(:,1),matrix(:,2));
[~,I] = max(abs(xc));figure;
plot(lags,xc);
legend(sprintf('Time delay %d',lags(I)));
title('Column 1 matrix 1 xcorr with Column 2 matrix 2')


texts: Is there a FFT algorithm with the circular buffering?

texts:
I have an audio signal and i want to process it sequentially. Is there an algorithm which allows to modify the spectrum of a signal when i remove a part of signal and add the other one?
For example, i obtain a set of 256 samples of input data. Then i do 256 FFT and i want to add the result spectrum to 16384 FFT which is the big buffer. So i have the  big FFT buffer with the good frequency resolution and i should sometimes modify to correspond the spectrum as if i do 16384 FFT.
The main objective is i need to avoid 16384 FFT because it's quite expensive. The way i want to optimize is (in time domain):

Shift the buffer to 256 samples;
Place the new samples to last 256 entries in buffer;
Make 16384 FFT.

So whenever the 256 sapmles come i need to recompute 16384 FFT whereas i've already have 16384 - 256 computed samples. Is there way perform such task (circular buffering FFT)?
My thoughs was (in frequency domain):

Multiply the 16384 FFT buffer by the complex sinusoid (that coresponds to circular shift the signal in time domain to 256 samples);
Calculate 256 FFT of new portion of signal;
Expand the spectrum to 16384 by inserting zeroes between bins (it's corresponds duplication of time domain signal);
Mix the signal with removing duplication of the expanded spectrum (by convolving with sinc function), and zero last 256 samples (by convolving with sinc function).

Is there way to avoid two comlementar convolving in the last step?


texts: Box ﬁlter approximations for oriented gabor ﬁlters

texts:
I'm currently working on a library to generate synthetic fingerprints using the SFinGe method (by Maltoni, Maio and Cappelli) link :http://biolab.csr.unibo.it/research.asp?organize=Activities&select=&selObj=12&pathSubj=111%7C%7C12&

One of the steps requires me to apply different gabor filters to an image, each pixel in the image have an orientation and a frequency associated, so the convolution is not done with one kernel over the entire image but the filter must change during the process depending on those attributes of the pixels, that way each pixel on the image is altered in a different way.

If you apply the filters this way, and convolve the image several times(you also have to binarize the image after each convolution) you obtain this:

A master fingerprint, I have currently implemented the standard gabor filter using 2d convolution. I would like it to be faster, Im trying to separate the filter but the problems with separable ﬁltering is convolution in the complex domain. Im working on C++ usin openCV and Im not sure how to perform such convolution.

Now reading this paper from a different answer, I saw an additional method, box ﬁlter approximation:

Now from reading the paper is not clear to me how to approximate the gabor filter using this method, what is B(x, y)?

Im not very knowledgeable in image processing (or mathematics). I would appreciate it a lot if someone could explain to me how this approximation works :).

Also an important note, can this be done in a contextual matter? The way in which I have to apply the filter in order to generate the master fingerprint is that I have to witch the kernel in each step of the convolution according to some local properties of the current pixel being convolved.

Thanks for reading!


texts: How to demodulation BPSK in real case(through a channel)?

texts:
I am configuring a channel sounding system. Using PRBS to correlate for time delay. Some details as follow.

Generate PRBS.
Using BPSK to modulate PRBS and Tx transmit.
Rx receive signal and demodulation.
Correlation the signal that demodulation from BPSK and original PRBS.

But now, I am stock at step 3. The following block diagram is the system that I set and the signal representation of each node.

I would get the wrong demodulated results that if (delta1) minus (k1R1 or k1R2 or k1R3) >= 90deg .
How to demodulation BPSK in real case with multipath??


texts: Does a collection of Gaussian random variables necessarily constitute a Gaussian Process?

texts:
If $\{X(t)\}$ is a Gaussian Process then the random variables $X(t_k)$ where $k = 1,2,3...n$, are jointly Gaussian.

If each random variable $X(t)$ is a Gaussian variable, then will the random process $\{X(t)\}$ be a Gaussian Process?

I am asking because while studying the Statistical proof of why we model noise in communication as Gaussian Process $N(t)$, I came across the fact that noise $N(t)$ at time $t$ is basically a sum of micro phenomenon which makes the random variable $N(t)$ a Gaussian random variable. Now that random variable $N(t)$ is Gaussian how do we prove that the random process formed by these random variables is a Gaussian Process?


texts:
I have been studying the following code but still unable to explain its output. I really hope someone could shed some light.

Let me begin with the baseline. Suppose I have a 20 seconds sinusoid with frequency components of [1,2,3,4....500Hz]. The signal has Fs=1024Hz, so it's FFT output looks as follows. The FFT amplitude is also scaled to preserve Parseval's Theorem. The frequency resolution is Fs/nfft = 1.0 Hz per bin in this case, as the result the FFT output smeared itself and make the output looks like a flat line, which is expected.

# Python code

import numpy as np
import matplotlib.pyplot as plt

dur = 20
fs = 1024
t1 = np.arange(0, dur, 1/fs)
x1 = np.zeros_like(t1)
for f in np.arange(1, 500):
    x1 += np.sin(2*np.pi*f*t1)

plt.figure(figsize=(8,4))
plt.title('Signal')
plt.plot(t1[:100], x1[:100])
plt.xlabel('sec')
plt.grid()
plt.show()

nfft = 1024
df = fs / nfft
X1 = np.fft.fft(x1, nfft)[:nfft//2]

# Scale FFT amplitude to preserve Parseval's Theorem
X1 = X1 * np.sqrt(2 / nfft)
f1 = np.arange(0, nfft//2) * df

plt.figure(figsize=(8,4))
plt.title(f'Baseline FFT, Freq Resolution={df}Hz')
plt.plot(f1, np.abs(X1))
plt.grid()
plt.xlim([0, nfft/2])
plt.xlabel('Hz')
plt.show()

p1 = np.sum(x1[:nfft]**2)
p2 = np.sum(np.abs(X1)**2)
print('Parseval check:', np.allclose(p1, p2))

Here comes the tricky part. The signal is resampled in the following way, x2 = x1[::10] (downsampling by a factor of 10?). If that's really downsampling then new effective Fs should become 1024/10 = 102.4Hz, then the Nyquist frequency becomes $102.4/2 = 51.2Hz$, so the FFT output should only contain frequency components up to 51Hz. At least that's from the theoretical perspective. However the following output suggests otherwise. I can see the frequency resolution has improved significantly from 1.0Hz to 0.1Hz, but then the output only contain frequency components up to 12Hz instead of 51Hz. Why is that happen? What do I miss? Can anyone shed some light? Thanks

step = 10
x2 = x1[::step]

plt.figure(figsize=(8,4))
plt.title('Downsampled Signal')
plt.plot(t1[:100], x2[:100])
plt.xlabel('sec')
plt.grid()
plt.show()

nfft = 1024
fs2 = fs / step
df2 = fs2 / nfft
X2 = np.fft.fft(x2, nfft)[:nfft//2]

# To preserve Parseval's Theorem
X2 = X2 * np.sqrt(2/nfft)
f2 = np.arange(0, nfft//2) * df2

plt.figure(figsize=(8,4))
plt.title(f'Downsampled FFT, Freq Resolution={df2}Hz')
plt.plot(f2, np.abs(X2))
plt.grid()
plt.xlim([0, 50])
plt.show()

p1 = np.sum(x2[:nfft]**2)
p2 = np.sum(np.abs(X2)**2)
print('Parseval check:', np.allclose(p1, p2))


texts: No-lag causal bandpass filter kernel formula

texts:
Can anyone point me towards a formula for a good (for any values of "good") no-lag causal bandpass filter?

I am doing sound processing and I need to identify starting points of certain sounds (target precision <10 ms.) I've been using a Lanczos filter followed by convolution with a sine wave to do spectral analysis, but that gives me skewed results because amplitudes at any given time are being affected by stuff that's happening at later times (the filter is non-causal).

Simply cutting off the right half of the Lanczos filter resolves that particular problem, but the resulting spectrum is badly messed up (I get heavy sideband spillover). All other simple filters I tried are even worse.

Half an hour of googling turned up several mentions of Butterworth filters, lots of transfer functions (which I could not care less about), but no formulas of actual kernels.

EDIT: Let me try to be more specific. I start with a narrow-band filter with impulse response:

$$ h(t) = \begin{cases}
 C \, e^{i 2 \pi t/T} \operatorname{sinc}^2\left(\tfrac{t}{NT}\right) \qquad & |t| \le NT  \\
 0 &  |t| > NT \\
\end{cases} $$

where $\operatorname{sinc}(u) \triangleq \frac{\sin(\pi u)}{\pi u}.$

With $T=$2 ms and $N=$10, that gives me 0 dB at 500 Hz and <-40 dB below 440 or above 560 Hz. The support is $2NT$=40 ms. If I use it on an input that goes from zero to a pure 500 Hz sine wave at time $t_0$, response goes from exactly zero to just over zero at $t_0$, to half the maximum level at $t_0+NT$, and stabilizes when the entire area under the filter is a sine wave. The "effective" lag (time to 50% level) is 20 ms.

What I'd like to do is to bring this "effective" lag as close to zero as possible, while maintaining the overall discrimination ability (being able to distinguish between a 440 Hz sine wave and a 500 Hz sine wave.)


texts:
I'm reading about PSC (primary synchronization code) in UMTS. It is a so-called Golay sequence.
It is said it has a good aperiodic auto correlation properties.
What does it mean (in terms of aperiodic)?
What is the difference between good periodic/aperiodic auto correlation properties?
How does it differ from good cross-correlation properties?
The questions ares not about PSC in UMTS in particular. It's rather about terminology
Thanks in advance.


texts: How to find out the noise level of a signal

texts:
Is it possible to find the amount of noise or SNR of a signal from its time series output? For example; let x = randn(100,1)  and y = awgn(signal,-5,'measured'); be two signals. Here signal is an unknown source signal of length 10,000 data points.

In Matlab this is what I did: snr(x) = 0.5177 and  snr(y) = -34.3797 dB Shouldn't I get -5 dB as SNR?

I am unable to understand what these values convey as they are a ratio: How is 0.5177 coming? Since I am generating y with an SNR of -5 dB, so snr(y) should be -5 dB.  Please help in clarifying if these answers are right and how to make sense of it.


texts: What is a good book reference on the definition of "subharmonic frequency"?

texts: Is there a good book reference on the definition of "subharmonic frequency" as integer fractions of the highest available frequency?


texts: As a new to DSP concepts  want to understand the meaning of impulse response and why we do $h(n-k)$ in convolution and $h(n+k)$ in correlation.


texts: Programming the IDWT for image processing

texts:
I want to program the 2D inverse discrete wavelet transform (only 1 level) in the case of image processing. In the matlab website there's this diagram:

now, I want to program the IDWT with haar wavelets. What are the kernels of the low pass filter and the high pass filters?


texts: What exactly is "sparse representation"?

texts:
I saw a recommended topic for the final project in my university called "dsp and dip applications using sparse representation techniques (MATLAB, C, C++)".

I consider  taking this topic as my final project, but there are not so many things about it on the internet to understand what this is exactly.

So, I would like to hear some clear answers about it and some opinions about how productive would that be for me considering that I want to learn more and more about dsp and dip.


texts: Processing digital sound in Audiokit to sound like the 1960's

texts:
I am looking to imitate 70's analog equipments sound. I was wondering is there a way to create tape saturation / distortion inside Audiokit? I want to add that and randomized tape crackles and artifacts to the mix which I know is possible in AudioKit luckily.

From what I understand Tape Simulator plugins normally distort the lower frequencies more than than the high which might be possible in AudioKit if I create multiple audio streams (Lower end, higher end) using equalization and distort the lower end more. Is there a way to make the curve of the sound more linear?

How to model Tape Saturation (Audio DSP)?

Analog tape sounds further back in the virtual 3d space than digital. To achieve that I typically slightly cut off the higher end with a curve in parametric equalizers and add a convolution reverb to imitate natural space. Tape also has a much smoother smaller range for low frequencies which maybe I will cut off with a low pass filter.

I was wondering if anyone has anymore knowledge to achieve an effective tape simulation and to get a very retro sound in Audiokit? I am making a video camera app full of retro filters with grain and artifacts. It would be awesome if I could also have the sound feel retro as well.

This is my ideal clip for a frame of reference how audio sounded back then. https://www.youtube.com/watch?v=hEpSVCYok8s


texts: Difference between acoustic emissions and ultrasonic waves

texts: I need some clarification about the difference between acoustic emissions and ultrasonic waves. I have the idea that both are high-frequency signals and informative frequency bands more or less overlap.  Is there any difference between the waveform or in other characteristics that may the propagation speed through a medium. Moreover, will the velocity of both the signals through a medium (carbon steel) will be the same or will there be any difference.


texts:
suppose I have a noisy sine wave (as R code):

# Sampling frequency
fs <- 1000
# Sampling period
ts <- 1/fs
# Length of signal
l <- 1500
# Time vector
t <- seq(0, l - 1)*ts
# Signal
s <- 0.4*sin(2*pi*50*t) + sin(2*pi*120*t)
# Noise
x <- s + 2*rnorm(length(t))
# Plot signal
qplot(t, x) + geom_line()

and want to determine the amplitudes of the signal with:

# FFT
y <- fft(x)
# Two-sided spectrum
p2 <- abs(y/l)
# Convert to one-sided spectrum
p1 <- p2[1:(l/2 + 1)]
p1[2:(length(p1) - 1)] <- 2*p1[2:(length(p1) - 1)]
# Define frequency
f <- fs*seq(0, l/2)/l
# Plot spectrum
qplot(f, p1) + geom_line()

How can I decide whether a peak in the spectrum is indeed the amplitude of the signal and not noise? Is there a formula to define a threshold?

Thanks!


texts: Is the first derivative operation on a signal a causal system?

texts:
Please help me sort this issue out.

Consider a system whose output $y(t)$ is the first derivative of the input signal $x(t)$.

We can write the first derivative of an input signal as follows:

$$y(t)=\frac{dx(t)}{dt}=\lim_{\Delta t\to 0}\frac{x(t+\Delta t)-x(t)}{\Delta t}$$

Thus, we can see that the system's output must remember an infinitesimal change in future input in order to calculate the derivative, thus system possesses memory. Perfectly fine till here.

The problem comes in checking if it's causal. By the definition of the derivative above, the system has to subtract an infinitesimal future value with the current value, to calculate the derivative. Thus by definition, this has to be a non-causal system, right?

However, most of the textbooks and websites mention the first derivative as a causal system.

Please let me know if I am making any mistakes. Thank you.


texts:
I'm reading up on Autocorrelation, but I'm not sure I understand exactly how it works and what output I should expect. Am I right in thinking that I should input my signal to the AC function and have a sliding window input. Each window (of 1024 samples, for example) would output a coefficient between -1 and 1. The sign simply states if the line is upwards or downwards and the value states how strong the correlation is. For simplicity, lets say I don't have an overlap and just move the window 1024 samples each time. In a sample of 44100, would I get 43 coefficients and do I need to keep all of them?

Lets say I perform this for a 200 second signal, giving me 8600 coefficients. How would I use these coefficients to detect repetition and, in turn, tempo? Should I create some sort of neural network to group them, or is that overkill?

Thanks for any help.


texts: MMSE - How to minimize a complex error with respect to a set of real parameters

texts:
Suppose there's a complex signal $X(k)$ (where $k \in \{0, 1, 2,...,N - 1\}$) corrupted by additive complex noise. Its estimate $\hat{X}(k)$ is a linear combination of a set of real parameters $A_r$ ($r \in \{1, 2, 3,..., R$})
$$\hat{X}(k) = \sum_{r=1}^R A_rZ_r(k)$$
where $Z_r(k)$ is complex (and known).

I wish to obtain the real values of $A_r$ for which the error is minimized. If I simply differentiate the mean squared error (MSE) $\frac{1}{N}\sum_{k = 0}^{N-1}(X(k) - \hat{X}(k))^2$ with respect to each $A_r$ and set the resulting derivatives equal to zero, the values of $A_r$ I'll obtain will be complex, so that's not a solution.

My question is: how do I obtain the optimal values of $A_r$ such that the MSE is minimized, under the constraint that each $A_r$ should be real?


texts: How to calculate the similarity of two electric waves

texts:
I have two electric signals from different sources and I think that are similar.
I have the FFT of the signals and I have their spectrums of frequency  in separate arrays.
How should I proceed to find the percentual of similarity the two signals?


texts: Is there an asymmetric crypto analog to direct sequence spread spectrum?

texts:
In section 12.1.1.2 of Digital Communications: Fundamentals and Applications, Bernard Sklar makes this claim about spread spectrum signals:

  Therefore, not only can the spread-spectrum signal be made difficult
  to jam, but additionally, the signal's very existence may be rendered
  difficult to perceive. To anyone who does not possess a synchronized
  replica of the spreading signal, the spread-spectrum signal will seem
  "buried in the noise."

This seems analogous to symmetric cryptography; there is a symmetric "key" (the spreading signal) that allows one to both transmit and receive messages, and without the key one can neither transmit nor receive a message. Furthermore, without the key one cannot even detect the presence of a message (or, at the very least, it makes it more difficult to do so).

This makes me wonder if there's an asymmetric cryptography analog to direct sequence spread spectrum. I assume there would be a public and private "key" pair, and one could use the public key to transmit a message, but only the private key could be used to receive that message. And without the private key, it would be very difficult to even detect the presence of a coded message (i.e. only the private key could be used to "de-spread" the message).

Does such a technique exist?


texts: *Hello I used k mean method to detect edges and that’s my results as shown in the pic but I want to get more accurate results so I need to enhance my edge detection more what methods would you suggest such AS in points M,P AND N


texts: Can an Delta Dirac function (not response) ever be sampled?

texts: This is purely theoretical question but can an actual delta-dirac function ever be sampled? To make it clear again, I'm not talking about an impulse response from the impulse function. Lets say we have a truly instantaneous impulse as defined with delta-Dirac that the impulse with infinite(or finite) height and no width, an instantaneous pulse. This pulse is to be detected by sampling, but to find this instant signal, wouldn't this require an infinitely fast sampling or infinite sampling rate? Thinking this in terms of frequency domain, I believe delta function is 1 or has infinite spectrum. So would delta function ever be sampled unless we have infinite sampling rate?


texts: Trouble with inverse Z-transform and calculating of samples

texts:
I have a little problem. I have to solve this task but I can't.

  Z-transform of sequence $\{x(k)\}$ describe by the formula:

  $$X(z) = \frac{2.5 -3.15z^{-1} + 1.2 z^{-2}}{1-2.3z^{-1} + 1.2z^{-2}}$$

  Calculate of samples $x(-2), x(-1), x(0), x(1)$ and $x(2)$
  double-sides sequence, which correspond this Z-transform.

Thank you very much.


texts:
I have a series of data (single array). If I take this data and plot it, I can see that there are multiple peaks. However If I zoom in to a section of data, I see that there is substantial noise.

I like to detect how many peaks in the data using little CPU time and energy. I thought of detecting the slopes of these peaks (could be sharp or smooth rising or falling) and look at the number of slopes to determine the number of peaks. (2 slopes positive and negative for every peak)

Any pointers how I can calculate the slopes in the presence of noise?

The code will go into an embedded system, memory is limited so preferably I like to implement something that doesn't require any significant data copy.


texts:
I got the time domain signal from experimental work as per below. The blue is x-direction,green is z-direction and red is y-direction.

When i convert the signal into FFT using MATLAB, I got as per below.

.

The problem is, the peak should happen in z-direction(Ch3) but happen in x-direction(ch1). The coding as per below:

clc;
y1 = Ch1;
y2= Ch2;
y3 = Ch3;
N=6678;
Fs=512;

%Plot FFT Moving Vehicle signal Ch1
NFFT1 = 2^nextpow2(N);
Y1 = fft(y1,NFFT1)/N; % the division by N is to scale the amplitude
f_plot = Fs/2*linspace(0,1,NFFT1/2+1); % vector of frequencies for the plot
subplot(311), plot(f_plot,2*abs(Y1(1:NFFT1/2+1))) % multiplicated by 2 to recover the energy related
% to the negative frequencies since in this way only half spectrum is plotted
title('FFT Moving Vehicle Ch1');
xlabel('Frequency [Hertz]');

%Plot FFT Moving Vehicle signal Ch2
NFFT2 = 2^nextpow2(N);
Y2 = fft(y2,NFFT2)/N; % the division by N is to scale the amplitude
f_plot2 = Fs/2*linspace(0,1,NFFT1/2+1); % vector of frequencies for the plot
subplot(312), plot(f_plot2,2*abs(Y2(1:NFFT2/2+1))) % multiplicated by 2 to recover the energy related
% to the negative frequencies since in this way only half spectrum is plotted
title('FFT Moving Vehicle Ch2');
xlabel('Frequency [Hertz]');

%Plot FFT Moving Vehicle signal Ch3
NFFT3 = 2^nextpow2(N);
Y3 = fft(y2,NFFT3)/N; % the division by N is to scale the amplitude
f_plot3 = Fs/2*linspace(0,1,NFFT3/2+1); % vector of frequencies for the plot
subplot(313), plot(f_plot3,2*abs(Y3(1:NFFT3/2+1))) % multiplicated by 2 to recover the energy related
% to the negative frequencies since in this way only half spectrum is plotted
title('FFT Moving Vehicle Ch3');
xlabel('Frequency [Hertz]');

For information,when the sensor do the measurement,the reading for y-direction will fall at-1 to 1 G due to the gravity.

Is there anything I need to do before convert it into FFT? Shift maybe?

Thank you in advance.


texts: How to find ROC of system when input is two sided and output is one-sided

texts:
Given the $\mathcal Z$-transform of input $x[n]$ and output $y[n]$, how can I find the ROC of the system function $H(z) = Y(z)/X(z)$? I have
$$X(z) = \frac{2z\left(z-\frac{10}{3}\right)}{\left(z-\frac 13\right)(z-3)},\quad\text{ROC}: \frac 13 < |z| < 3$$
and
$$Y(z) = \frac{-1}{4}\frac{z}{\left(z-\frac 13\right)\left(z-\frac 14\right)},\quad\text{ROC}: |z| > \frac 13$$

If I substitute $\hat{X}(z) = 1/X(z)$, then $ H(z) = Y(z)\hat{X}(z)$. Then, I know that the ROC of $H(z)$ contains the intersection of the ROCs of $Y(z)$ and $\hat{X}(z)$. But how can I find the ROC of $\hat{X}(z)$? I don't know if the corresponding sequence $\hat{x}[n]$ is causal or not.


texts: Help obtaining the variance of a signal from a power spectrum?

texts: I have been reading this. The paper presents a power spectrum (attached). The power spectrum shows the power spectrum of velocity perturbations, $v$, in the solar atmosphere. If I want to obtain the average velocity, $\langle v\rangle^2$, do I integrate the power spectrum over frequency space? I tried this but the value I obtained seemed  far too small $\sim 10^{-5}\,\text{km}^2\,\text{s}^{-2}$. We can see that the integral is less than $0.1\times0.01=10^{-3}\,\text{km}^2\,\text{s}^{-2}$ by looking at the x and y range in the graph below. The paper suggests they observed an amplitude of the order $1\,\text{km}\,\text{s}^{-1}$. Which means integrating over frequency space gives a value which is far too small. Additionally such an integral will give dimensions of $\text{km}^2\,\text{s}^{-3}$ which is not the same dimensions as $\langle v\rangle^2$. Does anyone know what I should do to obtain $\langle v\rangle$ from the power spectrum below? I would like to get an understanding of what the y-axis means in this graph and hopefully figure out how to use this graph to obtain $\langle\dot{v}\rangle$.


texts:
The general structure of Sub-Band Coding (SBC) is shown below:

The diagram shows that after each of the analysis filters, decimation is performed by a factor M which corresponds to the number of channels. Would this not result in aliasing when this decimation by N is performed for the higher frequency sub-bands and why is DDC not required in each path to translate each sub-band to baseband before performing the decimation and thus avoiding aliasing?


texts: Expressions for steady-state LMS coefficient covariance

texts:
I've seen in the literature a few transient analyses for the LMS algorithm.  These mostly focus on the mean estimation error for the LMS tap weights as well as the limiting signal estimation error covariance compared to the Weiner filter.
I am interested in simple expressions for the limiting coefficient error covariance.
In other words, if the LMS coefficient update is given by:
$$
 \hat{\bf w}_{k+1}=\hat{\bf w}_{k} + \mu {\bf x}_k e_k^*
$$
and the true coefficient vector is $\bar{\bf w}$, define the difference as:
$$
{\bf e}_k={\bf \hat{w}}_k-\bar{\bf w}.
$$.
What can be said about $\lim_{k\rightarrow \infty} E[{\bf e}_k{\bf e}_k^H]$?


texts: The frequency response function (FRF) fails to detect the antiresonance of a system

texts:
I am trying to identify a vibrational systems by computing the frequency response function (FRF) of the system when a chirp signal is applied to its input. After comparing the FRF computed and the bode plot of the transfer function I found that the resonance has been well detected but the anti-resonance has not.

The system has the following form:

$$\textrm{sys} =  \frac{s^2+\omega_{n_{z_1}}^2}{s^2+\omega_{n_{p_1}}^2}, \quad \textrm{where } \omega_{n_z} = 10\textrm{ Hz} \textrm{ and } \omega_{n_p} = 100\textrm{ Hz}.\\$$

I have a hypothesis:

The steady state response of a linear system to a sinusoidal signal is another sinusoidal signal with the same frequency with an amplitude and a phase that depend on the system. When I applied a chirp signal, the system never reached the steady state. So while the output of the system to $\sin(\omega_{n_{z_1}}\cdot 2\pi)$ should be zero, I am still obtaining some sinusoidal signal of frequency $100\textrm{ Hz}$ as a consequence of stimulating (weakly) that part of the spectrum with the transient response.

In my opinion the correct computation should be measuring only the components of the outputs that correspond to the frequency of the input at certain particular moment.

Do you know how to correct this problem?
Is it possible to apply a moving filter or window to correct this problem?

Here is my code:

% Identification of the Frequency Response of a Transfer Function with
% Resonance and Antiresonance peaks
%% Cleaning the house
clc;
clear;
close all;
%% Definition of the System
wn_z1 = 10*2*pi;           % Anti-Resonance Natural frequency 10 Hz
wn_p1 = 80*2*pi;           % Resonance Natural frequency 80 Hz
chi_z = 0;                 % Damping factor of the zeros
chi_p = 0;                % Damping factor of the poles
s     = tf('s');
sys_1 = (s^2+2*chi_z*wn_z1*s+ wn_z1^2)/(s^2+2*chi_p*wn_p1*s+wn_p1^2);
zpk(sys_1)
%% Definition of the Chirp Input
h                    = 0.001;
time                 = 0:0.001:60;
freq_ini             = 0.1*2*pi; % Initial frequency considered (0.1 Hz)
freq_final           = 240*2*pi; % Final frequency considered (240 Hz)
chirp_input          = 1.5*chirp(time,freq_ini,time(end),freq_final,'logarithmic');

frequency_evolution  = logspace(log10(freq_ini),log10(freq_final),length(time));
figure(1);
subplot(2,1,1);plot(time,chirp_input); xlabel('time (s)'); ylabel('Chirp    amplitude')
subplot(2,1,2);plot(time,frequency_evolution/(2*pi)); xlabel('time (s)'); ylabel('Frequency (Hz)');
%% Execution of the simulation
[sys_output,~]       = lsim(sys_1,chirp_input,time);
%% Computation of the FRF
fft_input   = fft(chirp_input,2^12);                                           %% Considering zero padding
fft_output  = fft(sys_output ,2^12);                                           %% Considering zero padding
fft_input   = fft_input(1:length(fft_input)/2);
fft_output  = fft_output(1:length(fft_output)/2);
fft_output  = fft_output(:);
fft_input   = fft_input(:);
FRF              = abs(fft_output)./abs(fft_input);
frequency_vector = linspace(0,1/(2*h),length(FRF));

%% Comparing the Bode Plot (Matlab) with the FRF (own implementation)
opts = bodeoptions('cstprefs');
opts.FreqUnits = 'Hz';
close all;
bodemag(sys_1,opts,'r'); hold on;
semilogx(frequency_vector, 20*log10(FRF ),'b')
xlim([0.1,frequency_vector(end)])
legend('Matlab Bode','FRF Identification')


texts: bandpass implementing with DWT in matlab

texts:
In matlab help in find the implementation of bandpass as following

  fs=8800;
  n = 0:159;
  x =sin(n)+sin(n)*pi; % whatever here say it music
  d = fdesign.bandpass('Fst1,Fp1,Fp2,Fst2,Ast1,Ap,Ast2',1/4,3/8,5/8,6/8,2800 /fs,200,3000/fs);
 Hd = design(d,'equiripple');
 %Apply the filter to the discrete-time signal.
 y = filter(Hd,x);
 freq = 0:(2*pi)/length(x):pi;

-> xdft = fft(x);
->ydft = fft(y);

 plot(freq,abs(xdft(1:length(x)/2+1)));
  hold on;
 plot(freq,abs(ydft(1:length(x)/2+1)),'r','linewidth',2);
 legend('Original Signal','Bandpass Signal');
 figure(2)
 plot(xdft(1:length(x)/2+1),'r');
 hold on
 plot(ydft(1:length(y)/2+1));
 legend('Original Signal','Bandpass Signal');

as you see the arrows refer to the dft function to make transform I want to make use the dwt function to transform to freqency domain and make it seen to me so I can easily use bandpass filter or if condition

thank you for helping


texts:
In the paper about speaker verification I can read that

Two quantities have
to be set: the length of the window and the shift between two
consecutive windows. For the length of the window, two values are most often used: 20 milliseconds and 30 milliseconds. [...] For the delay, the
value is chosen in order to have an overlap between two consecutive windows; 10 milliseconds is very often used.

and later

Once an FFT algorithm has been chosen, the only parameter to fix for the FFT calculation is the number of points for
the calculation itself. This number N is usually a power of 2
which is greater than the number of points in the window,
classically 512.

However, in the matplotlib's specgram I can only see NFFT and noverlap parameters which, according to the documentation, correspond to the second quotation.

How can I then set the window's length or what is the default value of it?


texts: All pass filter design for group delay compensation

texts:
I've designed an IIR bandpass elliptic filter in Scilab with the following parameters.

Sampling Frequency = 25000 Hz
Lower cutoff frequency = 100 Hz
Upper cutoff frequency = 150 Hz
Pass band ripple = stop band ripple = 0.001 ~ 0.005 dB

The transfer function comes out to be

(-0.0006050 + 0.0024134z - 0.0030118z^2 + 1.075D-18z^3 + 0.0030118z^4 - 0.0024134z^5 + 0.0006050z^6) / (0.9600937 - 5.7967985z + 14.586174z^2 - 19.578664z^3 + 14.785532z^4 - 5.9563371z^5 + z^6)

I need to make the group delay in the pass-band, i.e. 100 Hz - 150 Hz for which I'm looking at using an IIR all pass filter. There a function in Matlab called iirgrpdelay which accepts a vector of group delays corresponding to a range of frequencies between 2 frequency values (the pass band range) & generates an all-pass filter. The vector contains the delays needed for the frequency range in order to make it linear. An example for this can be found here.

This function isn't available in Scilab or Octave & I'm not sure how to design such an all pass filter which accepts a vector of the required group delays as input. What would be the technique used for this?


texts: I have a receiver which is left to open for sensing. I'm using the periodogram method for estimating the PSD. I need to calculate the SNR of the incoming for which I need to calculate the variance. Can anyone suggest me ways to calculate the variance from the PSD?


texts: downsample signal to a specific number of samples in Python

texts:
The Python Scipy library provides several functions to downsample signals, but they all have limitations:

The resample function is based on Fourier method, which means it assumes periodic signals.
The resample_poly function assumes "values beyond the boundary of the signal to be zero."
The decimate function can only downsample signals by some integer factors, not to  specific number of samples.
The upfirdn function requires FIR filter coefficients as inputs, which I am not sure how to get.

I have a signal that's NOT periodic and the values beyound the boundary of the signal are NOT zeros, and I want to downsample the singal from 611 samples to 100 samples.

Is there a simple way to do that in Python?

Thanks.

EDIT: The reason why I need to downsample signals to a fixed number of samples is that I have signals with different length, and I want to use those signals to train a 1D CNN for classification/regression tasks, which requires inputs of fixed size; I prefer not to use RNN/LSTM, although those models can accept inputs with various sizes.


texts: How to fix frequency drift when producing an arbitrary length signal by repeatedly looping over a single-period waveform?

texts:
I have an audio application that produces an arbitrary length signal by repeatedly looping over a single-period waveform. In pseudo-code:

// Build a single period of a waveform at the given frequency
samples_per_waveform = int(sampling_rate / frequency)
waveform = make_sin_waveform(samples_per_waveform)

...

// Repeat the single period waveform to produce a continuous signal
count = int(duration * sampling_frequency / samples_per_waveform)
repeat `count` times:
     write(waveform)

Since the number of samples per waveform is obviously an integer, I think this is subject to rounding error:

As an example, when rounding to the nearest integer, at a 5kHz sampling rate, a single-period waveform for a 300Hz signal would require 17 samples--and thus leads to an actual frequency of 294Hz. I understand the higher will be the sampling rate, the lower will be the rounding error. However, it will still exist.

Is there a way to fix that, or is this approach doomed from the start?

EDIT: As a possible workaround, I envision to keep track of the rounding error in the output loop, and duplicating one sample when the error exceeds a sample duration. I don't know if this is clear enough, but it is somewhat similar in my mind to what the Bresenham line drawing algorithm is doing.


texts:
I am very new in signal processing and using digital filters. I have to use a low-pass filter to analyze my data in LabVIEW and have a question about it. Any help and advice is appreciated.

I am trying to simplify my problem here:

Let’s say there is a digital sine wave (made by LabVIEW) with $V_{offset}=1 \ \mathrm{V}$, $V_{peak}=0.1 \ \mathrm{V}$, $f=10 \ \mathrm{kHz}$, $N=2000$ (number of samples), and sampling rate $f_s=200 \ \mathrm{kHz}$. Now, if I pass this signal through a low-pass filter with cutoff frequency $f_c=1 \ \mathrm{kHz}$, then the output should be a constant number equals the DC offset (here $1 \ \mathrm{V}$), is it true?

Another question is the concept of “cutoff freq” and “sampling freq” as the inputs of the filters in LabVIEW. Cutoff frequency as an input of a filter makes sense to me but what is that “sampling freq” ? Can anyone explain to me please? I am very confused. Is it the same rate at which the sine wave is created? For example in the attached code, what is the real cutoff frequency (with $f_l=200000$ and $f_l=1000$)?

I have attached the screenshots of the Front panel and Block diagram of my simple vi.


texts: Is sparse signal in time domain is sparse in Frequency domain

texts:
In OFDM system, if we have a multi-path environment, we supposed to have a sparse received signal in time domain, then  $FFT$ is performed to convert it into frequency domain.

My question, does the results of FFT in frequency domain must be sparse since its equivalent time signal is sparse too ?

thank you.


texts:
I try to develop a guitar game in android platform.

And I need to do the real-time pitch detection of guitar chord.

I will get the input from the microphone, and then analyze the input (the input playing which kind of guitar string/chord).

I find two kinds of method that I can use, one is YIN, another one is FFT.

Which method can get better performance and exact result?


texts: Kalman Filter for estimating position with nonconstant velocity & acceleration

texts:
I am trying to estimate the position & head direction of a rodent going through a 2D environment (a circular surface of 1m radius). Above his head is an overhead camera which records 4 LEDs mounted on his head. The position is calculated as the centroid of the LEDs positions and the head direction as the angle of the LEDs (see images).

The issue is that wires attached to his head occasionally occlude 1-4 of the LEDs and thus prevent/complicate measurement. Furthermore, the head direction estimation is particularly subject to lots of noise due to nonhorizontal head movements which add a 3rd dimension to the LED movements which are assumed to be parrelel to the surface (thus introducing a sort of parallax).

I thought a kalman filter would be a great way too smooth this data. However, in Matlab it seems that to implement this I would need to assume either constant acceleration or velocity which is not the case since the rodent is freely moving.

Is a standard Kalman filter still possible? Maybe a nonlinear (extended) KF? Are there any other better solutions?

Thanks
Roman

Rodent in a square environment. I included this image since here the overhead wires and LEDs are more easily visible

Rodent in the circular environment. Notice that 1 of the LEDs (green) is occluded by the wire:


texts:
A basic theorem in communications is the matched filter maximizes the SNR at sampling. I'm a little confused on how this relates to discrete time systems and sampling rate.

Normally if you sample at Nyquist rate it's "good enough" but how does this interact with the discrete time correlator? Does sampling at a higher rate benefit the system? Does up-sampling the signal and passing it through the correlator yield the same benefit?

I understand the idea of Noise being White and the output of white noised passed through the correlator limits its PSD but lets say your analog signal is pre-filtered to remove bands outside of some region and your sample rate meets the Nyquist rate for this bandwidth. Is there any benefit to over-sampling and then discrete time matched filtering?


texts: Filter for Reducing Background Noise in Image

texts:
I have images where there is a lot of "black" background (few shades of black). (In many images at least half the pixels are background). I need to get interest points from the image, but because the background is not smooth (although black) , I get a lot of "false" key points in the background.
What filter/s will fit best to keep the details in the image but will smooth the background as much as possible ?
I tried median filter, but the size needed to reduce background noise ([5 5]) was to strong for the small details in the image I wanted to keep sharp, same for Gaussian with very small sigma.

Thanks


texts:
I wrote an algorithm to do bicubic interpolation of an image.
I used the method desribed in the wikipedia page.
On simple images, the result looks good, but on more complex ones, I got strange artifacts in non-smooth zones.

What could be a probable source of those artifacts, and what can I try to get rid of them ?

My code is on github, but I checked it several times, so I don't think the problem comes from it but rather from the simple implementation suggested in wikipedia

I'm using this particular equation (3rd-order Hermite polynomial interpolation) to do the interpolation (in one dimension):

$$x(n+t)=\frac{1}{2}
\begin{bmatrix}
1 & t & t^2 & t^3 \\
\end{bmatrix}
\begin{bmatrix}
0 & 2 & 0 & 0 \\
-1 & 0 & 1 & 0 \\
2 & -5 & 4 & -1 \\
-1 & 3 & -3 & 1 \\
\end{bmatrix}
\begin{bmatrix}
x[n-1] \\
x[n] \\
x[n+1] \\
x[n+2] \\
\end{bmatrix}$$

where $x[n-1]$, ... $x[n+2]$, .. are the values of the surrounding pixels (in one dimension), and $t$ is the fractional coordinate (where $0 \le t < 1$) of the interpolated pixel between the adjacent pixels $x[n]$, $x[n+1]$. $n$ is an integer and $x(n)=x[n]$ for any integer $n$.

The expected output:

My result, look closely at the bird head,a lot of white pixels have appeared.

---EDIT---
To show that my code does produce coherent results, here is the example of my scaling up on simple image with 16 red pixels.

And here is the ouput of the library (which also use 3rd-order Hermite polynomial interpolation):

And here is my results:

In this case the difference are IMO totally acceptable, but it shows that my interpolation does something legit.


texts: Has this Principal Component Analysis (PCA) been done correctly?

texts:
I have a set of 3D data points, indicated by the blue color in the picture below.

I then project them onto the x-y plane, i.e. setting z values of all the points to 0, shown by the yellow color below.

Then, I perform PCA on the 2D data points (the yellow ones).

I end up with the 1st Principal Component being the green line, and the 2nd being the yellow line.

As I am a newbie in PCA, I am not sure whether I have done something wrong or not.

Because in my intuitive, the 1st Principal Component should be like in the picture below:

I mean unlike the above picture given in wiki, my 1st component is NOT obvious and seems to make no sense. Did I do something wrong?

Edit:

To avoid the confusion, the 2D view from the above is as follows:


texts: How are Object Recognition and Content-based Image Retrieval interlinked?

texts:
I need help on understanding the linkage between Object Recognition and Content-based Image Retrieval. Does OR use CBIR to find similar objects with the help of algorithms like SIFT etc or does CBIR use OR to improve its results on finding images showing the same object?

CBIR could use color histograms to find similar images but that has not too much in common with or so CBIR can stand alone i guess. OR needs to retrieve information which could come from learned images thus I suppose OR needs CBIR in some way. That would lead to the conclusion that OR is dependent and builds on CBIR, is that correct?


texts: Why/How a Gaussian Channel is both band-limited and Power limited?

texts:
All physical channels are band-limited, such that,

  X(f) = 0 for |f| > W

A Gaussian Channel is limited in both power and bandwidth. For this reason we encounter phrases, such as, "A band-limited Gaussian Channel" (for e.g., in the opening statement of Channel Capacity theorem).

I have two questions:

Why/How a Gaussian Channel is band-limited and Power limited?
What are physical implications of using such a channel?

Please help/guide me to understand this.


texts:
I'm trying to quantize an image 8 bits to 4 or 2-bits uniformly. I searched internet, interestingly I could not find what I want exactly. Then I wrote an simple code for it myself. I'm curious about whether there is a build-in function in MATLAB which convert 8-bit image to 4-bits uniformly . My results using the methods from internet is not good. Am I doing something wrong when I use these methods? Thank you in advance.

Method from internet giving strange results

reducedImage = uint8((single(monalisa)/256)*2^4);

Method writen by me. I'm putting code to clarify what I want to do.

monalisa_2= monalisa;
figure
monalisa_2(monalisa<=63)               = 32  ; %% 00
monalisa_2(monalisa<=127&monalisa>63)  = 85  ; %% 01
monalisa_2(monalisa<=191&monalisa>127) = 159 ; %% 10
monalisa_2(monalisa<=255&monalisa>191) = 223 ; %% 11
imshow(monalisa_2)


texts: Measure of intermittency/continuousness of a signal

texts:
I have three signals (below) each having the same standard deviation, however, are clearly very different temporally.

Is there some such metric that could be calculated for each of these signals to give an indication of how they vary over time, such as a measure of continuousness or intermittency if such a metric exists.

If possible, I would like to do this by calculating a single number for each of these signals, rather than windowing into smaller chunks.

Thanks


texts: How to initialize shift registers in polyphase filter based resamplers to avoid or compensate for group delay?

texts:
I have a signal a = [a(t0),a(t1),a(t2),..,a(tN-2),a(tN-1)]

I'm using a polyphase filter based resampler with filter length $L$ of $8$, fixed interpolation rate $I$ of $10000$ and decimation rate of $D$ to resample the signal $a$. After resampling, I'm seeing a group delay of $N = \frac{L \times I}{2 D}$ samples.

I'm able to partially compensate for the group delay by initializing shift register S with the samples of the signal a in the following configuration,

S = [a(t0),a(t1),a(t2),a(t3),a(t0),a(t1),a(t2),a(t3)]

The above configuration compensates for the group delay, but does not reconstruct the signal from t0 to t4 accurately. What would be the ideal way to initialize the shift register to reconstruct the signal with better accuracy?


texts: Using MIMO with DS-CDMA and spread spectrum

texts:
I have a question about using MIMO with DS-CDMA. suppose we are using 4 TX antennas, and 4 RX antennas, then the Reyleigh channel generated randomly  of MIMO system will be of dimension  H = [4x4],  assuming also the transmitted signal using QAM modulation is x done as below:

x =
[0.7 + 0.7i;
0.7 - 0.7i;
-0.7 + 0.7i;
-0.7 -0.7i];

So the received signal r is supposed to be r = H*x , which is of dimension of [4x1].

Now, suppose we are spreading that signal with PN/walsh code of dimension [4x1], by using the function of "kron" in MATLAB, so the dimension of the new transmitted signal will be [16x1], equivalent that x_1 = kron(x,c); where c is the code used to spread the signal before transmission.

My question is, how can we transmit that signal after spreading it x_1 over the noted channel H of dimension H with dimension [4x4]?

Thank you


texts: How to verify the obtained SNR of the noisy signal after adding noise?

texts:
I defined the SNR as follows:
$$\rm SNR = \frac{\text{Signal power density}}{\text{Noise power density}}$$
Given noise-free signal, I added the noise to this signal by using the formula above. Now, I am looking for a way to verify if the obtained SNR of the noisy signal is equal to the SNR value that I set when I generated the noise.
Thanks in advance!


texts: How to apply Hamming Window and compute FFT?

texts:
I'm borrowing an EEG device whose accompanying software also computes the FFT of the data it collects. The software provides the user with a file of EEG values and a file of FFT values computed from the EEG file, and I wanted to learn how to compute the FFT values on my own.

Here's the manufacturer's description of how they do it: "Each path contains 129 decimal values with a range of roughly -40.0 to 20.0. Each array represents FFT coefficients (expressed as Power Spectral Density) for each channel, for a frequency range from 0hz-110Hz divided into 129 bins. We use a Hamming window of 256 samples(at 220Hz), then for the next FFT we slide the window 22 samples over(1/10th of a second). This gives a 90% overlap from one window to the next. These values are emitted at 10Hz."

Going through MATLAB's documentation I ended up with this code:

X = data(:, 1);
h = hamming(256);
buf = buffer(X, 256, 234);
hWin = buf.*repmat(h(:), [1, size(buf, 2)]);
Y = fft(hWin);
P2 = (abs(Y))';
P1 = P2(:, 1:256/2 + 1);
P1(:, 2:end-1) = 2*P1(:, 2:end-1);
P1 = 10 * log10(P1);

But it doesn't produce the sames values that I see in the four FFT files given to me by the manufacturer. First, let me start by saying that their EEG file contains a 193011 X 4 matrix of data, and each column corresponds to the data collected by one of the four channels/sensors on the EEG device. And with each column of data, they were able to compute four separate files, which they named raw_fft0, raw_fft1, raw_fft2, and raw_fft3. Each raw_fft file contains a 8607 X 129 matrix whose values should be somewhere between -40 and -20. The actual observed values for their FFT files are somewhere between -30 to -20 to the low-to-mid 40s.

On the other hand, my FFT files are of the dimension 8774 X 129, and the range of values in my files are bigger than what the manufacturer says they should be. Would anyone mind pointing out what I've done wrong or pointing me to information that could clarify what I've done wrong?


texts: ADC sampling rate in OFDM baseband processor

texts:
I have optical OFDM baseband TX and RX implemented in FPGA. The design is tested using one FPGA board in loopback connection from DAC to ADC via coaxial cable. It works fine, the BER is about 1.10^-8 for QAM-16 modulation. The ADC output is captured, and looks like this:

When I implement the TX and RX in two board (one board for TX and another for RX), there is problem. The ADC output (captured from FPGA) looks like this:

The amplitude of the OFDM signal is reduced, so the synchronizer (using training symbol) fails to detect the start of symbol.

System specification:

DAC sampling rate: 50 MHz, resolution 14-bit
Analog LPF (after DAC): 50 MHz
Analog LPF (before ADC): 50 MHz
ADC sampling rate: 50 MHz, resolution 14-bit

Question:

What is the mistake in design, is it because of the ADC sampling rate, which is equal to the OFDM bandwidth?
When sampling OFDM signal, should I use ADC sampling rate larger than 2x bandwidth (nyquist)?


texts: How to recover data/information from real complex baseband signal?

texts:
I have recorded real complex baseband signal. I want to recover data from the same. The complex baseband signal is basically a composite CDMA signal. How can i recover that information ?
Original Tx Composite signal

Real Rx Baseband signal

I know fundamentally I just have to multiply the composite CDMA signal with the relevant PN code in-order to get the original signal back. code length is 1024.
p_recovered = np.array([], dtype = float)
for i in range(p_len):
    p_recovered = np.append(p_recovered, 1.0*sum(composite[i*codelength:i*codelength+codelength])/codelength)
p_recovered = np.repeat(p_recovered, codelength)

r_recovered = despread(composite, r_code, codelength)
q_recovered = despread(composite, q_code, codelength)

I did that in simulation(below), and recover the original messages back.

My question is how can i recover information from actual modem received signal ?


texts: How do I plot the square of the amplitude response?

texts:
I have calculated the transfer function of an FIR filter

$$
y[n] = x[n] + α · x[n − R]
$$
This is what I have
$$
H(z) = 1 + αz^{-R}
$$

Now I should plot the square of the amplitude response. So I calculated the frequency response by inserting e^jωTs for z. I should choose α = 1/2 and R = 3. So I get

$$
H(e^{jωTs}) = 1 + αe^{-RjωTs} = 1 +  \frac{1}{2} e^{-3jωTs}
$$

$$
|H(e^{jωTs})|^{2} = (1 + αe^{-RjωTs})^{²} = 1 +  \frac{1}{2} e^{-6jωTs}
$$

I am also not sure here about the algebra of the square of the magnitude.

Now I should plot from 0 to ωs and I think I am missing Ts. I should probably just choose a value, right? How is the usual approach?


texts: How to calculate the delay generated by a 2nd order all-pass filter

texts:
I'm trying to simulate the dispersion effect of a stiff string in a digital waveguide system, by using second-order all-pass filters.

I use a second order all-pass filter with this transfer function:

H(z) =    a2 + a1 * z^-1 + z^-2
        --------------------------
         1 + a1 * z^-1 + a2 * z^-2

Let's suppose the filter coefficients are:

a1 = -1.86377279
a2 =  0.86982215

With the above coefficients, if I insert the 2nd order All-Pass filter inside a digital waveguide loop with a delay line of 168 samples (=> 262.5 Hz if the playback rate is 44100), I get a nice dispersion effect (a little exaggerated to study the phenomenon) but I get a big delay too in the output signal, so now my note is detuned.

A little simplified schema:

impulse --> (+) -------------------------------+-----> output
             ^                                 |
             |                                 |
         very much                             v
          detuned!                       |DELAY LINE|
             |                           |168 smpls |
             |                                 |
             +---------|2nd ord.|<-------------+
                       |ALL-PASS|

I noticed that if I manually modify the delay line length (by subtracting 41 samples from it) to compensate for the delay induced by the All-Pass filter, I could get the right tuning of the original note (about 262.5 Hz).

Now my question is:

is there a formula to know exactly how much I should modify the delay line length to compensate for the delay generated from the 2nd order All-Pass filter with arbitrary a1, a2 coefficients?


texts:
In a quantization process what happens when the signal range is more than the quantizer range?
Suppose i have a sinusoidal signal which ranges from -5 to +5 but the quantizer which has four levels ranges from -3 to 3.now for the samples which are greater than 3 and less than -3 what will happen?


texts: Is ideal LPF a LTI system? And is it time invariant?

texts: Ideal filter time-domain impulse response is $$h(t) = 2B \operatorname{sinc}(2Bt)$$ It is not causal nor BIBO stable. Is it time invariant? And how can i prove it?


texts: Fractional/Interpolating Delay Line still sounding glitchy

texts:
I'm trying to implement a simple digital delay line. I want it so that when the user changes the delay amount, it sounds smooth and not glitchy. Currently I'm implementing a fractional delay line but I'm still getting a very glitchy sound when I change the delay amount.

What else do I need to implement besides a fractional implementation?
I can post code if need be.

EDIT:
Here's some code-

- (void)setDelay:(Float32)dly {
    // get delay amount
    delay = dly*srate;
    if (delay != 0) delay -= 1;

    // Get Delay difference
    Float32 outPointer = w - delay;
    // wrap write
    while (outPointer < 0) outPointer += DELAY_MAX;

    // Integer conversion
    r = (int)outPointer;
    // Wrap read
    if (r == DELAY_MAX) r = 0;
    // Fractional Calculation
    alpha = outPointer - r;
    // Get interpolation values
    omAlpha = 1.f-alpha;
}

// interpolation method
- (Float32)nextOut {
    // If we can get the next output
    if (doNextOut) {
        // First 1/2 of interpolation
        nextOutput = cBuf[r]*omAlpha;
        // Second 1/2 of interpolation
        if (r+1 < DELAY_MAX) {
            nextOutput += (cBuf[r+1]*alpha);
        }
        else {
            nextOutput += (cBuf[0]*alpha);
        }
        // set false
        doNextOut = false;
     }
     // return interpolation
     return nextOutput;
}


texts: Cyclic convolution between two matrices in matlab

texts: How can I perform cyclic (periodic) convolution between two matrices in Matlab? And also, what would be the normalization term?


texts:
Precisely how does the performance of an FIR filter depend on its design criteria?
To illustrate, here is a specific question that probably does not have a convenient answer:

The Parks-McClellan algorithm converges to a FIR filter optimally aligned to design parameters. With it one can design a low-pass filter with minimum pass- and stop-band ripple:

Take $\mathcal{S}_n$ to be the collection of symmetric FIR filters of order $n$:
\begin{align}
\mathcal{S}_n := \left\{(h_1,\dots,h_n)\in \mathbb{R}^n: h_k=h_{n-k+1} \text{ for each }k\right\}.
\end{align}
Write the frequency response of $h\in \mathcal{S}_n$ as $\hat{h}(\omega)=\sum_k h_ke^{-j2\pi \omega k}.$
Specify a "stopband weight" $\alpha > 0$, cutoff $\omega>0$ and transition width $\varepsilon>0$. Define filter mismatch as follows:
\begin{align}
    \gamma_n(h; \alpha,\omega,\epsilon)\!
:= \min_{f\in [0,\omega]} \left|1-|\hat{h}(f)|\right| + \alpha\cdot \max_{f\in[\omega+\epsilon, 2\pi]} |\hat{h}(f)|.
\end{align}
The Parks-McClellan algorithm with parameters $(\alpha, \omega,\varepsilon)$ converges to an $h_{best}$ which minimizes $\gamma_n(h; \alpha, \omega, \varepsilon)$ over all $\mathcal{S}^n$.

What estimating expressions depending on $(n, \alpha,\omega, \varepsilon)$ provide tight upper and lower bounds for $\min_h \gamma_n(h;\alpha, \omega,\varepsilon)$ ?
Alternatively, fixing $n,\ (\alpha,\omega, \varepsilon)$  and $\Delta > 0$, are there bounds for $\min_h \max_{f\in[\omega+\epsilon, 2\pi]} |\hat{h}(f)|$ among those $h\in \mathcal{S}_n$ that have $\min_{f\in [0,\omega]} \left|1-|\hat{h}(f)|\right| \leq \delta$


texts: Ideal low-pass filter in the Frequency Domain

texts:
I try to create a low-pass filter to apply it for an image.

$$H(U, V)=1 \text{, if }D(U,V) ≤ D_0 \text{, and } 0 \text{ otherwise}$$

and $D(U, V) = \sqrt{U^2 + V^2}$

U,V = np.shape(padded)
H = np.zeros((U,V))
D = np.zeros((U,V))
#cut off point
D0 = 5
for u in range(U):
    for v in range(V):
        D[u,v] = np.sqrt(u**2+v**2)
        if D[u,v] <= D0:
            H[u,v] = 1
        else:
            H[u,v] = 0
plt.imshow(H,cmap='gray')
plt.show()

After I Run it I expect this image

But I get this circle at the left top of the image as bellow :

so I think I should get the center of $\frac{U}{2}$ and $\frac{V}{2}$ but I do not know how.

How can I get the expected filter image?!


texts:
I'm working on a 10-second sound, sampled at 44.1 khz.

I want to do filtering, and have a desired EQ (equalization) curve that varies over time, as suggested here (here $f0=250\ Hz$)

How to implement a time-varying filter?

How to go from this bad idea of zeroing FFT bins

[0, 0, ..., 0, 1, 0, ..., 0, 1, 0, ..., 0, 1, 0, ... 0]
   |           |             |             |
   |          freq bins for f0, 2*f0, 3*f0  are kept
   |
 most freq bins are zeroed

to a good filtering?
Then, how to make this vary over time? Should I use overlap-add, etc.? Then if so, can't I do it inside a STFT framework?
Would you have a pseudo code or Python demo?

PS: I've read http://recherche.ircam.fr/anasyn/roebel/amt_audiosignale/VL3.pdf, but I can't go from here to a working algorithm.


texts: How to analyze a pattern of discrete events belonging to a 24 hour and ~90 minute cycle?

texts:
I'm a sleep researcher and since last december I've been collecting data on mydreaming and awakening in the middle of the night.

Dreaming occurs within a well-known 80-120 minute "cycles", with about 5-6 cycles within a typical night. Each cycle is defined by 4 stages, from light sleep to deep sleep to dreaming. I'm lucky enough that most of the time I wake up after or close to the dreaming sleep stage.

The start of the sleep cycles depends on the choice of bedtime, which in turn is governed by a larger, ~24.5 hour circadian rhythm, which is entrained by the light/dark cycle.
The circadian rhythm is exhibited by all cells within a body as well as major organ systems.

What I'm trying to do is better understand the circadian rhythm by analyzing discrete awakening events that I experience over night. On the attached screenshots;
Each row represents one day of data. The events are plotted in real time.

Each discrete awakening from a  dream is marked in green

Each non-dream awakening is marked in red

Each lucid - conscious dream is marked in cyan.

Each morning "get out of bed" time is marked in Orange

The onset of the first sleep cycle is marked in black.

Each row has 2 rulers:

The purple ruler's ticks are 90 minutes apart - the approximate sleep cycle length

The cyan ruler represents real time starting with 0 midnight, going to noon.

Together, the markers form vertical patterns, with dream awakenings happening at approximately the same times on subsequent days, with the pattern sloping one way or the other.

How can I analyze hundreds of discrete data points like that to see if they do form some kind of a waveform that persists from day to day?

I would be grateful for the keywords related to the transforms or algorithms that can help make more sense of such data, as well as any engineering articles dealing with a similar problem, or any other input.

Thank you!

Full 90 day sleep history


texts: How exactly does 8PSK work on each I and Q leg?

texts:
I am working on making a modulator in VHDL and am implementing 8PSK. Currently, I have BPSK and QPSK looking good, and I have verified the spectrums. However, now that I am moving to 8PSK, I am stuck on how the data mux to the IQ legs should work?

Currently, my implementation is two-fold for QPSK. My first method splits the bits into the I and Q legs for even/odd bits (I am not pulse shaping, purely PCM), then I basically take this I or Q bit and xor it with the MSB of the signed DCO value, in-phase for even bits and quadrature for odd bits (I am using a digitally controlled oscillator as my carrier). The second method I have is I just generate the symbols for QPSK and feed my DCO with a phase shift at 45, 135, 225, and 315 (this saves logic as I only have to use 1 DCO). I know I can use the phase request for 8PSK as well, however, this method produces a lot of spurs.

If I want to perform the following block diagram for 8PSK, how does this work? Am I accumulating 3 bits, or one symbol, for each leg, and performing a phase shift? I'm really struggling to figure out how to do this. I know I can just accumulate one symbol and send a phase shift to the DCO, but this doesn't use I and Q and feels like a crappy way of doing this, especially considering it doesn't work well when I tried it for QPSK. Here is the fundamental question, when I go from serial to parallel through the data MUX, what exactly is going onto the I and the Q leg? If it is still 1 bit, how is this any different than QPSK?


texts: Calculating covariance matrix for MVDR beamforming

texts:
I am trying to calculate the covariance matrix that is required for the calculation of an MVDR beamformer. I am getting confused as to how to calculate it. I have an array of 3 microphones each with a vector of speech signals with N samples in the time domain. I apply an STFT to the signals giving a N x M matrix, where N represents the time domain or a single frame and M is the frequency bins. The covariance matrix can be described as:

$$
R_{xx} =E[xx^{H}]
$$

Thus, for a 3 microphone system we get:

$$      R_{xx} = E\left(
        \begin{matrix}
        x_1x_1^* & x_1x_2^* & x_1x_3^* \\
        x_2x_1^* & x_2x_2^* & x_2x_3^* \\
        x_2x_1^* & x_3x_2^* & x_3x_3^* \\
        \end{matrix} \right)
$$

Keeping in mind that the system will be done in real-time, would the expectation function apply a time-average across past and current signals (i.e. for the first frame divide by one, second divide by two and so on) or just current signals (i.e. divided by one)? And so, would each frequency bin in each frame have its own $R_{xx}$ giving N x M $R_{xx}$ values?

Thanks!


texts:
according to
Estimating acquisition noise P. 2

Estimating the mean and the standard deviation for each pixel is calculated  as below:

$$
     \begin{equation*}
      i, j = 0,\dots,N-1 \\
      \bar{I}(i, j) = \frac{1}{N} \sum_{k=0}^{N-1}  I_k(i,j)\\
      \sigma(i, j) = \left [ \frac{1}{N-1} \sum_{k=0}^{N-1}(I_k(i, j) - \bar{I}(i,j)^2) \right]^{1/2}
     \end{equation*}
$$

My Question

What is $I_k$ ?
What should mean of each pixel be? Is it the mean of whole Image?

Thanks


texts: How to compare imaginary coherence values?

texts:
High, I am currently writing my Master's thesis on susceptibility for hypnosis. I am comparing EGG signals of lowly and highly susceptible individuals.

The imaginary coherence is given by:

$$ \text{ImCoh} = \frac{\mathcal{Im}\{E\{S_{xy}\}\}}{ \sqrt{ E\{S_{xx}\} E\{S_{yy}\}} }$$

I found out that I can calculate the coherence like the scipy.signal.coherence function with the following code:

sfarrA, Pxx = scipy.signal.welch(ex_sig1[0], fs=128)
sfarrB, Pyy = scipy.signal.welch(ex_sig1[9], fs=128)
sfarrAB, Pxy = scipy.signal.csd(ex_sig1[0],ex_sig1[9],fs=128) ```

coherence = abs(Pxy)**2/(Pxx*Pyy)

Next I calculate the imaginary coherence using following code:

Icoh = np.imag(Pxy)/np.sqrt(Pxx*Pyy)

As a result I get an array of negative and non-negative values.

array([ 0.        ,  0.15618352,  0.09958203,  0.00483643,  0.02767622,
        0.07590722,  0.06308123,  0.05571067,  0.0312004 , -0.06364436,
       -0.15178242, -0.0313257 ,  0.01023862, -0.02718701, -0.03659143,
        0.04386298, -0.05511231, -0.05728281, -0.03789811, -0.07106073,
       -0.34574961, -0.27614587, -0.27608529, -0.37483751, -0.34951999,

How could I compare these values in a meaningful way?

Any help would be greatly appreciated.

Many thanks


texts: Can addition of two band limited signals create aliasing?

texts: Can mixing/adding two band-limited signals create any frequencies above Nyquist?


texts:
Suppose that we have some peaks of a signal and we know the time lag between dominant peaks. How is it possible to remove the extra peaks by applying time lag ?
Lets say it in MATLAB

[corr]=xcorr(signal,signal);
[val1 ind1]=max(corr);
corr(ind1:length(corr))=0; % symetric properties of correlation
[val2 ind2]=max(corr);
lag=ind1-ind2;
[pks ind]=findpeaks(signal);

Actually, I would like to remove the spurious peaks by using time-lag but maybe some dominant peaks a little less or more than exact time-lag.


texts: Separation of many signals all containing 1 common component

texts:
I have been struggling with a research question and would appreciate any pointers (or full solutions).
I have a data set of many (+/- 50k) observed - correlated - signals. These signals are known to contain 3 main components: S_total = S1 + S2 + noise.
S1 is common to all signals, S_total is observed.
For one signal, this could be visualized as:

A few signal together could look like this:

Note that the variation between these - simulated - signals is only coming from S2.
What method / algorithm would you recommend to obtain S1, which is common to all signals observed?
Many thanks in advance,
Bart


texts:
I'm making a guitar tuner for iOS.Due to the fact I'm Beginner I'm struggling a bit to gather all the resources and information about it. I know the theory like ( correct me if I'm wrong ) :-

First I need to get the input from microphone.
Then need to apply apply FFT algorithm to get the frequency
Then compare the frequency with the fundamental frequency of notes.

according to my research the fundamental frequency of the guitar are measured as:-

String  Frequency
1 (E)   329.63 Hz
2 (B)   246.94 Hz
3 (G)   196.00 Hz
4 (D)   146.83 Hz
5 (A)   110.00 Hz
6 (E)   82.41  Hz

and from an test app example I'm getting the frequency:-

String  Frequency
1 (E)   333.02  Hz
2 (B)   247.60  Hz
3 (G)   398-193 Hz    (398 when start and 193 when end)
4 (D)   290-150 Hz    (290 when start and 150 when end)
5 (A)   333-215 Hz    (333 when start and 215 when end)
6 (E)   247-161 Hz    (247 when start and 161 when end)

  One thing to notice the example was showing the Max Frequency.

So what I'm asking here is there anybody who have implemented this before
and give me the direction and some detail information about the topic so In future nobody have to research more and more to gather all the resources.I need help with the following topics:-

What is the best possible way to get the accurate frequency of guitar string.
Which type of frequency need to tune a guitar like maxFrequency, is there a role of altitude, magnitude (i have less information about this topic).
What to do after getting the right frequency.

Any help would be truly appreciated and my apologies in advance if i'm asking something stupid .


texts: Kalman Filter - Non Linear Measurement Model

texts:
I'm new to Kalman filters and estimation in general. I'm running a simple test of an EKF to check my understanding, but I'm getting some odd results with a particular case.

Given a state vector:
$$ \textbf{x} =
\left(\begin{matrix}
p_x \\
p_y \\
p_z \\
\end{matrix}\right)
$$

If I have a sensor that measures these values directly, i.e.:
$$ \textbf{z} =
\left(\begin{matrix}
p_x \\
p_y \\
p_z \\
\end{matrix}\right)
$$

Then I get a converging solution. However, if I decide that my sensor actually measures the square of each value, i.e.
$$ \textbf{z} =
\left(\begin{matrix}
p_x^2 \\
p_y^2 \\
p_z^2 \\
\end{matrix}\right)
$$

Then following the process on the EKF wiki page causes my solution to diverge very rapidly. This seems to be irrespective of the dynamic model of the system, but entirely dependent on whether what the sensor is measuring is linear. Is there a standard approach to follow in this situation?

EDIT: Turns out it was a silly mistake on my part. I just implemented the filter incorrectly. The EKF was valid in this case.


texts: Frequency Domain Equalization [Sources request]

texts:
When I read papers, I often see the term "Frequency Domain Equalizer".

I do not understand what a "frequency domain equalizer" is.

Do they mean equalizers like MMSE, LMMSE?

Or is there an equalizer called FDE?


texts: What's the advantage/disadvantage of oversampling followed by decimation, verses sampling at the correct rate to begin with?

texts: I read somewhere that oversampling somehow reduces sampling noise?  How to quantify the reduction in noise?  and why would it reduce noise if you are decimating the oversampled signal.  What's the advantage of oversampling followed by decimation, verses just sampling the signal at the decimation rate to begin with?  Or is the advantage moot, and its generally better to just sample at the decimation rate to begin with and forgo the oversampling.  would there be a case where one is better than the other, and vica versa?


texts: Bit-reversal equivalence on IFFT (radix-2 Cooley-Tukey)

texts:
From everything that I read on the internet, I understood that when performing FFT or IFFT, one should perform bit-reversal to obtain a result in natural order (Assuming input is also in natural order). I also understood that this bit-reversal can be performed before or after the FFT or IFFT processing.

However, when interverting the order of operations, from :

bit_reversal();
IFFT();

to :

IFFT();
bit_reversal();

I don't get the same results. I would expect some difference due to rounding (I am working in fixed-precision here). But I get really completely different results.

Unfortunately, I am not allowed to share the code here, so my question is not about the implementation of the bit-reversal or IFFT, but more about fundamentals that I may have missed : can we operate bit-reversal equivalently before or after the radix-2 IFFT, without changing the radix-2 algorithm ?

Thanks in advance.


texts: How to deal with signal not equally spaced in time when performing FFT?

texts:
I wonder to know what is the best way to handle not equally spaced in time signal when performing FFT ?

I guess it depends on the signal itself. I work with signal with about 1000 - 5000 samples and they can contain between 0 to 20% missing samples. And I am interested in both main frequencies AND amplitudes.

Currently I do a basic linear interpolation to "fill" missing samples. I heard somewhere that we could actually "fill" missing samples with NaN (I didn't try yet).

What do you think about these two options ? Advantages and disadvantages ?


texts:
I have the following block diagram

I know that it's solution is:

$$\dfrac{10 K_p}{Js^2+(B+2K_pK_v)s+10K_p} $$

But when I try to derive the solution with reduction rules I get different result.

Assumptions:

$5$ and $K_vs$ are in parallel: $(5+K_vs)$

$K_p \cdot \dfrac{2}{(Js+B)s}$

result: $$\dfrac{K_p \cdot \dfrac{2}{(Js+B)s}}{1 + K_p \cdot \dfrac{2}{(Js+B)s} + (5+K_vs)}$$

What am I doing wrong?


texts: How to compose a Discrete Prolate Spheroid (DPSS) dictionary?

texts:
I have a model of signal as
$$
Y=AX + N
$$
where $Y$ is received data in a linear array, $A$ is steering matrix, $X$ is data of sources and $N$ is noise. If $A$ has the form of $A=\exp(\alpha \sin(\theta))$ where $\alpha$ can be consider as a constant, $\theta =[0,180]$. How can I compose a dictionary of $A$ based on DPSS-to solve an estimation problem?


texts:
As part of an audio peaking filter calculation, I need to implement hyperbolic sine ("sinh") function in a fixed-point arithmetic (fixed-point DSP processor).

What is the proper way to do this?


texts: Why are zero values added in the FFT of a concatenated noise signal?

texts:
After concatenating an 11-second *.WAV (sample rate 44100 Hz) file 27 times (just by gluing them head-to-tail), I obtained a ~300-sec noise signal. I used MATLAB(2015b) and the FFT function on both the original file and the concatenated signal. The FFT of the concatenated signal was, as expected, pretty similar to the original. However, for every y-axis point (spectral power) in the concatenated FFT, there were 2 zeros (not precisely zero, but values in the order of x*10-19, so about 18 orders of magnitude below the 'real' FFT values). The original noise file was specifically designed to allow for looping. Below a close up of the concatenated FFT showing the added zeros:

Hence, the spectral power plot of the concatenated signal looked like the original, but every frequency power point was preceded and followed by two zeros. How can this be? Does it have to do with the fact that lengthening a signal increases the resolution of the FFT and that frequencies are analyzed in the concatenated signal that are not present in the original signal due to resolution restrictions? Or is this rather an error in my script? Other spectra showed up fine, though.


texts: Fast Normalized Cross-Correlation for Stream-Data

texts:
I read some papers about fast normalized cross-correlation techniques. Most of them are used to identify a template in an original image.

I was not able to find information about input-data in a stream-based fashion:
Assume, that there are two audio signals a and b both of the same length. Now, I want to find the position where both signals are best aligned.

It should possible to use normalized cross-correlation to calculate this.

Example (1):

a = {1,2,3,4}; b = {1,2,3,4}

Using normalized cross-correlation the result is given as:

{-1.8, 1.2, 1, 4, 1, 1.2, -1.8}

Now, let us assume that the signals will be extended, but we have a fixed queue-size for a and b.

Example (2):

a = {2,3,4,7}; b = {2,3,4,4}

Again, we use normalized cross-correlation to get our results:

{-0.967, -1.450, -0.161, 3.223, 2.256, -0.483, -2.417}

Normalized cross-correlation makes use of mean and variance of a and b and for me it seems, that there is no option to get a stream-based version of the normalized cross-correlation without recalculate the complete normalized cross-correlation. For sure, variance and mean can be calculated in a stream-based fashion.
$$
\frac{1}{n} \sum \frac{(f(x,y)- \overline{f})(t(x,y)- \overline{t}) }{\sigma_{f}\sigma_{t}}
$$
Is it correct, that because of each calculation consists of mean and variance, it is not possible to have a stream-based normalized cross-correlation, because they change if a and b change?!

Would be happy to have a great discussion about this.


texts:
I have trained a machine learning model to recognize a specific audio event from mel spectrograms, and it works pretty well on a hold-out test dataset. I now want to run the model on a device and feed it with real audio but I have a problem that I don't now how to solve.

The problem is that when I'm testing my audio recognizer by playing audio clips it fails to recognize most of them, even if it can correctly recognize the audio in the original file. When I compare the spectrogram from the original audio file with the spectrogram computed by the recording of my device, I can see some differences:

Some of the differences has to do with the fact that the real time algorithm is doing the time-domain segmentation a little bit differently so that it "misses" the first few milliseconds of the event. But I have verified that this is not a big problem for the machine learning model. So the reason for why it fails seems to be because the loudspeaker and/or the microphone messes with the audio so that some frequencies are attenuated. The audio file has a 50kHz sample rate, which I think my high-end loudspeaker should be able to handle. My device samples the audio at 44.1kHz and the mel spectrogram uses frequencies up to 22.05kHz.

How can I make my feature extractor more robust? So that the spectrogram looks more similar irrespective of which microphone is used for recording or if I use a loudspeaker or if the audio event is actually happening in front of the microphone.


texts: Why do the in-phase and quadrature components of a narrowband signal each have half the bandwidth of the signal?


texts:
I'm using a very simple 1st order Butterworth Filter shown in Matlab code:

order = 1;
cutOff = 0.1;
[b, a] = butter(order, (2*cutOff)/SampleRate, 'high');

So from the Matlab documentation of butter, it is not hard to know the transfer function is:

$$H(z) = \frac{b_1 + b_2z^{-1}}{a_1 + a_2z^{-1}}$$

The question is, if I have a signal $x$ with variance $v$, after the signal $x$ is filtered by the above high-pass filter, what is the variance of the output signal? Will it change? How can I calculate it? And is the filter linear?


texts: Finding the impulse response given response to another signal

texts:
I was trying to solve this question :

I respresented $x(t) = u(t+1)-u(t-1)$
writing the convolution as $[u(t+1)-u(t-1)]*h(t) = y(t)$
I then used the property of differentiation to convert from the step to the impulse function :
$[\delta(t+1) -\delta(t-1)]*h(t) = y'(t)$
$\implies h(t+1) - h(t-1) = y'(t)$
However , I am now stuck here. Moreover , the fact that the width of $y(t)$ is the same as that of $x(t)$ , bugs me a little.
How can I find the impulse response $h(t)$ using the given info ?


texts:
I'm quite a newbie in dsp domain, and there's something I can't understand clearly.

I'm doing a IFDT of a spectrum regularly spaced with complex values : the spectrum is a complex conjugate symetric frequency domain sequence so the output of the IFDT should be real.
It is when using a sampling rate of 50MHz with 2048 samples.
But when I change the number of samples (for exemple 1024), it does not work anymore.

Can someone explain why it won't work ? I'm clearly missing something...


texts: Is the Channel shortening similar to Time domain Equalizer (TEQ)

texts:
I used to read about Time domain equalizer (TEQ) in many lectures, i.e here , But in many other papers,  that's called  Channel shortening i.e here .

However, in the a two above links, the cost functions is almost similar, but in first paper it's called Time domain equalizer, and in the second, it's channel shortening.

What's the difference between them (I mean between Time domain equalizer TEQ and Channel Shortening)?  And why do we call it sometimes TEQ and another time channel shortening?

Thank you.


texts: Self study question about Padé Approximation and transfer functions

texts:
Given a signal $x[n] = [1, 2, 3, 4, 5]$. How many transfer functions can be found with Padé approximation which have a causal impulse response and start with those five samples. How many of them are stable?

Found this question and I'm not sure if I understood it right. I'd say none of them is stable since Padé can't guarantee stable solutions. The first five samples will be perfectly recreated in the impulse response. Therefor $h(z)$ between $0-4$ looks like $x[n]$. But what about how many? Are there multiple solutions? Is this question badly formulated?

A further question is then: Find the Padé model with exaclty one pole. One pole? Isn't that an all pole model? Or is there some kind of rule like always the same amount of poles and zeros? Can't solve it with the Padé euqations.


texts: How to Draw a Line at Non Integer Pixel Coordinates on an Image

texts:
I have lines (or segments really) I want to draw into a new image I.
Let's define the intensity of the pixel at row = r and col = c as I(r,c).  However, the lines are defined by 2 points whose coordinates are not integers.

For example, one line might be from point A (with row = 2.2 and col = 3.5) and to point B with (row = 10.8 and col = 30.4).

What is the algorithm or method for filling in the intensity values of image I at the integer rows and cols?


texts: What is the difference between closed form and opened form in z transform

texts:
I am new in DSP and I need to have a clear udrastanding about below example.
Consider the below $h[n]$ as the impulse response of a system:

\begin{equation}
h[n] = a^n(u[n] - u[n-N])
\end{equation}

so the opened z transform can be expressed as :
\begin{equation}
H(z) = 1 + a^1 z^{-1} + a^2 z^{-2} + \dots + a^{N-1} z^{-(N-1)}
\end{equation}
and its closed form is expressed like :
\begin{equation}
H(z) = \frac{1}{z^{N-1}} \frac{z^N-a^N}{z-a}
\end{equation}
defenetly, closed form is better than the opened form for hardware implementation because as below, the closed form can be written in the time domain like this:
\begin{equation}
y[n] = x[n] - a^N x[n-N] + a y[n-1]
\end{equation}

Now the only difference between the closed and opened form is a zero and a pole at
$z=a$. Actually, it seems to me, by adding zero and pole, we obtained a better representation. But my question is,

What is the drawback of adding zero and pole to the system at $z=a$.

How can we analyze the benefits and drawbacks between closed and opened form?

I would appreciate any help you can provide.


texts: Gaussian wavelet generation of a given frequency

texts:
How can I generate a Gaussian wavelet (time domain) with a given central frequency.

I mean, if I take the Fourier Transform then its spectrum should be around that given central frequency. For example the peak of spectrum is 20Hz and its side lobes becomes nearly zero around 20±10 .

I carried out a few coding exercise but the fourier transform of the gaussian wavelet was always centered around 0 Hz regardless of shape of gaussian spectrum.

Code and explanation are here Link


texts:
I have a certain signal that has a linear trend that should be removed:

The problem is that the signal needs to be processed in real time(sampling rate is 20Hz).
If I only apply the detrending per window of say, 20 samples and then concatenate detrending results, the signal is disrupted.

I wonder if there is some efficient technique to detrend the data per window in real time and then concatenate it so the final signal will make sense? Or should I just detrend every new incoming data with the older data which is far from being efficient.

I would like to note that taking the entire signal and processing it offline using MATLAB's detrend works correctly. The problem occurs when breaking the signal into certain fixed size windows and then applying detrend on each.

Thanks,
Mike


texts: Confirmation on how to calculate phase differences at every frequency point between two wideband signals?

texts:
I have a problem that I thought was going to be simple, but it has become surprisingly stubborn and I question my method...

I have used the method (described below), but I wanted to confirm that this was the right approach, and see if there was another method I could use.

What I basically have is two wideband signals, (on the order of 100 Mhz BW), in quadrature. So I have the I and Q of those two signals. What I am trying to find is the phase difference between those two signals as a function of frequency.

The method I am currently trying is the following:

I take my quadrature complex signals, and window each of them with a hamming window 3 times over.
I take the FFT of the windowed complex signals, at some arbitrary length, say 2048. Call the results $S_1(w)$ and $S_2(w)$.
I multiply $S_1(w)$ by the conjugate of $S_2(w)$.
I take the argument out of the complex exponential vector from the above step (inverse tangent of imag over real). This gives me the phase differences between the two signals as a function of frequency.

My questions on this method:

Is this the right approach for this simple type of problem?
When I tried this method for a simulation with 1, 2 or a couple tones, I can get good results. I find the index of where my tones exist, (look at the abs of the FFT), I take that index, and look at the $\Delta$phase value from step 4 above. This gives me the right phase.

However, when I look the $\Delta$phase values for frequency indicies where said those frequencies dont exist, I get wild answers. Why is this? I realize those frequencies dont exist but shouldn't the phases just be 0 - 0 = 0 in that case?

Based on this (with a couple tones), can I readily extend this to a wideband signal and also expect good results?
My last question regards the nature of the phase delay here - physically what is happening is that one of the signals (optical) is going through a device that delays different frequencies differently. Of course the delays are manifested as phase offsets between all frequencies... but at the same time there is a time delay associated with this as well. Even if we were to measure the $\Delta$phase between the two signals properly, wouldnt it still be impossible to measure the time delay between them at those frequencies? I can have a time delay of 10 wavelenghts, but a phase offset of 0.

Thoughts?

Thanks in advance!

EDIT FROM COMMENT FEEDBACKS:

I think I am clearer now as to what my question really is:  For a case where you have well separated tones and high SNR, my method can work.

However, what to do in the wideband case? In this sense, suppose $s_1[n]$ and $s_2[n]$ both contain two tones that are very close to each other in the frequency space. But your FFT can only give you so much frequency resolution. So instead of applying my method to two peaks, you only have one peak to go on. Now what?

What relation would the $\Delta$phase of that peak have to the actual $\Delta$phase of the  'true' peaks? Is it their average?
One method I am thinking of is correlating the two signals in the time domain, and working backwards to get the phase offsets per frequency taking advantage of the fourier translation property, although part of me thinks it might not work since there is no new information, just looking at the problem from a different domain.
Another thought I had for the wideband case would be to put the whole thing through many narrow band filters, and then apply my original method. Is there a name for something like this?

Thanks again


texts: Confusion regarding sinusoidal synthesis formula,

texts:
I am reading signal processing first by McClellan

As shown in attached snapshot, we have formula of sinusoidal synthesis

I am confused ,why in expression of x(t), there is a '2' in first line immediately after summation sign but that'2' disappear in 2nd line?


texts: Filtering signal with added noise (frequencies known)

texts:
I have been assigned a project to take a given noisy sound file, with added broadband and narrowband noise. The goal is to reach a certain mean-squared error between the filtered signal and the original, non-noisy signal.

Here is my code for calculating the mean squared error (in MATLAB)

%  let x be the original signal, xn be the noisy signal, and xnf be the filtered
%   noisy signal. Let the prefix 'G' indicate a power spectral density,
%   a lowercase 'x' indicate time domain, and uppercase 'X' indicate FFT.
%   Compute the cross correlation between the original signal and the
%   filtered noisy signal
Rxxnf = N*real(ifft(ifftshift(Gxxnf))) ;
%   Find the number of samples of delay, I, for peak cross correlation
    [Rmax,I] = max(Rxxnf) ;
%   Circularly shift the original signal by I samples to make it have a
%   maximum correlation to the filtered noisy signal.  Call the shifted
%   signal "xs"
 xs = [x(end-I:end);x(1:end-I-1)] ;
%   Compute the mean squared error between the original signal (delayed)
%   and the filtered noisy signal.  Compute only for the samples that occur
%   after the delay.
   mse = mean((xnf(I:end)-xs(I:end)).^2) ;

The sound clip is about 10 seconds long, and is of a man speaking.

The locations of the narrowband noise are specified as ranges from

967.884 Hz to 1069.95 Hz
715.365 Hz to 810.994 Hz

Looking at the FFT of the original signal (left) and the noisy signal (right), we can see that there is indeed a significant amount of broadband noise.

To inspect this further, I look at the power spectral density.

I can confirm that the  broadband noise is the outermost section of the noisy PSD and that narrowband noise has indeed been added at the indicated frequencies.

I have applied a lowpass filter (6th order butterworth, $f_c=700\textrm{ Hz}$ (the broadband noise is actually introduced closer to $1000\textrm{ Hz}$, but $700\textrm{ Hz}$ gives the best mean-squared error), and 2 bandstops at the given frequency ranges (chebyshev type II, both second order, one with a passband ripple of $3\textrm{ dB}$ (low ripple gives best MSE), and one with a passband ripple of 18 (best MSE).

As of right now, my MSE is not as low as my professor's claimed MSE for this particular signal (he indicates a mean-squared error of .00143, I have only achieved .0023)

My questions are:

are there any other filter types I should be using rather than butterworth/chebyshev?
Should I place any other filters anywhere to perhaps mitigate ripple effects?

I very much appreciate the help.


texts: best similarity measure for audio classification

texts:
I'm doing a project on audio classification, mainly focused on Bird sounds.

I would like to know which is the best similarity measure that I can adopt to check test signal with training data.

I intent to do the project in matlab.


texts: Can we broadcast two messages using NOMA?

texts:
Is it possible for a transmitter to broadcast two messages to a set of receivers at the same time using superposition coding (Non-Orthogonal Multiple Access - NOMA)?

A transmitter has two messages $m_1$ and $m_2$ that wants to broadcast to a set of $n$ receivers. I think this not possible even for $n=2$. If the transmitter broadcasts $P_1m_1 + P_2m_2$ with transmit power factors $P_1$ and $P_2$, then only the receiver with higher channel gain can decode both $m_1$ and $m_2$. Is this true?


texts: What kind of effect is this and what is it caused by?

texts:
In the screenshot below I have two audio signals. The upper one is the audio signal that is created  with python and then played, whilst the lower one is the "same" signal but recorded via microphone on an android device. Tried my best at aligning them in audacity but I am unsure where the "fade-in" comes from. The same thing happens on the other end of the signal, but in reverse of course.

Also, there is no other noise in the recording, as I am recording from the audio output relayed to the input and not directly from the speakers (so no background noise).
Don't know which terms to search for. Any help would be much appreciated.


texts: How to solve this question by using convolution?

texts:
I had previously posted this question but my method was wrong in it.
Two signals are given:

$$x[n]=2^nu[n]$$
$$h[n]=u[n]$$

I have to find $y[n]$. I have spent hours on this question but I am unable to solve it. Kindly please tell me the way to solve it.


texts: Strange results with all pass compensation

texts:
I am new to DSP and, as a learning exercise, I am trying to design an all-pass filter in order to compensate the phase of an input impulse response (heqm.mat, which you can download from here).
The MATLAB code is the following; for the phase compensation I am using the time-reversed technique, also called matched-filter:
load heqm.mat; % Load the heqm vector
HEQM      = fft(heqm);
% Design the all-pass filter
HEQMP     = exp(1j.*angle(HEQM));
tmp       = ifft(HEQMP,'symmetric');
% Time-reversed technique
eqpr      = tmp(end:-1:1);
% Figure on the left
figure; plot(unwrap(angle(fft(heqm))));
hold on; plot(unwrap(angle(fft(eqpr))));
% Convolution of the original impulse response with the all-pass filter
o = conv(heqm,eqpr);
% Figure on the right
figure; plot(unwrap(angle(fft(o))));

As you can see from the left picture, the phase of the all-pass filter is the opposite of the original one. But if I convolve the two vectors, why the resulting phase (right figure) isn't flat as it should be? What am I missing here?


texts: What does MATLAB display when plotting a sound signal?

texts:
When  you load the speech file 'mtlb' in MATLAB (which is by the way the voice of a woman saying the word 'MATLAB'), and plot the data, you get the following plot:

My question is this: what does the y-axis represent? What is MATLAB actually plotting here? Thanks!