texts: What kind of discrete convolution is this?

texts:
Let x[n] and y[n] be defined for n=0,1, ..., N-1.

This (for example with Python's scipy fft / ifft)

a = ifft (fft(x) * fft(y))

should give the convolution

$$a[k] = (x * y)[k] = \sum_{n=0}^{N-1} x[n] y[k-n]$$

But is it done with

y[i] = 0 if $i<0$  ?

or

y[i] = y[N+i] if $i<0$ ?  (ex : y[-3] = y[N-3] )     (in this case do we call it circular convolution?)

Note : More generally, I have good books about signal processing / Fourier transform in general, but they are not really handy handbooks for discrete finite signals (x[n] for  n=0,1, ..., N-1). Do you have a good reference of handbook for this case (discrete finite signals) ?


texts: How to interpret correlation angle from frequency domain

texts:
Since the coefficient of one frequency could be seen the length of the signal projects on that frequency basis.

And the length is
$$|a|\cos\theta= \frac{a\cdot b}{|b|}   $$

If I divide the coefficient by $|a|$, how to interpret the correlation cosine angle from frequency domain?


texts: What is the phase of the Fourier Transform of $h(n)={\frac13,\frac13,\frac13}$ (zero position is at 2nd element)

texts:
I am using Digital Signal Processing Principles, Algorithm and Applications,  4th edition written by Proakis and Manolakis. Here is s a part that I don't understand.

$h(n)={\frac13,\frac13,\frac13}$. The zero position is at the 2nd element.
Hence, its Discrete-Time Fourier Transform (DTFT) is: $$H(\omega)= \frac13(e^{j\omega}+1+e^{-j\omega}) = \frac13(1+2\cos{\omega})\,.$$

Since only the real part exists, I assume that the phase $\tan^{-1}\frac{H_I(\omega)}{H_R(\omega)} = 0$, where $H_I(\omega)$ and $H_R(\omega)$ are the imaginary and real part, respectively. However the book says that the phase is $0$ for $0\leq\omega\leq\frac{2\pi}{3}$ and $\pi$ for $\frac{2\pi}{3}\leq\omega\leq\pi$.

Sorry that I'm not good at Latex.


texts: How can I properly detect a ball moving in a film capture?

texts:
I have problem in detecting a ball moving in a flume (image below).

since converting the frames into grayscale format didn't work I converted it into rgb and make it more saturated. problems are:
1- in the first frames which ball is in the right position, light disturbed the ball detection (image subtraction didn't work well)

So I did it without subtraction.

this is my problem in detecting the ball in the first frames. in the rest frames the problem is its shadow and the white dots on the ball. I can't see the exact ball to find it's exact position in each frame (its shadow lead to some uncertainties)

anyway, I'm using matlab. I shared the code below (subtraction is not working well because of light, note that ball moves since water make it move in the flume).

Any Ideas about how can I detect the ball more precisely? handle the light in the first frames? handle the ball's shadow?

thanks in advance.

% the code starts

I=imread('frame3860.png');
I=imadjust(I,stretchlim(I));
I=lensdistort(I,-0.3);
I2=imcrop(I, [263.5 85.5 1468 940]);
I2=imrotate(I2,-1,'bilinear','crop');
I2=rgb2hsv(I2);
I2(:,:,2)=I2(:,:,2)*1.9;
I3=I2(:,:,1)>0.57 & I2(:,:,1)<0.75;
I4=bwareaopen(I3,200);
I4=medfilt2(I4);
I4=imfill(I4,'holes');
str=strel('disk',1);
I5=imdilate(I4,str);
I5=I5-I4;
e = edge(I5, 'canny');
radii = 30:1:70;
h = circle_hough(e, radii, 'same', 'normalise');
[~,maxIndex] = max(h(:)); %find the global maximum of your accumulator array
[i,j,k] = ind2sub(size(h), maxIndex);
radius = radii(k);
center.x = j;
center.y = i;
N = 200;
theta=linspace(0,2*pi,N);
rho=ones(1,N)*radius;
[X,Y] = pol2cart(theta,rho);
cent(1,:)=[j,i]; % centroid of the ball in each frame

% code finished

If you have any idea for a better ball detection in my situation I will be  grateful.

P.S. I removed the images since I didn't have the reputation needed to post more than 2 images

thanks in advance


texts:
I have series of amplitudes of a signal recieved at center frequency of 433.9MHz, sampling at 250000 Hz. I can see the information received is the packet I am after, just sampled too high. I believe the transmitter transmits at 3918bps. What is the best method to downsample from 250000 Hz to 3918 bps?

I have tried using a butterworth filter, which doesn't seem to work and I don't understand why. I've tried picking every 63 samples (because 250000/3918 = ~63), but didn't seem to work either. Please bare in mind I am a bit of a newbie regarding DSP.


texts: What is the successor to the GNURadio liveCD?

texts:
I used to burn and run the GNU Radio livecd - it was a very nice package that gave me a working gnuradio environment - with drivers and tools - without doing a lot of sysadmin or package management, etc.
Now, however, the GNU Radio liveCD is no longer supported nor recommended and, while still available from some places, the last release is about four years old.
What is the next best option if I want to quickly get a gnuradio environment up and running with tools like gqrx and I don't want to spend a week working through dependencies ?
Thank you.


texts: Design of discrete time second order all pass filter

texts:
I want to design a second order  digital all pass filter with the transfer function given by
$$H(z)=\frac{a_2 +a_1 z^{-1}+z^{-2}}{1+a_1 z^{-1}+a_2 z^{-2}}$$
The input to this filter has frequencies ranging from $1$ to $8\ \rm kHz$. The sampling frequency is $20\ \rm kHz$.
How to design the filter coefficients $a_1, a_2$ such that for a particular input frequency, I have to apply 180 degree phase shift? The frequency which is given a 180 degree phase shift will vary. So, how to vary the filter coefficients based on which frequency has to be phase shifted by 180 degree?


texts: Equivalence of difference equation with system of difference equations

texts:
I have this system and tried to write down its difference equation as follows

$$    y[n] = b_0v[n] + b_1v[n-1]\\
    v[n] = x[n] - a_1v[n-1]   $$

I was asked to prove that the system can be described by the following difference equation

$$    y[n] + a_1y[n − 1] = b_0x[n] + b_1x[n − 1]$$

Please help me how to prove it.


texts: Gaussian pyramid: why needs the image to be downsampled

texts:
In the context of a gaussian pyramid, why is the image downsampled separately although the numbers of pixels are decreased through smoothing already?
In my understanding, the number of pixels will decrease due to the smoothing process already, if a gaussian filter is used within the image boundaries. An image with 8 x 5 pixels smoothed with a gaussian filter of 3 x 2 pixels will create a 6 x 4 pixel image if the filter can only be slided within the boundaries of the original image without crossing the edges.

So if the constraint of not crossing the edges of the image is applied, why would one use additional downsampling?


texts: How to store the signal in a variable and plot the signal in MATLAB

texts:
I am trying to store the signal I have created in a single variable and graph the signal in two subplots in one single graph, one for left channel, one for right channel.

Fs = 44100;
N = 400
t = (1:N)*(1/Fs);
freq = 500;

sound_array = zeros(N,2);
%starts N x 2 matrix
sound_array(:,1) = sin(2*pi*freq*t);
sound_array(:,2) = sin(2*pi*freq*t);

%Left and Right channels
sound_L = sound_array;
sound_L(:,2) = 0;
sound_R = sound_array;
sound_R(:,1) = 0;
plot(sound_array);

R = 400;
beep = [sound_array; zeros(R, 2)];

%Sound 1
sound(sound_L, Fs);
pause(1);
%Sound 2
sound(sound_R, Fs);
%Play Consecutive
sound(sound_array, Fs);


texts: Design narrow bandpass filter for signal with high sampling rate

texts:
I have to bandpass-filter a signal which has been sampled with 4000 Hz.

Only frequencies around 15 Hz shall remain after filtering (let's say in a band between 10 Hz and 20 Hz; the narrower the better).

I have several questions here:

What is the recommended way to perform tasks like this?
What bandpass filter is suitable?
Do I need to downsample the signal before filtering. If yes, what type of additional low pass filter should I use to avoid aliasing?

In my case, it is important that the phase of the original signal will not be distorted or shifted in the whole process. Furthermore, the filter step response should have close to zero overshoot. The settling time could therefore be a little larger.

I do not want to perform a sophisticated filter design which delivers optimal results for exactly this specific signal. I am more interested in a general solution that delivers solid results using standard IIR or FIR filters (e. g. as available in Python's SciPy library) which could be reused afterwards for similar tasks.

UPDATE:

After the answer of @MarcusMüller and the provided links, I basically tried every FIR filter design method available in Python's SciPy library and deciced to go with remez. I developed the following code which may be used for further discussion:

import math
import matplotlib.pyplot as plt
import numpy as np
from scipy.signal import lfilter, remez

F_test = 20.0
duration = 10.0
fs = 8000
samples = int(fs*duration)
t = np.arange(samples) / fs
signal_test = (5.0 * t * np.sin(2.0*np.pi*F_test*t)) + (0.5 * np.sin(2.0*np.pi*5.0*t)) + (0.5 * np.sin(2.0*np.pi*100.0*t))

#design filter

ntaps = 5000
edges = [0, F_test - 5.0, F_test - 2.5, F_test + 2.5, F_test + 5.0, 0.5 * fs]
taps = remez(ntaps, edges, [0, 1, 0], Hz=fs, maxiter=2500)

#apply filter
signal_test_filtered = lfilter(taps, 1, signal_test)

#create plot
fig = plt.figure()
ax0 = fig.add_subplot(111)
ax0.plot(t, signal_test_filtered, label='signal_test_filtered')
ax0.set_xlabel("time [s]")
ax0.legend()
fig.show()


texts:
My book with answers to proakis and manolakis book says that $F(x,n)=x(-n + 2)$ is time invariant while my calculation shows that it is time variant. I have repeated them many times and I am mostly convinced that book is wrong. Can somebody please check my calculations?

First we shift output signal:

$$F(x, n - k) = x(-(n-k) + 2) = x(-n+k+2)$$

now lets put $x'(n) = x(n-k)$ and check what happen when we input shifted signal:

$$F(x', n) = x'(-n+2)=x(-n+2-k)=x(-n-k+2)$$

since $F(x, n - k) \neq F(x', n)$ this system is time variant. Is my method correct?


texts: Applying a mathematical function to a 16-bit Python array

texts:
I am using scipy.io.wavfile to read a WAV file and playing it with pyaudio. Before the numpy array is written to a string and output through pyaudio, I would like to be able to apply some sort of function (sin, cos, etc.) to the 16-bit numpy array.

Does anyone have any experience doing something like this? Should I convert to float, apply the function to the float data, then convert back? I am also aware of fixed point arithmetic, but not sure if I can apply it here.


texts:
Is there a function, in continuous or (primarily) discrete time, such that convolving with it preserves the input's energy? For $x$ that is finite-valued,  finitely supported in time, and:

Real or complex?
Real?
Real and non-negative?

where "energy" refers to L2^2:
$$
\text{E}(f(x)) = \int_{a}^{b} |f(x)|^2 dx \ \ \text{or}\
 \ \sum_{n=a}^{b} |f[n]|^2
$$
and following must hold
$$
\text{E}(h \star x) = \text{E}(x)
$$
where $\star$ is linear or circular convolution (either works for me). If  there's no such kernel, can this be proven?

I'm aware of the trivial case of the delta; $h$ should be such that $x$ is transformed non-trivially  (not just a copy, conjugate, rescaling, sign flip, etc).
Idea: we can satisfy Plancherel's formulation per Envidia's answer with a distribution, example being the Dirac delta, but unsure there'd be something one could implement finitely.


texts:
I understood that SNR is the ratio of signal power to the noise power. In terms of images, how the original image is affected by the added noise. In PSNR, we take the square of the peak value in the image (in case of an 8 bit image, the peak value is 255) and divide it by the mean square error. The SNR and PSNR are used to measure the quality of an image after the reconstruction. I understand that higher the SNR or PSNR, the reconstruction is good. What I don't understand is how SNR and PSNR differs in terms of their conclusion about the reconstructed image.

What the PSNR of an image concludes that the SNR of the same image can't conclude ?
Simply how the conclusion of PSNR differs from the conclusion of SNR?


texts:
I need to sample a signal of 40MHz with a sampling frequency of 120MHz at least with a 8/12 bits ADC and send them through the USB serial port of a PC. That means the serial communication bandwith should be at least 960/1440 megabits per second.

What kind of DSP should I use to satisfy this bit-rate, and what signal-processing-methods should I learn to optimize this bit-stream to use the same data simultaneously with a mainstream PC. (The bitrate seemed a bit huge for me.)

I have no experience with DSPs before. I was seeking for answers for a while but not even my teachers, or their assistants could help me for a little long period of time.

I just need a beacon to find the answers of the questions, thanks in advance.

Edit: I wrote 40Hz instead 40MHz, I am really sorry.


texts: Using unlabeled EEG data for Machine Learning

texts:
I am working on a project that is basically a game for motor-paralyzed people.

It should take an EEG signal from FP1 channel of brain and then after processing it should generate command for the game (like left,right, stop, fire, jump etc).

In the start stage of the project, I want to train my ANN on data that has already been gathered.

I am new to the signal processing field. The dataset I found over internet is unlabeled data (i.e. it is just a signal data; it does not tell about output of the specific data row).

So can somebody, who has working experience with EEG signal, tell me what I should do? More specifically, I am unable to understand the signal.


texts: Question about delayed sampled sinusoid math expression

texts:
I have been studying the digital audio processing by using the book <Designing Audio Effect Plugins in C++>.
For analog Sinusoid:
Complex Sinusoid = $e^{jωt}$
Delayed Sinusoid = $e^{jω(t−n)} = e^{jwt} * e^{-jwn}$, a delay of n seconds
For digital sampled version:
sampled complex sinusoid = $e^{jωnT}$, T is interval for each sample, n is the index of sample
I understand all above, but I got confused about the delayed sampled sinusoid which described as: $e^{jω ( nT −M )}$, M = samples of delay
But I think it should be described as $e^{jωT( n − M )}$, since the T is a constant for a fixed sample rate, n and M has the same unit.
At first I thought that maybe a typo, but the following computation parts of the book are all using the $e^{jω ( nT −M )}$ as basis.
Anyone can explain it for me?


texts:
I'm doing a work on the acoustic emission of a leak. I am using a 1.80 m cylindrical tube with water, and the leak is simulated by means of a tap that opens after I have applied negative pressure inside the tube.
I am taking this spectrogram as a test on a tube filled with water and simulating a leak. I notice that in the spectrogram there are some horizontal lines visible at some specific frequencies throughout the spectrogram, and I have no idea how to interpret those lines at some specific frequencies.
The sharp difference in the middle of the image is when the tap is turned on.


texts: How to choose a fixed adaptation step for decision feedback equalizer

texts:
How do I choose $\mu$ - a value of a step size for adaptation in a decision feedback equalizer (DFE) with adaptive reference control (ARC)? For a regular adaptive FIR filter an adaptation step depends on a variance of input, but DFE gets symbol values as input. Another problem is that DFE employs a sign-error LMS adaptation algorithm, while ARC - a simple LMS. So do $\mu$ value for $ARC$ and $h_i(k+1)$ can be the same?
The algorithm for the DFE and ARC is the following (based on this M.Sc. thesis by Mirna Hage as a reference).
With $x^{in}_k$ being a $k^{th}$ sample:
\begin{align}
x_k &= x^{in}_k - \sum_{i=1}^N h_i(k) \cdot a_{k-i}\\
a_k &=
\begin{cases}
+3 & \text{if} & x_k > 2 \cdot ARC_{k-1}\\
+1 & \text{if} & 0 \leq x_k \leq 2 \cdot ARC_{k-1}\\
-1 & \text{if} & -2 \cdot ARC_{k-1} \leq x_k < 0\\
-3 & \text{if} & x_k < 2 \cdot ARC_{k-1}
\end{cases}\\
e_k &= x_k - ARC_{k-1} \cdot a_k\\
ARC_k &= ARC_{k-1} + \mu\cdot e_k \cdot a_k\\
h_i(k+1) &= h_i(k) + \mu\cdot \operatorname{sgn}(e_k) \cdot a_{k-i}
\end{align}
where

$N$: DFE length
$x_k$: DFE output
$ARC_k$: Signal level estimate
$a_k$: Symbol estimate for a 4-level PA
$e_k$: estimate error
$h(k)$: DFE weights.

I tried $\mu=0,125 * 2/N$, but it seems to be wrong in some cases.


texts: Modal analysis and frequency response function - interpretation needed

texts:
For couple of weeks now, I have been trying to measure cutting forces during machining processes. The system which I am currently examining is highly distorted at high frequencies due to its structural dynamics. My set-up consists of force sensor, workpiece holder and workpiece itself - they are screwed all together with M12 screw.

The problem is that I struggle to interpret FRF of the setup correctly. For better understanding, please see attached figure.

What is the physical explaination of the peaks marked in the figure? Does it mean that forces measured at these frequencies are damped?

Thanks in advance for any kind of help.

Best,
Daniel


texts:
Image registration algorithms can be classified into intensity-based and feature-based algorithms, I want to know When and Why should i choose intensity or feature based algorithm?

another question Is there a difference between image registration and optical flow?


texts: How to convert a spatial frequency in a 2D-DFT into the units radians per pixel?

texts:
Let's say I have a 2D image, and I take the discrete Fourier Transform (via FFT) of that image.  In the frequency domain, I get the following image:

In this image, let's just assume all the spatial frequencies are zero, except for the selected point at (608, 513).  Now, given the x and y-coordinate of the spatial frequency (note this has been centered via fftshift), how can I convert the spatial frequency for this selected pixel into radians per pixel?

My Guess: Since the 2D-DFT is centered, our X spatial frequency is 608-513 = 95, and our Y spatial frequency is 513-513 = 0.  Then, we need to multiply by 95 by 2 * pi and divide by N to get radians per pixel?


texts: What transformations have been applied to this image to yield four derivatives?

texts:
In a few days I have an exam of image processing. I'm currently having a look at the example exam and I came across this question. "Which transformation of this image is being done here. Explain the components in the quadrants NA (second image)."

Na

Can anyone help me with this ?


texts: Purely theoretical question about idelal filters and infinite oscillations

texts:
I am asking perhaps a naive question, but still, it would be nice to have this formally stated one time:

in theory, if we could do it (of course, we cannot, but imagine we could), if we could convolve an upsampled image with an infinite sinc function (resulting from the ideal low-pass filter in freqency domain, i.e., the brick-wall) would the interpolation really be perfect, i.e., no Gibbs ripples ?


texts: Find out when an audio file is actually recorded?

texts:
Is there any technology to help find out when an audio file is actually recorded?

Here is the problem:

Given an audio file, how can we find out when it is recorded? File creation time is useless because it can be altered.


texts: minimum-phase phase via Hilbert transform returned values

texts:
Following my previous question: HRIR Minimum phase I managed to compute the minimum-phase phase of a FIR filter (in my particular case, HRTF filters).
However I am not sure of the phase values returned by my function.
The Python function is the following. It expects an arbitrary HRIR as input. I am using a 44100 long FFT to have the freqency bins perfectly aligned at integer frequencies at 44100Hz sampling rate (so for example a sinusoid at 1000Hz will have a magnitude of exactly 1 with no spectral leakage)
import numpy as np
import matplotlib.pyplot as plt
from scipy.signal import hilbert
from scipy.fft import fft, ifft, fftfreq

def min_phase_coversion(HRIR):

    HRIR_fft = fft(HRIR,44100)

    xf = fftfreq(len(HRIR_fft), 1/44100)

    phase=np.angle(HRIR_fft)

    minimum_phase = np.imag(-hilbert(np.log(np.abs(HRIR_fft))))

    plt.figure()
    plt.plot(xf, phase)
    plt.plot(xf,minimum_phase)
    plt.plot()
    plt.grid()
    plt.show()

For the following plots I am using the KEMAR MIT dataset, in particular the plot have been computed for the HRIR at [azimuth=144°, elevation=30°, distance=1.4m]. The selected HRIR is the following (one per ear):

Once I perform my minimum-phase computation using the above Python function, this is what I get (phase for a single ear):

However, looking at the plot, I can't understand the values returned by
np.imag(-hilbert(np.log(np.abs(HRIR_fft))))

Is the returned phase unwrapped? Or am i doing something wrong? The original phase "wraps" as expected once it gets over +3.14 or under -3.14 (+/- pi), while the minimum one doesn't.


texts:
kernel A = [ 0 -1 0 -1 4 -1   0 -1 0 ]
and B = [-1 0 -1  0 4 0 -1 0 -1 ]

i have found the frequency spectrums with freqz2 function of Matlab.

A's spectrum has the highest magnitudes in the corners (1,1),(-1,-1),(1,-1),(-1,1)
and B's in (1,0),(0,1),(-1,0),(0,-1).

I know that the center of the frequency spectrum has the low frequency information and the corners the high frequency information such as edges. It appears that B filter produces a sharper image although A has higher magnitude in the corners. Why is this happening ?


texts: Can anyone help me with good reference books for Discrete Wavelet Transform (DWT)

texts: I am working on EEG signal. To analyze the signal I need to use DWT. Therefore I need good reference books. If the book contains MATLAB implementation of DWT that will be more helpful.


texts: Error analysis for uniform quantization of uniform input

texts:
I was told that

  error energy $$\delta^2_e=\Delta^2/12,\quad \Delta=(B-A)/2^R$$
  signal energy $$\delta^2_x = (B-A)^2/12$$
  Here signal bounded:$A\le x[n]\le B.$ $[A,B]$ is split into $2^R$ subintervals.

But i don't understand how the signal energy was derived. More precisely, i don't understand why signal energy equals to variance.

EDIT: Ah, probably because we center signal around $\frac{B-A}{2}$.


texts: How to properly smooth the phase of a spectrum (or any unit-complex function)

texts:
I want to smooth the phase of a measured (transfer) spectrum without destroying unit-complexity of the phase factor. Suppose
$$f:\mathbb{R}\to \mathbb{C}\qquad , \qquad f(\omega)=r(\omega)\cdot {\rm e}^{i\phi(\omega)}$$
is the spectrum, then it is well-known how to smooth the power spectrum (i.e. $f^*f=r^2$): just multiply the autocorrelation function (i.e. inverse Fourier transform of $f^*f$) with a window function $w(t)$ with limited width and use the convolution theorem to show that this is the same as filtering $f^*f$ with the Fourier transform of the window. For example, if $w(t)$ is a Gaussian, $w(\omega)$ will also be Gaussian (with reciprocal width with respect to the width of $w(t)$).

But what about smoothing the phase? If one simply smoothes the phase factor (which has unit modulus) in the same way, the result usually will not have unit modulus anymore. Moreover, phases that change rapidly but linearly will cause vast cancellation. But fast phases in a transfer function correspond to big translations/lags of wave groups, which I want to be represented faithfully in the transfer function. If I just cut them off (which is what happens if I filter the phase factor), bigger translations/lags will just disappear.

For example take
$$f(\omega)=e^{i\omega t_0}$$
as a pure phase function, then averaging this over a frequency interval from $\omega=\omega_0-\pi/t_0\dots \omega_0+\pi/t_0$, the result will be zero (in complex) independent of $\omega_0$ because it always integrates over a whole unit circle:
$$\bar{f}(\omega_0)=\frac{1}{2\pi/t_0}\int_{\omega_0-\pi/t_0}^{\omega_0+\pi/t_0} e^{i\omega t_0}d\omega = 0$$
By contrast, the average phase angle will be linearly increasing, i.e.:
$$\bar{\phi}(\omega_0)=\frac{1}{2\pi/t_0}\int_{\omega_0-\pi/t_0}^{\omega_0+\pi/t_0} \omega t_0 d\omega = \omega_0 t_0 = \phi(\omega_0)$$

So what I actually want to smooth is the phase angle $\phi(\omega)$ itself! But joining phases continuously is very challenging for measurement data, so the ATAN2 function is not an option. What I have so far come up with is the following. Suppose we had a pure phase function
$$f(\omega)={\rm e}^{i\phi(\omega)}$$
Then
$$f^\prime(\omega) f^*(\omega)=i\phi^\prime(\omega)$$
So if $\phi(\omega)$ is a linearly increasing phase angle, $\phi^\prime(\omega)$ is constant, independent of how fast the angle changes. Hence, $\phi^\prime(\omega)$ is perfectly destined to be smoothed by convolution, in my opinion. I don't know yet how to properly reconstruct a (smoothed) phase factor from that, other than directly integrating it.

Does this make sense? I want make sure it does before I calculate any further. I haven't seen anything like that in the literature. Maybe this method is already known. Of course, what I would most enjoy in the end is a simple smoothing formula that does not depend on whether $f$ has unit modulus or not.

Edit: Upon further thinking, I have recognized that the computation of the phase gradient (aka group delay) via $f^\prime(\omega) f^*(\omega)=i\phi^\prime(\omega)$ is exact only for a continuous spectrum, while it can get very inaccurate when applied to a discrete spectrum. The derivative for the discrete spectrum becomes a finite difference and if the phase changes very rapidly, the angles between adjacent spectral points become very blunt, making the finite differences inaccurate. So there is no way around computing the phase explicitely and unwrapping it, like Dan Boschen has proposed in his solution.


texts: How does reconstruction function relate to watershed segmentation?

texts:
I'm trying to understand how watershed work by reading The Morphological Approach to Segmentation: The Watershed Transformation S. Beucher and F. Meyer

This document mentions the reconstruction function. But I don't understand how reconstruction function relates to the watershed algorithm?

More detail, in reconstruction algorithm implementation, it is implemented by ordered queue (below image)

And input of algorithm :

First, I understand parcel and film in function but I have no idea about it in image. Can you give me some example about it on image and how can we get film known parcel ?

Finally, result of this algorithm :

I know that is not result of watershed, so how does it relate to watershed's output ?


texts: Does decoding jpeg images uncompresses them?

texts:
I train a convolutional neural network on compressed jpeg images. I noticed that my code takes pretty much the same amount of space in memory whether I use high quality jpeg images or the same jpeg images but compressed. I know that processing a jpeg image starts by decoding it, so my question is: does 'decoding' a jpeg image uncompresses it?
In other words, does a compressed image regain its original size (in bytes) when decoded?


texts: Imposing sparse structure in Gabor space

texts:
I have a 2-D signal collected in the frequency vs time domain. There is a feature of interest in this signal (a chirp) appearing at some arbitrary times which I have successfully extracted using a properly orientated 2-D Gabor filter.

What I would like to do is exploit the detection of this feature to represent the data in a sparse domain, so that Compressed Sensing can be used to capture the data regions which contain the chirp trait without sampling the entire data set a-la Nyquist.

I have searched for a way to do this online but haven't found very much. The literature I've read says that the basis must be orthogonal for perfect reconstruction, and that Gabor domain is most times not, but given that I'm only interested in a faster detection (yes contains/no, does not contain chirp) and not really a perfect reconstruction, is such a thing possible? If the approach above is not the best way to go about doing this, recommendations would be welcome as well.


texts: How does Fourier Transform manipulates image technically?

texts:
I need a little bit more explanation on technical details regarding Fourier Transform in Image Processing, mainly on the frequency aspect.

I see that we can use Discrete Fourier Transform to change an image into a signal. From spatial domain to frequency. But from this definition, suppose I want to manipulate, let's say, for example, overall intensity of an image, once my image is in frequency domain, how can I know which frequency is what in the spatial domain? How could I know which frequency belongs to what part of an image?

If there is anything seems incorrect, please correct me.


texts: How can I generate an ideal 2D low-pass filter in MATLAB?

texts:
I've been tasked with creating a 32 x 32 half-band low-pass image filter in MATLAB. My thinking is to generate the ideal filter mask in the frequency domain and compute the corresponding convolution mask using the inverse FFT. I first generate the filter in the frequency domain.
filter = zeros(32);
filter (1:8, 1:8) = 1;
filter (1:8, 24:32) = 1;
filter (24:32, 1:8) = 1;
filter (24:32, 24:32) = 1;

This filter turns out to be the following in the frequency domain. I've confirmed my MATLAB code produces this pattern which is symmetric.

Note, I'm assuming I need to define the mask with the frequency ranges from 0 -> 2pi rad/sec, hence putting the "ones" in the corners.
I then generate the convolution mask using the "iift2" MATLAB function. However, this mask turns out to be complex which raises some red flags. I would expect it to be a real filter.
mask = ifft2(filter)

When I convolve this filter mask with my image however, I get some weird artifacts.
image_filtered = imfilter(image_orig, mask);

I'm curious if there is some flaw in my thinking here.


texts: Algorithm to Count Zeros Outside Unit Circle for FIR Filter

texts:
As detailed in this post Can I set a constraint on the first tap of an FIR filter such that its inverse is stable? I show how Cauchy's Argument Principle can be used to easily confirm if an FIR filter is minimum phase, meaning all zeros are inside or on the unit circle, or the converse, maximum phase, simply by plotting the frequency response on a complex plane.

However, the approach I use is completely graphical, in that I have to plot the frequency response and then manually count encirclements of the origin. Each encirclement of the origin indicates a zero that is outside the unit circle.

My question what would be  a compact and efficient algorithmic approach (not graphical) based on the Cauchy Argument Principal that would be more efficient than solving for the roots and determining that $|z|>1$ for all roots, or any other approach to count the number of zeros outside the unit circle?

The graphical approach is as follows, using the example from the linked post above with a maximum-phase FIR filter with coefficients [1 -3 -3 2 5].   The frequency response (magnitude and phase versus frequency) is plotted below on a complex plane. The graphical way to count encirclements is to draw a vector from the origin out toward infinity at any angle and count how many crossings of the frequency response take place. Below you would count 4 encirlements indicating all four zeros are outside the unit circle. (A minimum phase filter would have no encirclements).  The easy way to count encirclements is to note a direction on the frequency response with a forward direction consistent with increasing ω, and then draw a vector from the origin out toward infinity at any angle and count how many crossings of the frequency response take place: if the cross is of a forward direction the count increases, and if of a negative direction the count decreases. Note given the crossing locations it is not as simple as just counting crossings on the Real or Imaginary axis (here the Real can't work but the Imaginary could, but that is not a general solution).


texts: How is DJI OcuSync wireless link physical layer implemented?

texts:
I'm looking for documentation on the wireless command link of the DJI OcuSync. I found this website (in Chinese). It seems like an OFDM link that uses Alamouti codes?

I'm specifically interested in the receiver signal processing involved in synchronization of this link. Can anyone point me to useful papers or block diagrams?


texts: Why does the impulse response determine the transfer function of a system?

texts:
why does it describe the transfer function.. How come? Especially for LTI systems.
I thinking about the theory about how come an impulse input can provide information about a complete system.. As it doesn't operate fully of all configuration.


texts: Calculating MTF: negative values in the LSF

texts:
I want to implement an MTF measurement. The goal is to measure a tilted slit, like this one:

However the first half of the LSF is positive but on the second half (white-black edge) it obviously has negative values:

Applying a FFT, I get a curve that almost always looks the same and does not really depend on sharpness alterations. In my understanding this curve is meaningless.

However, when I set all negative LSF values to zero, I get a curve that really looks like an MTF curve:

Is "cutting off" negative values a viable solution when measuring a slit?

I hope this is the right place to ask.

Regards


texts: Aliasing in the Short time Fourier Transform of a pulse

texts:
When attempting to take the Short Time Fourier Transform of a pulse, at the end of the pulse I'm running into problems.

The signal looks like this at the end (it is a simple $sin^2$ pulse envelope, with a carrier frequency $\omega$.  The only frequency that I care about in the STFT is $\omega$):

And when I perform a STFT on it using a flattop window of $6 * 2\pi/\omega$ width (and only look at the carrier frequency):

The error here for the STFT seems abnormally large towards the end.  When I look at the difference between the envelope ($A * sin^2(\frac{ t * \pi }{ T} )$ where A is fitted by hand to minimize the error  in the center of the pulse) and the absolute value of the STFT produced signal, I get something that looks like this:

The two lobes at the end are what the problem is.  Does anyone have any sources that discuss this source of error?  Properties, how to mitigate them, etc?


texts: Rake Receiver vs OFDM and multi-path channel

texts:
My Question is, Can we use the Rake Receiver with OFDM modulation? As you know, Rake receiver is usually used with DSSS modulation, so can we use it with OFDM?

in other words, Does the Rake Receiver depends on type of modulation or on the channel environment (I mean, If the channel has high multi-path, the Rake Receiver will be better to be used regardless the type of modulation)

Thank you


texts: Convolution and Modulation Theorem with Octave Example

texts:
As far as I know, convolution theorem stated that
$$
x[n] * y[n] = X(e^{j \omega}) Y(e^{j \omega})
$$
$$
X(e^{j \omega}) * Y(e^{j \omega}) = x[n] y[n]
$$

I tried a simple experiment using GNU Octave which code is provided in the following

x = zeros(1,31);
h = zeros(1,31);
x(14:18) = [20 -30 40 -30 20];
h(14:18) = [1/5 1/5 1/5 1/5 1/5];

So now, we have

x =

 Columns 1 through 20:

    0    0    0    0    0    0    0    0    0    0    0    0    0   20  -30   40  -30   20    0    0

 Columns 21 through 31:

    0    0    0    0    0    0    0    0    0    0    0

h =

 Columns 1 through 10:

   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000

 Columns 11 through 20:

   0.00000   0.00000   0.00000   0.20000   0.20000   0.20000   0.20000   0.20000   0.00000   0.00000

 Columns 21 through 30:

   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000   0.00000

 Column 31:

   0.00000

Next, the following code demonstrates whether $ x[n] * y[n] = X(e^{j \omega}) Y(e^{j \omega}) $

conv_result = round( conv(x,h,shape="same") );
freq_product_result = round( ifft( fft(x) .* fft(h) ) );

Then the results are

conv_result =

 Columns 1 through 25:

   0   0   0   0   0   0   0   0   0   0   0   4  -2   6   0   4   0   6  -2   4   0   0   0   0   0

 Columns 26 through 31:

   0   0   0   0   0   0

freq_product_result =

 Columns 1 through 25:

   0   6  -2   4   0  -0   0  -0   0  -0   0  -0  -0   0  -0  -0   0   0  -0   0  -0   0  -0   0  -0

 Columns 26 through 31:

  -0   4  -2   6   0   4

Indeed both have the $\{4 , -2  , 6 ,  0  , 4 ,  0 ,  6, -2,   4\}$ pattern, but the pattern doesn't appear at the same column index. Is my experiment correct? Should both the result be identical? I think I still don't understand the concept completely.


texts:
I am building an optical OOK (On off Key Ring) system without a carrier frequency. [However I have a guard time in between symbols, so a consecutive "1" messages will result in a pulse train as opposed to DC, see image].  Essentially, the presence of signal indicates one and lack of it indicates zero. I have a precise clock that is syncs the receiver with transmitter. The system operates with a low SNR and I like to improve SNR using DSP techniques.

I have a few questions:

I do selective sampling in my hardware, in other words, I do not continuously sample the channel but only sample when the probability of seeing the signal is maximum (i.e. this is a light pulse, I time the ADC such that ADC samples at the end of the pulse where I know the whole analog chain is stabilized). See image.

Naturally, this drawing doesn't show the noise but it is there. This is particularly a  low signal system and primary noise sources are shot noise, johnson noise and amplifiers internal noise. (optical system so there are no other interferer except Sun). My observation of noise indicates, it is similar in all frequencies. (At least what I see in Scope)

I now use a simple threshold comparison in software to determine if the data is one or zero.  Is there a better way? I have thought of some options but I like to hear from the experts.

So far I have considered the following options:

Do continuous ADC and try to integrate during rise time:Not entirely
sure on the benefit (There may be other benefits, I don't know).
Matched filter in software: Don't really understand the math but
based on what I read, a possibility
Sample during the guard time and subtract this from the signal ADC
value (This may    provide some further details but also not so sure,
guard time would be the noise measurement)
Change hardware to a synchronous decoder, costly, time consuming and
may not work well    since my data rate is fast and getting a
synchronous demodulator would mean expensive board since I have to
built a multi-MHz carrier frequency system.


texts: Is it possible to define Fourier or wavelet transforms on DNA sequences?

texts:
I am wondering how and if it is possible to define a Fourier transform or Wavelet transform on DNA sequences which are basically arrays with the values $\{T,C,G,A\}$ in them.
I have found a paper which uses a "Voss 4D binary indicator representation":

Yin, C., Chen, Y. and Yau, S.S.T., 2014. A measure of DNA sequence
similarity by Fourier Transform with applications on hierarchical
clustering. Journal of theoretical biology, 359, pp.18-28.

Their code is supposed to be available online FFTDNA4D.m, but the link is broken.
The general idea is to define an indicator mapping per nucleotide type ($\alpha \in \{T,C,G,A\}$) for the DNA sequence $s(0), s(1), \dots, s(N-1)$
$$
\begin{equation}
  u_{\alpha}(n)=\begin{cases}
    1, & \text{if $s(n)=\alpha$}\\
    0, & \text{otherwise}
  \end{cases}
\end{equation}
$$
Then the DFT is given by:
$$
U_x(k) = \sum_{n=0}^{N-1} u_x(n) e^{-2\pi i \frac{kn}{N}}
$$
the power spectrum is then given by:
$$
\Phi(k) = \sum_{x\in \{ A, T, G, C\}} |U_x(k)|^2 \quad k = 0,1, \dots, N-1
$$
Is there a better way to do it? This method masks out the differences between the different types of DNA bases.


texts:
I am doing research on images using gabor filter.
I have read the following post:

Apply a Gabor filter to an input image

There must be some formula to determine the value for gamma,sigma,lambda,psi...can i know how to determine it?


texts:
I have to find the Fourier transform of $$ x(t)= \frac{1}{T}e^{-\frac{t-T}{T}}u(t-T) $$

First I applied traslation property , so $$ F[x(t-T)] = X(f) e^{-i 2 \pi f T} $$

after I applied time scaling property , and I obtained
$$ F[x(\frac{t}{T})] = |T| X(Tf) e^{-i 2 \pi f T} $$

Now I calculate $$ \int_{T}^{+\infty } \frac{1}{T}Te^{-\frac{t-T}{T}}u(t-T) e^{-i 2 \pi f t } $$ that’s results $$ \frac{T}{1+ 2 i\pi f t } $$

So I obtained that the final Fourier transformation should be $$ X_f = |T| \frac{T}{1 + 2 i \pi f T} e^{- 2 i \pi f T } $$

On my book the result is the same except for T|T| but I don’t know what’s wrong. Thank you so much


texts: Finger print skeletonization without using matlab functions

texts:
I am trying to thin foreground pixels of a finger print until they are one pixel wide without using the built in bwmorph and thin operations.

Can anyone explain to me how to write them in matlab without actually using the built in function. Thanks.


texts: How would you describe this signal? It's like a sinusoid but as if its frequency was constantly increasing: could you write down a mathematical description? Thanks. And no, this is not homework...


texts: ESPRIT signal subspace and rotational invariance

texts:
so in ESPRIT[1,2] we have some signal of the form:

$$y[n] = \sum_{k=1}^{K} a_i(k) e^{j\omega_k n} + N[n]$$

where N is noise, then we get the autocovariance matrix

$$R_{yy}=E[y[n]y[n]^*]$$

take the eigen-decomposition

$$\hat{R}_{yy} = U \Lambda U^*$$

and this is broken up into our signal subspace/noise subspace

$$U = [S \; N]$$

ESPRIT really loses me at the next step where it selects the first M and last M columns of the signal subspace:

$$S_1 = [I_{M-1} \; 0]S $$

$$S_2 = [0 \; I_{M-1}]S$$

then states:
$$S_1 \phi = S_2 $$

not, sure what that statement means.

it then computes
$$\phi = (S_2^* S_2)^{-1} S_2^* S_1$$

I'm having a hard time wrapping my head around how selection of the first and last M columns relate to the frequencies embedded within the signal calculated as the eigenvalues of $\phi$. The rotational invariance statement $S_1 \phi = S_2 $ is pretty non-intuitive.

[1] https://ieeexplore.ieee.org/document/32276

[2] http://citeseerx.ist.psu.edu/viewdoc/download;jsessionid=67D758D42882A6927A1947F292E3A27D?doi=10.1.1.13.1219&rep=rep1&type=pdf


texts: When should the sum of all elements of a gaussian kernel be zero?

texts:
I found an approximation of a 5x5 2D convolution kernel like this :

Here, the sum of the elements is zero and this one was used for Laplacian of Gaussian!

Another one here :

This one has all positive and the sum is not zero and is being used for smoothing.

Could anyone please tell me what is with the sum being zero in the previous case ?

References:

Link1 and  Link2


texts: Meaning of the frequencies vector using fft

texts:
I'm so confused, and it would be nice if someone could clear my head, or at least direct me to suitable references.
Assume that I have two gaussian pulses with the same bandwidth in time, but for 2 different wavelengths. I'm going to handle them separatly- 2 separate graphs of one gaussian in time.
Both graphs are the same because they have the same bandwidth.
Then, I use fft on each graph.
Of course I'm going to have the same graphs after the fft- because I had the same graphs before.
Now, I want to creat the frequencies axis for each signal. In order to do that, I use the exact direction from matlab. df=1/T,  T=N·dt.
So the frequencies axis are also the same!
But how can it be if f=c/lambda?  I'm talking about two different wave length, how is it possible that the both have the same spectrum of frequencies?
Further more, I don't understand what is the different between F and df and why is F=N·df. I do understand that T is the time in which I sample, dt=T/N. but that's it.


texts: Effect of aliasing when sampling at higher rate

texts:
Considering a signal $f(t)=\sin(2\pi t)$ with frequency $10\textrm{ Hz}$:

Middle and bottom plots show the results of sampling at different rates.

I'm trying to understand the case shown at the bottom.

Is the sampling rate too high that it's capturing noise?
What's the name of this effect? May I call it a kind of aliasing caused by high sampling rate?


texts: All window functions, which go from zero amplitude, up, then back to zero amplitude, are all base-band signals. Why is this? Are all functions that have this behavior always band-band aka low-pass?


texts: Digtial FIR Impulse response & transfer function

texts:
I am currently working the figure through below.
as it is an FIR Filter i have worked out using convolution that the output is 4,2,4,6,0,0.
i am trying to obtain the 'z' domain transfer function of the filter below using h(n)
using the $h(n)$ 2,-1,3
I have started to obtain the $z$ domain transfer function with
$y(n)=2x(n)-1x(n-1)+3x(n-2)$
any help?


texts:
I'm trying to smooth a Sawtooth wave coming out of my LFO.   Effectively the transition from 1.0 to -1.0 in a 1 sample step is causing clicks and pops if I modulate my audio signals amplitude with this LFO output.

What is the best way to smooth this basic waveform to eliminate the clicks without adversely affecting performance (significantly)?


texts: How does sensor fusion work? I want to understand the math/physics/algorithms

texts:
Sensor fusion algorithms can give a more precise 3D orientation (and possibly postion?) of a device by combining readings from an accelerometer, gyroscope, and magnetometer.

Can anyone explain, or provide links, explaining the details behind this? I want to understand the math and physics, such that if I have raw 9-DOF sensor data, I might be able to implement sensor fusion. Assume I have sufficient background of linear algebra, calculus, etc.

Cheers!


texts: SNR for complex valued signal and when to use SNR and when to use Eb/N0?

texts:
This question is in continuation to the one asked earlier Help in understanding the formula of Signal-to-Noise-Ratio (SNR) - Part 1

Consider a model
$$y[n] = x[n] + w[n]$$

Signal means the desired information which in this case is x and y is the output noisy observation.

Question A

Case 1: When noise is AWGN $ \sim C\mathcal N(0,2\sigma^2_w)$; $y$ and $x$ is also complex valued i.e.,
$w \sim CN(0,2\sigma^2_w), x \sim CN(0,2\sigma^2_x)$

Then would SNR = $\frac{E[x^2]}{2\sigma^2_w} = \frac{2\sigma^2_x}{2\sigma^2_w} $ ?

Case 2: AWGN $ \sim \mathcal N(0,\sigma^2_w)$; $y$ and $x$ are real valued. Then would SNR = $\frac{E[x^2]}{\sigma^2_w}$?

These topics are not clearly mentioned in the Answers given for the other Question. The Answers explain what is SNR but my confusion remains that if the signal and the noise are complex valued, then would there be any change in the SNR due to the real and imaginary component? In the denominator of the formula would there be a 2? That is why I have asked this here.

Question B: https://en.wikipedia.org/wiki/Signal-to-noise_ratio shows 2 formulae for SNR. In many research articles, the formula of SNR is not used, instead EbN0 is used. I am using Least squares estimator to estimate channel coefficients and then trying to plot the mean square error (in db) of the channel estimates with varying SNR. NExt, I want to plot the bit-error rate (using simulation) for bpsk and say QAM modulation. I used SNR for the plot of mean square error. So, should I use SNR for Bit error rate? OR is there no such rule that both the plots should be using SNR and not EbNo formula?


texts: Problem with calculating recursive equation

texts:
I'm facing with this error when I'm trying to solve the equation :

  Attempted to access x(101); index out of bounds because
  numel(x)=100.

actually the original equation is

so I used "i+1"instead of"i" to avoid following error:

  Attempted to access y(0); index must be a positive
  integer or logical.

n=[1:100];
yy=zeros(1,99);
xx=zeros(1,99);
y=[1 yy];
x=(.5).^n;
for i=n

  y(i)=.5*y(i-1)+2*x(i)+4*x(i-1)
end

error resolved by initializing i by  2

the improved code is :

yy=zeros(1,99);
xx=zeros(1,99);
y=[1 yy];
x=(.5).^n;
 for n=2:100

 y(n)=.5*y(n-1)+2*x(n)+4*x(n-1);
end

but the results of **y(n)**is different from what I get on paper.


texts:
So I'm used to DFT being on an equation like this:

x(n) = cos(n), ∀n , with window n = 0,1,2

I can clearly plot discrete values and take a DFT to find Xk

But I found this equation on a past exam that's like this

x(n) = δ(n+2) +δ(n+1)−δ(n−1) with window of length 5.

This appears to be an impulse response to me, and I don't understand how one would take a DFT of it?


texts: Colored noise and adjusting covariance matrix with change of frequency

texts:
Assume a signal sampled at 2Hz.
Assume it is a colored noise with known PSD.
The measurement equation to be used in a Kalman filter requires a covariance matrix.
If I want to sample at 9Hz (for example), how should I modify the covariance matrix?


texts: RF direct conversion down to DC and signal recovery from IQ

texts:
I've previously worked with super-heterodyne receivers where the IF spectral content has been above DC.

Reading the below tutorial; direct down conversion receivers are used to mix a signal of interest down to base-band. Figure 7 shows the output of the direct conversion process where the modulated signal has its frequency content either side of DC i.e. half as negative frequencies.
https://www.microwavejournal.com/articles/3226-on-the-direct-conversion-receiver-a-tutorial

A similar question was asked in the link below which explains that there will be indeed overlap and corruption because of the negative frequencies.
Downconversion / Demodulation of RF signal to DC

If that is the case, how is the signal recovered and what do the resultant I and Q samples look like in terms of magnitude and frequency content?

If specifics help; if we assume an already modulated carrier 1500MHz with +/-1MHz FSK signal is mixed with a 1500MHz LO. The output (and image) will be centered around DC but half the signal bandwidth will be negative, the missing half will be on-top. Its this with which I'm confused.


texts: New to audio processing, why is my data so lopsided when I try to break into frequency ranges?

texts:
I have a wav file audio data, I broke it up into 1024-length windows (no overlap), and performed fft on each one.

If I visualize this data it actually looks pretty good, but the problem is that the data is really lopsided. It seems like buckets 1-3 get good activity, but they're generally much larger than buckets 4-8, so when I visualize the data I have to have a weird conditional multiplier on the higher frequency buckets so I see some activity.

So then, what is the proper way to break my fft into frequency buckets? A simple explanation would be best. Thank you!


texts:
I am trying to check whether this black circle is exist or not on my images.

My current approach is like that;

Hough Circle detection
Crop the outer roi
Adaptive threshold

count zero pixels
if ( zero pixels >500 )
black hole exist
else
black hole does not.

I just want to know is this method useful if this algorithm go for real time check ? if not what kind of methods you could suggest me ?

Here is a image with blackhole

Here is a image without a blackhole.


texts: Why can't music recognition services like Shazam or SoundHound recognize music if its speed is changed by more than 5 percent?

texts: Music services like Shazam and SoundHound take an audio sample from your microphone and compare it to its database to recognize music. This algorithm works very well under high noise, unequal frequency bins, and even a simple hum or whistle to the music. However, if you change the song's pitch or speed by more than 5%, the algorithm cannot detect the song at all. A 5% change is very small and not even that noticeable, yet it completely fools the music detection algorithm. Why is this? Is it due to the different speed/pitch missing the targets in the audio fingerprint?


texts:
I'm new to SP Stack Exchange . I am doing a project on implementing a handwriting recognition machine using Neural Networks in Real Time .

The Image Processing Part involves convolving a 4*5 kernel over a 28*28 image. I want to know the pros and cons involved with using a DSP vs FPGA .

The FPGA I plan to use belongs to Xilinx's Spartan 6 family which has 90 DSP slices. Each DSP Slice has 18*18 bit MAC units which can operate in parallel at 200MHz.

My college has a DSP starter kit TMS320C6713 DSK from Spectrum DIGITAL. I'm not sure how MAC operations are performed in DSP .Do we have parallel MAC units? The datasheet isn't very clear on the MAC units .

Are there any advantages on using the DSP over the FPGA which I might have overlooked?

Any help is appreciated .


texts:
I was asked a question, as posted here, and the answer given is (A) i.e. $\frac{3}{2} A^2 N_0$.
My solution steps was:

Finding the mean of the output process: Since input is gaussian the output will be also Gaussian, as per property of Gaussian variable.

Finding the mean square value (MSV) of the process i.e. $\mathbb{E}(Y^2)$.

Variance = MSV$-\mathbb{E}(Y)^2$

But I am confused do mean of output process will be zero? Since the input is zero-mean white Gaussian noise if yes how?
Is there any other way to solve the question?


texts:
My objective is to create a Gardner loop to remove the timing offset present in my signal. To achieve best timing SNR, I downsample my signal to 2 samples per symbol. It is depicted below:

where blue and red vertical lines indicate the symbol and the sample duration respectively.
Next, I plug my signal into Gardner loop, here is my code:
def Gardner_loop(time_support, sig):

    # Initialization of the parameters
    sample_rate = 1/time_support[1]
    sps = int(sample_rate/f_symb) # samples per symbol
    sig_n_1 = int(np.ceil(sps/2)) # x(n-1)
    sig_n_2 = int(sps) # x(n-2)
    loop_filter = f_symb//20 # BW of the loop filter
    loop_num = 100
    error_symbol = np.zeros(loop_num)

    error_list = []
    sig_wo_time_err = []
    step = sps

    # Gardner loop
    for i in range(1, loop_num):
        # Calculate the error
        error_symbol[i] = sig.real[sig_n_1 + step + int(error_symbol[i-1]*loop_filter)]\
        *(sig.real[sig_n_2 + step + int(error_symbol[i-1]*loop_filter)]\
        - sig.real[step+int(error_symbol[i-1]*loop_filter)])\

        + sig.imag[sig_n_1 + step + int(error_symbol[i-1]*loop_filter)]\
        *(sig.imag[sig_n_2 + step + int(error_symbol[i-1]*loop_filter)]\
        - sig.imag[step+int(error_symbol[i-1]*loop_filter)])

        # Make a step
        step+=sps

        # Save the error
        error_list.append(error_symbol[i])

        # Save the signal with reduced error
        sig_wo_time_err.append(sig.real[sig_n_2+step+int(error_symbol[i-1]*loop_filter)] \
        + 1j*sig.imag[sig_n_2+step+int(error_symbol[i-1]*loop_filter)])

    error_to_plot = np.array(error_symbol)
    # Plotting
    plt.plot(np.arange(1, loop_num+1), error_to_plot)
    plt.xlabel('loop number')
    plt.ylabel('error')
    return np.array(sig_wo_time_err)

the plot of the error is represented in the figure below:

the constellation looks like this:

indicating that there is still the timing error.
My first specific question is: can somebody check if my implementation of Gardner loop is correct? Because I am not sure if the code I wrote does what is needed and I cannot find any software implementation of Gardner loop in internet.
In addition to that, I know that I need to interpolate the signal at specific points. These points are indicated by the error from Gardner loop. My understanding of the process in the following (see first figure): if the error is positive, then the sampling is too late and the interpolation must happen between 0 and the red line (sample duration). For negative error, interpolation happens at the last sample of the symbol (between red and blue lines). Then we use this interpolated point for further processing inside the loop and eventually drive the error to 0.
To implement this, in the code mentioned above, inside the for loop I use the function scipy.interpolate.interp1d to get the value of the signal at the time instant indicated by the error. However, with an interpolation length consisting of two points (0 and red line), I receive an error and cannot get the point in between (if needed, I can post my code as well).
So the last question is how do I implement the interpolation inside the Gardner loop.
Thank you!
Update:
The eye diagram of the signal before matched filtering and upsampling is depicted below:

The eye diagram after matched filtering (MF)

After interpolation (upsampling with symbol_rate * 5) and without MF, the eye-diagram looks like this:

After upsampling symbol_rate*5 and MF:


texts: Blocking effect after reconstructuring the image. Jpeg compression

texts: The image is splitted into 8x8 blocks each of which are processed through DCT, quantization, huffman encoding/decoding, dequantization and inverse dct and used to recontruct the image. But there is a horrible blocking effect after this compression process. What would be the reason why this is happening? I simply create a new matrix for storing the data as big as the original image. And each block after processed is copied into new image. For example first 8x8 block is put from (0,0) to (7,7). The second block is put on the indices between (0,8) and (7,15)


texts: Synchronization of BPSK with variable-length symbols (ISO 14443)

texts:
Does anyone have an idea how one would perform BPSK symbol synchronization if symbol length is allowed to vary?

As dumb as that sounds, this is not a theretical question.The "Proximity" Type B NFC protocol defined in ISO/IEC 14443 has such a scheme.

Brief Background of ISO 14443 Type B

For communication from a Proximity tag/card (e.g. a wireless smart card) to a reader, BPSK modulation is used*. The symbol period is (usually) equal to 16 carrier cycles. To begin communication, a start-of-frame synchronization sequence is sent, allowing the reader to A) identify the carrier frequency, and B) determine the phase for a logic 0. Then the card sends data out in characters. Each character contains a logic 0 start bit, 8 bits of data, then a logic 1 stop bit. (Think UART.) This means that between bytes you always see bit 7 -> 1 -> 0 -> bit 0.

For historical reasons, the specification allows the card to include "Extra Guard Time" (EGT) between each character. Before proceeding to the logic 0 start bit, the card is allowed to arbitrarily extend the logic 1 stop bit up to 32 carrier cycles, in increments of 0.5 carrier cycles. (So phase changes still always occur at 0° or 180° of the carrier.) This is something that can be done dynamically; the reader does not know in advance what the card will do.

Since the bit following the EGT is always the logic 0 start bit, you can of course pick back up symbol synchronization, but it complicates the whole process. It's hard to find information online on such a convoluted scheme, and the textbooks I have do not cover this. I'm wondering if anyone else has run into a similar situation and has some ideas on how decoding might work.

*The actual details are more complicated. The card performs load modulation of the 13.56 MHz carrier frequency using an 848 KHz sub-carrier. However the 13.56 MHz carrier signal can be easily filtered out by the reader, leaving just the sub-carrier information from the card.


texts: Is there any method/algorithm to estimate the magnitude of non-stationarity in a signal?

texts:
e.g. the global Lyapunov exponent can give sense of the level of chaos in the signal. Is there any reliable numerical technique to estimate "how" non-stationary (or how predictable) a signal is?
Also, is the Hurst exponent appropriate for this purpose? The time signals I'm interested in come from natural processes like fluid flow. I want to study signals based on how non-stationary they are.


texts: Looking for Fastest 2D Convolution in Python on a CPU

texts:
Convolutions are essential components of many algorithms in neural networks, image processing, computer vision ... but these are also a bottleneck in terms of computations...  In the python ecosystem, there are different existing solutions using numpy, scipy or tensorflow, but which is the fastest?

Just to set the problem, the convolution should operate on two 2-D matrices. We will here always consider the case which is most typical in computer vision:

a first matrix $A$ is the input and is typically large ($N \times N$ where $N$ is typically larger than $2^{10}=1024$),
a second matrix $B$ is the template and is typically smaller (say $M=128$),
the result of the convolution $C = A \ast B$ is padded such that it is of the same size as $A$.

Often, you need to do that with many images on many kernels, so a method that does it on many has a bonus point.

Thanks for any hint!

PS: the goal is at some point to sum up these results in this blog post which already contains some examples...


texts: Python equivalent code for Matlab bandpass function

texts:
I have used a bandpass filter function in Matlab as follows
Fs = 128;  sampling rate
Fd = [1 4]; passband freq
[data_1_bandpassed, digital_filter] = bandpass(signal, Fd, Fs);

How can we design this bandpass filter in Python to recreate the exact same output?
When using butter & filtfilt commands in Matlab & Python, I am able to get the exact same filtered signals between the two programming languages. However, I explicilty need to recreate this output of Matlabs bandpass function in Python.
I have also attached the output of the digital filter properties and values created using Matlab bandpass function.

A similar question is already asked by another user few months ago in different domain but with no solution, so I am adding that here for reference as well!


texts: Increasing SNR and Dynamic Range using Oversampling

texts:
How much gain in dynamic range and SNR can be expected if we are to oversample a signal with fixed analog input bandwidth. For Example if I have a analog filter at the input which limits the bandwidth to 100 MHz. The noise power will also be limited to the 100 MHz bandwidth ($kTB$ formula for noise).
Now if I were to oversample the data by a factor of 2 how much gain should I expect and why? I know we get 3 dB gain due to spreading of quantization noise, will a similar gain be achieved also due to limiting the noise bandwidth at the input?


texts: Cross-Correlation : Comparing two signals or finding the location of a target

texts:
I read across places that cross-correlation of two signal (a radar sending a signal and receiving it back) one can detect the location (distance) of a target from a point and that cross-correlation can be used for comparing signals.

I solved a question to find cross-correlation of two finite length sequences where:

x(l)={1,2,1,1}
y(l)={1,1,2,1}

The cross-correlation sequence came as:

{1,4,6,6,5,2,1}

The answer was also correct as per my verification.

But I am not able to understand how using the above cross-correlation sequence value I can state it to be used for -

  a)location of a target? b)comparing signals x and y?


texts: Drawing a Spectrogram off of FFT - what are the decibels relative to?

texts:
I am trying to draw a spectrogram from FFT data I get from audio files.
The problem is that I can't figure out is what the yielded power (decibels) are relative to. I'm getting some arbitrary -100db-ish numbers, but I don't understand how I can have them represent a finite intensity-scale.

If I have a gradient which will represent the intensity of these numbers, where silence is black and the maximum is yellow I need to know what the maximum is. I'm really not good in the DSP-area so help would be appreciated.


texts: Gradient of transfer function (z-transform) with respect to coefficients/parameters?

texts:
my problem is perhaps very simple but I just can not find the answer, even though this is used (though not explained) in several books:
What is $\frac{\partial G (z, \theta)}{\partial \theta}$ when $G$ is some transfer function, and $\theta$ the coefficients? (at least I call them coefficients). I tried using "regular" calculus on the following, but the result is very strange.
Example:
$$G(z, \theta) = \frac{b_1 z^{-1} + b_2 z^{-2}}{1 + a_1 z^{-1} + a_2 z^{-2} + a_3 z^{-3}} = \frac{B}{A}$$
$$ \theta = [b_1, b_2, a_1, a_2, a_3] $$
Using the chain rule, and treating this as a regular rational function in $z$, I get:
$$ \frac{\partial G (z, \theta)}{\partial \theta} = [\frac{\partial G (z, \theta)}{\partial b_1}, ..., \frac{\partial G (z, \theta)}{\partial a_3}] \\[19pt]
= [\frac{z^{-1}}{A}, \frac{z^{-2}}{A}, (-1)z^{-1}\frac{ B}{A^2}, (-1)z^{-2}\frac{ B}{A^2 }, (-1)z^{-3}\frac{ B}{A^2}] \\[19pt] $$
I agree, this looks weird.
EDIT: What I am looking for is the gradient, $\nabla_\theta G$.
Perheps there are some rules of differentiation in this domain that I do not know.
The problem is related to parameter estimation.


texts: Suggested Preprocessing methods for OCR on Circular Images

texts:
Hello this is my sample image I am going to do real time character detection on images like that.

I've tried SURF, SIFT, MSER and template matching on original image without any preprocessing. I can detect characters but when the image changes this feature detection methods does not work because there isn't a lot of interested points in the characters on the image. I am not exactly sure what kind of preprocessing I should do for higher accuracy.

so my question is what methods can you guys suggest me to do for character detection with higher accuracy.

thanks in advance.


texts: Playing/Recording acoustic signal in matlab to tell difference between 2 sounds

texts: I will be using one central node/latop to play a tone and will have a second node/laptop to record and read the sound and be able to say this node is on the left or right. So trying to figure out the best way to distinguish between the two tones. I want to have a tone set for left and one set for the right. I want the sound to have 400 Samples at 44100 for 10ms between 18-22freq just outside the human hearing range so that the sound is not annoying etc. I want the sound to be different slightly but be able to tell the diff between them for the system to say the node that recorded and processed the sound is on the left of right. How would I go about doing this in matlab because I don't dont really code in this language much and only time I did was for image processing just doing this so I can graduate.


texts: Applying DFT twice does not actually reverse an array. Instead, the first element stays in place while the rest of the array is reversed. Why?

texts:
I've heard a million times that applying DFT twice will result in a reversed array, but that is not what actually happens. Instead, the first element remains where it was, and the rest of the elements are reversed. Is there an intuitive reason for that?
Octave code:
N = 16;
x = linspace(1, N, N);

fft(fft(x)) / N

Result:
1   16   15   14   13   12   11   10    9    8    7    6    5    4    3    2


texts:
I am having difficulty with the case when a is real.

But after this, I am not going anywhere.


texts:
I have information related to the wired network topology in a room with the dimensions 50x10x20m. This includes the virtual link information containing source, target and datarate on a particular link.
Based on this information, I calculated the overall bandwidth as the sum of all the datarates on all the links in the network. Let's say this value is $1000 MHz$.

As a part converting the above wired topology to wireless using Zigbee, I intend to calculate the bandwidth when using Zigbee in the network. For this, I have removed all the multiple point to point links from the wired setup and ensured that using just one link from each device in the Zigbee network (since wireless is inherently multicast).

For the bandwidth calculation of the Zigbee network, I have used all the datarates of the wireless links and multiply them with the coding factor (which is 8 in 802.15.4's case) and sum all the datarates after multiplication with the coding factor. Doing this gives me $500 Mhz$. (Note that this is without implementing any scheduling which would be done later. It is just to get a raw idea on how much reduction in bandwidth I may expect)

I would like to know if the above method is correct when calculating the bandwidth for the wireless case.

If this is right, in case I intend to use 802.11g instead of Zigbee, how should I go about in calculating the coding factor for the 802.11g case?


texts: Resampling of signal with non uniform sampling frequency

texts:
I have a non uniform sampling frequency signal and I have to convert it in a constant sampling frequency.
I tried to interpolate it with an Hermite spline interpolation but it make a lot of wrong peaks, like in the figure:

For example at 14887433 there is a peak too big.

A cubic (4 point) interpolation is too sharp for me.

Which method can I use to resampling this signal?
Is there a way to filter a not constant sampled signal in order to use a sinc or a similar window?


texts: Design a filter which passes all frequencies except $\omega=\pm\frac{\pi}{2}$ and plot its pole-zero diagram

texts:
Also draw its normalized frequency response.

What is the ROC?

This has to be done in z-plane so there must be two poles at $+i$ and $-i$ since they cannot be included in region of convergence. Is my assumption correct?

And what is meant by normalized frequency?


texts: Hardware selection for acquiring signal in digital form

texts: I want to it clear that I am a mechanical engineer and so my knowledge on electronics and signal processing is limited. I am working in a hydrophone system. My job is to find the position of an underwater pinger with an AUV. The pinger sends the sound at a frequency of 45kHz and my hydrophone send signals in range of 0-5mV. I have doubts regarding what implementation of the anti-alisiang filter, signal conditioning, ADC conversion and digital filters for removing noise in terms of hardware. I do know that the sampling frequency should be minimal double that of the signal but that us not enough and aiming at 1MHz I am struggling to find a microcontroller board with such specs. Also, I need simultaneously sampling of two hydrophones depending on the difference in the time the signal takes to reach the two hydrophones, I would calculate the angle at which the pinger is present. Basically, I need help figuring out the hardware required to get such results.


texts: pre-processing to improve feature detector before tracking

texts:
I am trying to make tracking for soccer player, I need to detect features from this player and then estimate the distance difference of these pixels over number of frames.

first I have to detect the main features of the player and then match them to those in the next frames, then calculate the distance difference.

My problem, there is no good feature detector at MATLAB to my case, these sample from the existing detectors:

FAST

Harris

SURF

Q: I see that I can improve feature detector by doing preprocessing such like:
adaptive thresholding, opening & closing and sharping to highlight the player edges.
I have do the first one to make adaptive threshloding, and then make the feature detector for the output logical image, such like the following

After adaptive thresholding Harris, for logical image

I tried to make the feature detector on gray image after thresholding but I couldn't change logical image to gray image.

Also any suggestion or alternative approach to do?


texts: Getting the amplitude and phase of harmonics in C

texts:
I have a periodic series of 360 values. How can I use the fftw c libraries to get the amplitude and phase of the first, second, third and fourth harmonic?

If I do this for N=360

fftw_plan_r2r_1d(N,input_array,output_array,FFTW_R2HC,FFTW_ESTIMATE);

I can get the 360-point Discrete Fourier Transform (DFT) where output_array is complex Hermitian in half-complex format.

Where k$^{th}$ component of output_array is the real part and (N-k)$^{th}$ component of output_array is the imaginary part.

So if I want the amplitude or magnitude(??) of the first harmonic should I do this?

ampl_1sth=sqrt(output_array[1]*output_array[1]+output_array[N-1]*output_array[N-1])

Is that correct? then the amplitude of the second harmonic will be the same but with output_array[2] and output_array[N-2] and so on right?

And then how can I get the phase of the first, second, third and fourth harmonic with the info from the DFT?

For a reference to fftw see here www.fftw.org/fftw3.pdf

Thank you.


texts: Finding crosses and their orientation on a given image

texts:
I want to recognize table structure from some arbitrary photo image and store it in some formal notation (Let it be HTML table notation).

For example I have this blurry image as input:

After binarization I achieved this:

Now to detect structure I want to detect intersections between lines, their type (T-type, X-type, or simple corner), their orientations and their location on the image. After this I'm going to use all gained information for further joining of adjacent crosses in some structure, and translate this structure into formal representation.

The problem is in detecting crosses and their orientation:

In general these crosses may be scaled and/or rotated.
Maybe someone could help with method of solving this problem? Or maybe recommend different methods for this whole task?

That's what I have written so far:

 # -*- coding: utf-8 -*-

 import cv2
 import cv
 import numpy

 original = cv2.imread("/home/user/my_photo-1.jpg")
 grayscale = cv2.cvtColor(original, cv2.COLOR_BGR2GRAY)
 smoothed = cv2.GaussianBlur(grayscale, (5,5), 0)
 cv2.imshow("original", original)
 cv2.imshow("grayscale", grayscale)

 binarized = cv2.adaptiveThreshold(grayscale, 255, cv2.ADAPTIVE_THRESH_MEAN_C, cv2.THRESH_BINARY_INV, 7, 8)
 binarized = cv2.Canny(grayscale, 50, 200)
 cv2.imshow("binarized", binarized)
 cv2.waitKey(0)

Thanks in advance for any response/idea.


texts: What should be the cut off frequency of a Filter operating on an OFDM signal?

texts: I have a time domain OFDM signal and to this I need to apply different linear filter however I am not sure on the cut off frequency. I am simulating for multiple symbols where each symbol has 1200 central subcarriers in a 2048 FFT/IFFT bin. The samp freq is 30.72MHz. While designing filters in MATLAB using commands like "fir1", I am not sure how to set cut off freq.


texts: How to properly filter a BPSK demodulated signal?

texts:
I am new to this field and yesterday I began a project for BPSK modulation. So this is what i did:

$\ x(t) = 0 $ or$\ 1 $ binary data

$\ Fs = 15 $ sampling frequency kHz

$\ Fc = 10 $ carrier frequency kHz

BPSK modulation

$\ BPSK(t) = A.sin(2.\pi.Fc.t + \pi.x(t) ) $

BPSK demodulation

$\ S(t) = BPSK(t) . (Ac.sin(2.\pi.Fc.t))$

low pass filter

void lowPassFrequency(float input[], float fs, float output[], int points)
{
  float  CUTOFF = 2 * Fc;
  float RC = 1.0/(CUTOFF*2*M_PI);
  float dt = 1.0/fs;
  float alpha = dt/(RC+dt);
  output[0] = input[0];

  for(int i = 1; i < points; ++i)
  {
      output[i] = output[i-1] + (alpha*(input[i] - output[i-1]));
  }
}

it seems to me that the filtering part is not as it should be ... it supposed to be smoothed so i can get values like -1 and 1 ... Any ideas please?


texts:
I am going through digital filter design. But I really get confused with so many

Digital filters :

IIR
FIR
Moving Average
Linear phase filter
Allpass filter
Comb filter

Analog filter:

Buterworth
Chebychev
Elliptic

Is their any specific application where each filter is used ?

I suppose biomedical signal sensing or temperature sensor value filtering use digital filter (?).

Please suggest.


texts: Is DFE equalizer always blind or can we use reference signal?

texts:
I have made a blind DFE Equalization, I mean decision directed .
But I am not sure how to make an implementation based on reference signal known by receiver and emitter. In the back forward filter adaptation how do I take this into consideration?

My second equation is to know if the updating rule for the RLS is the rigth one ?

Here is my rls tracking DFE matlab code

`

    clearvars;
    close all; clc;
    MQAM = 2;
    N = 1e2;
    L = 3;
    LFF = 8;
    LFB = 8;
    target = 1e-6;
    error_min = 30;
    trans_min = 1e6;
    trans_max = error_min/target;
    BER=0;
    SNR=0;
    iSNR=1;
    endSNR = 16;
    stepSNR = 2;
    h = [1 0.5*exp(1i*pi/6)  0.1*exp(-1i*pi/8)]';
    eq = comm.DecisionFeedbackEqualizer('Algorithm', 'RLS', ...
               "AdaptAfterTraining",false, 'NumFeedbackTaps', LFB, ...
               'NumForwardTaps', LFF); eq.ReferenceTap = 1;
    BER_mtlb = 0;
    handle = waitbar(0, 'Please wait...');

    p=LFF; lamda=0.98; sigma=1; P1=(sigma^-1)*eye(LFF);

    P2=(sigma^-1)*eye(LFB);

    while SNR(iSNR)<= endSNR
         error_num=0;trans_num=0;  error_mtlb = 0;
         while  (error_num <= error_min) && (trans_num <= trans_min)
             MemFF = zeros(LFF,1); MemBW = MemFF; hFF = zeros(LFF,1);
             hFB=hFF; hFF(ceil((LFF-1)/2)) = 1;

             s = randi([0 1],N,1);
             smod = qammod(s,MQAM,'gray','InputType','bit');
             xin = awgn(filter(h, 1, smod), SNR(iSNR), 'measured');

             out = zeros(N,1); dec = zeros(size(s)); sdec = zeros(size(s));
             smodCplx = complex(smod);
             [y, err, weigths] = eq(xin, smodCplx);
             rxx = qamdemod(y, MQAM, 'gray', 'OutputType', 'bit');
             error_mtlb = error_mtlb + sum( rxx ~= s);
             %- RLS DFE Equalization loop block
             for n = 1: length(xin)
                 MemFF = [xin(n); MemFF(1:LFF-1)];
                 out(n) = transpose(hFF)*MemFF + transpose(hFB)*MemBW;
                 sdec(n) = qamdemod(out(n),MQAM,'gray','OutputType','bit');
                 dec(n) = qammod(sdec(n),MQAM,'gray','InputType','bit');

                 alpha(n) = dec(n) - out(n);% error of RLS

                 g1(:,n)=P1*MemFF*((lamda + MemFF'*P1*MemFF).^-1); % Gain
                 P1=(lamda^-1)*P1 - g1(:,n)*MemFF'*(lamda^-1)*P1;

                 g2(:,n)=P2*MemBW*((lamda + MemBW'*P2*MemBW).^-1); % Gain
                 P2=(lamda^-1)*P2 - g2(:,n)*MemBW'*(lamda^-1)*P2;

                 hFF = hFF  + alpha(n)*g1(:,n);
                 hFB = hFB + alpha(n)*g2(:,n);
                 MemBW = [dec(n); MemBW(1:LFB-1)];

             end

             M = min(length(s), length(sdec));
             Ryx = xcorr(1-2*sdec(1:M),1-2*s(1:M),'unbiased');
             [valMax,Delay] = max(abs(Ryx(M:M+99)));
             Delay = Delay-1;

             if (Ryx(M+Delay) < 0)
                 sdec = qamdemod(-out,MQAM,'gray','OutputType','bit');
             end
             sdec_synchronized = sdec(1+Delay:end);
             s_synchronized = s(1:length(sdec_synchronized));

             error_num = sum(s_synchronized~= sdec_synchronized) + error_num;
             trans_num = trans_num  + length(sdec_synchronized)  ;
         end % end this number of transmitted symbols

         BER(iSNR) = error_num/trans_num;
         BER_mtlb(iSNR) = error_mtlb/trans_num;
         waitbar(log10(BER(iSNR))/log10(target),handle);
         iSNR = iSNR+1;
        SNR(iSNR) = SNR(iSNR-1) + stepSNR;
 end
    SNR=SNR (1:length(BER)); berTheory =
    berawgn(SNR/log2(MQAM),'psk',2,'nondiff');  SNR =
    SNR(1:min(length(BER), length(BER_mtlb))); h = figure;
    semilogy(SNR,BER, 'b--'); hold on; semilogy(SNR, BER_mtlb, 'g--'); hold on;
    semilogy(SNR, berTheory, 'kx-'); legend('Simulation DFE', 'Theory');
    xh=xlabel('E_b/N_0 (dB)'); %yh = ylabel('Probability of bit error');
    xlabel('SNR (dB)'); ylabel('BER'); grid on;
    XMIN= 0; XMAX = 30; YMIN= 1e-7; YMAX= 1; axis([XMIN XMAX YMIN YMAX])


texts: Get autocorrelation function from the power-spectrum (python)

texts:
I am trying to compute the autocorrelation function of a signal for which I only know the power-spectrum.

In order to test my approach I wanted to try it out on the spectrum of $1/f^2$ noise for which the autocorrelation is known (decreasing slowly). here is my test code

import scipy as sp
import colorednoise as cn

def plot_acf(length, avg):
    fig, ax1 = plt.subplots(1,1)

# get 1/f^2 noise.
    data = cn.powerlaw_psd_gaussian(2, length)

# getting the power-spectrum. scipy.welch returns the power-spectrum, with
# positive frequencies first and then negative frequencies.
    f, pxx = sp.welch(data, nperseg = avg, window = "hanning", return_onesided = False)

# remove the spurious imaginary part of the inverse fourier transform.
    acf = np.real_if_close(np.fft.ifft(pxx))

    ax1.plot(acf[:int(avg/2)])

But it produces plots that look like that :

Which is obviously wrong : the strengh of the correlation is way too weak and it decreases too fast for a $1/f^2$ noise signal.

Do you know where the poblem might come from ?


texts: Real and complex low pass equivalent of band-pass signal

texts:
A real band-pass signal $s(t)$ with center frequency $f_c$ can always be represented as

$$ s(t) = A(t) \cos(2\pi f_c t) $$

where $A(t)$ is real and can be viewed as amplitude.

Also

$$ s(t) = A(t) \cos(2\pi f_ct) = \Re\left\{(A(t) e^{j2\pi f_ct})\right\}$$

Then we consider

$$s_l(t)= A(t) e^{j2\pi f_ct}$$

to be the corresponding complex band-pass signal of $s(t)$.

My question is how can I find the real low pass equivalent of $s(t)$ and the complex low pass equivalent of $s_l(t)$?


texts: What is the Technique to Find Variance of Estimation Error

texts:
Given an $n$-vector $y$ (responses) and a design matrix $X$, I wish to fit them with a simple linear regression model
$$y=X\beta+e,$$ or,

$y_t = x_t'\beta_0 + e_t$
where $e\sim\mathcal{N}(0, \sigma^2I)$ and $\beta_0$ is the true parameter. Then, we have
$$y\sim\mathcal{N}(X\beta, \sigma^2I).$$
Then the maximum likelihood estimations (MLE) of $\beta$ which is the Ordinary Least Square estimator is
$$\hat\beta_n=(X_n^TX)^{-1}X_n^TY_n$$. Then how to calculate the Expectation of the variance of the estimation error viz. $E[(\hat{\beta_n} - \beta_0){(\hat{\beta_n} - \beta_0)}^T]$ . I am struggling to understand what to substitute for $\beta_0$ in this calculation and how to find the expression. Any help will be extremely useful. Thank you.


texts: Denoising short, non-uniformly spaced, bandlmited sequences

texts:
Consider a bandlimited signal $x(t)$ with bandwidth $BW$.  Samples of this signal are observed at non-uniformly spaced sample points in the presence of noise:

$y(t_i)=x(t_i) + n_i \qquad i\in\{0,1,\cdots,N-1\}$

where $N$ is on the order of 10-30 samples.

The samples are "densely-spaced" in the sense that:

$\frac{1}{t_i - t_{i-1}} > 2\times BW$

The noise samples can be modeled as a realization of independent, identically distributed random variables of zero mean and some variance $\sigma^2_n$.  A variant of the noise model would be the inclusion of occasional outliers.

My goal is to recover estimates $\hat{x}(t_i)$ of the desired sequence which explicitly exploit the bandlimited nature of $x(t)$.

Savitzky-Golay filters can be generalized for non-uniformly spaced samples, but it's not clear how the bandwidth can be taken into account.  (For uniformly spaced samples, such an interpretation is available, see "What Is a Savitzky-Golay Filter?" by R. Schafer).  I guess I could use the frequency response interpretation of the SG filter for a worst case effective, constant sample frequency, but I'm wondering if other approaches may be better.

An extension would be interpolation/resampling of the data on a uniform grid.

Any suggestions?


texts:
I have a question about finding the FFT duration for 802.11a preamble.

According to the standards, when the bandwith is $20 \textrm{MHz}$ and for $N=64$ we have $\Delta f = 20\textrm{MHz}/64 =312.5 \textrm{kHz}$ and it says FFT duration should be $1/\Delta f=1/(312.5 \textrm{kHz}) = 3.2 \textrm{us}$.

My question is that why FFT duration should be "1 over subcarrier spacing"?


texts: Identifying sounds patterns from wav file

texts:
I have an audio recording stored as a .wav file. My goal is to identify exact sounds patters from it using Python.
The problem I'm facing is there's a clear "ting" sound in the audio file, but it also has "other" noise apart from my expected sound. Also, the sounds I'm expecting is either "ting" or can be "ting ting" (lower intensity) or something similar.

What's the best way to identify exact peaks from the signal?
How do I distinguish between my "ting" and "ting ting" (lower intensity/amplitude compared to normal "ting")?
How do I filter out "noise" from the signal that's not what I expect?

I tried some basic Python code using librosa and https://github.com/MonsieurV/py-findpeaks#scipysignalfind_peaks, but it wasn't reliable/scalable beyond a point.


texts: Lane detection and tracking using B-snake?

texts:
I am truly sorry if my question sounds little bit stupid, but I am getting stuck at this part for a long time. I am currently working for my final year project of lane detection for curved road using B-snake. I am using Matlab software, and I have done the image pre processing part by converting the RGB value to grayscale. After the converting I simply use the Canny edge detection but it seems that there is too many noise.

This is the first part that I obtained :

Then, I tried the region props function for which the Matlab code I obtained from this website. The noise has reduced a little bit.

The only region that I want is the road region. But I don't manage to get the ROI that I need. So I just proceed to Hough Transform to obtain the straight line. I try to run the coding of Hough Transform that I found on the internet, but the result is so confusing. I don't know how to change the code so that my result will turn out better.

The image after Hough Transform is applied :

This is my full code:

I=imread('C:\Users\LENOVO\Documents\MATLAB\PSM2\road\draw\test1.jpg'); %         Input image file here
dim = size(I); % Get size of the image and store in an array called 'dim'.%size array matrix; 1 is height, 2 is width

I2 = rgb2gray(I); % Convert to grayscale
BW = edge(I2,'canny',[]); % Edge detection using Canny
cc = bwconncomp(BW);
I4 = labelmatrix(cc);

a_rp = regionprops(cc,'Area','MajorAxisLength','MinorAxislength','Orientation','PixelList','Eccentricity');
idx = ([a_rp.Eccentricity] > 0.99 & [a_rp.Area] > 100 & [a_rp.Orientation] < 60 & [a_rp.Orientation] > -90);

BW2 = ismember(I4,find(idx));
[H,T,R] = hough(BW2);

%  figure, imshow(H,[], 'XData', T, 'YData', R, 'InitialMagnification', 'fit');
%  xlabel('\theta'), ylabel('\rho');
axis on, axis normal, hold on;
P  = houghpeaks(H,50,'threshold',ceil(0.1*max(H(:))));

% Set houghpeaks parameters, threshold unsure
x = T(P(:,2));
y = R(P(:,1));
plot(x,y,'s','color','white');

% Apply median filtering
I3 = medfilt2(I2);

% Find lines and plot them
lines = houghlines(BW,T,R,P,'FillGap',20,'MinLength',10);
figure, imshow(I3),imagesc(I3), hold on
max_len = 0;

for k = 1:length(lines)
    xy = [lines(k).point1; lines(k).point2];
    plot(xy(:,1),xy(:,2),'LineWidth',2,'Color','green');

    % plot beginnings and ends of lines
    plot(xy(1,1),xy(1,2),'x','LineWidth',2,'Color','yellow');
    plot(xy(2,1),xy(2,2),'x','LineWidth',2,'Color','red');
end

showlines = struct(lines);
cellData = struct2cell(showlines);

% X-coordinates are for width
% Y-coordinates are for height
%point1(x y) etc

for i = 1:280
    % 'A' stores all 'x' coordinates of point 1
    A([i,i+1])= [cellData{1,i}];
    % 'B' stores all 'x' coordinates of point 2
    B([i,i+1])= [cellData{2,i}];
    % 'C' stores all 'y' coordinates of point 1
    C([i,i])= [cellData{1,i}];
    % 'D' stores all 'y' coordinates of point 2
    D([i,i])= [cellData{2,i}];
end

I do really hope someone can help me. At least general idea for me to proceed to the next step is enough. Thank you in advance.


texts: How do you interpret FFT of an FFT of a discrete signal?

texts:
I am new to DSP & Radar domain. I have one question.
Say I have a discrete time signal $x[n]$ and I take an FFT of it say the output is $X[k]$. The result what i get is frequencies present in $x[n]$.

Now what will happen if I take fft on $X[k]$ ? How would the spectrum look and what would you interpret from it?

thanks in advance !


texts: What does frequency domain denote in case of images?

texts:
I was just learning about the frequency domain in images.

I can understand the frequency spectrum in case of waves. It denotes what frequencies are present in a wave. If we draw the frequency spectrum of $\cos(2\pi f t)$, we get an impulse signal at $-f$ and $+f$. And we can use corresponding filters to extract particular information.

But what does frequency spectrum means in case of images? When we take the FFT of a image in OpenCV, we get a weird picture. What does this image denote? And what is its application?

I read some books, but they give a lot of mathematical equations rather than the physical implication. So can anyone provide a simple explanation of the frequency domain in images with a simple application of it in image processing?


texts: Differentials - Differences: Non causality in the system

texts:
I'm still learning DSP and referring to Oppenheim video lectures.

In that lectures, differential difference equation is obtained for IIR filter design, in Lecture 14.

$$\mathcal{L}[\frac{\mathrm d}{\mathrm dx}y_a(t)] = sY_a(S)$$
$$\mathcal{Z}[\frac{y[n+1] - y[n]}{T}] = \frac{z-1}{T}Y(z)$$

And using these two equations, difference equation $z-1=sT$ is obtained.

To calculate the z-transform of the differential, $y(n+1)$ is used. By using that term, are we making our system non-causal..? (As our equation depends on future values)

Thank you for your time..


texts: Find the autocorrelation function of signal $x(t) = u(t) - u(t-1)$

texts:
I have used the energy-type signal autocorrelation function:
$$\mathcal{R}_{xx}(\tau)=\int_{-\infty}^{\infty}x(t)x^*(t+\tau)dt$$
I have rewritten the equation as:
$$\begin{align}
\int_{-\infty}^{\infty}\big[u(t)-u(t-1)\big]\big[u(t+\tau )-u(t+\tau-1)\big]dt \\
\end{align}$$
How do I simplify this equation?


texts:
I have created a fairly simple TDOA system that uses ultrasonic signals emitted from two speakers to geolocate (relative to the speakers) mobile phones.  The two signals are separated by frequency.

The system has the following constraints:

The signals must be inaudible.  To that end we stick to frequencies above 17 kHz.  A few people can still hear that, but most can't.
Sample rate is 44.1 kHz.
Music will typically be playing, so there is lots of interference at the lower frequencies.
We don't have control over how well the speakers and microphones work at the upper frequencies, so we've kept our upper limit at around 20 kHz.

The particular signal that I am using is BPSK modulated 13-bit Barker codes because of their good autocorrelation properties.  The autocorrelation looks like the following-

When I cross-correlate the expected signal against the received signal in real life, though, what I get typically looks like this-

The blue is the cross correlation with the speaker 1 signal, and the red is the cross-correlation with the speaker 2 signal.  It appears that the echoes are significant and, unfortunately, often stronger than the direct path signal due to the directional gain of the microphone.

I tried simply detecting the earliest appearance of the signal as that is likely to be the direct path.  This approach is very sensitive to the threshold that I use for deciding when the signal is present and so is not robust at all.

I would like a robust approach for determining the "true" arrival time of the signal- i.e. the arrival time of the direct path signal.  Perhaps some form of channel estimation and deconvolution?  If so, how would that work?

Data/Code: I want to make it clear that I am not expecting anyone to analyze the data or inspect my code.  I have made them available in case you want to do so.  I am mostly interested in ideas.

I made the raw received signal and modulated expected signals available for download.  They are all sampled at 44.1 kHz.  Correlating the received signal with the expected signals will produce something similar but not identical to the picture above because I move the received signals to baseband and decimate before correlating with the expected signals.

Received signal

Expected signal #1

Expected signal #2

Matlab scripts
The Matlab scripts has both the signal generation script (genLocationSig.m) and my receive/processing script (calcTimingOffset.m).


texts: How to make the error of the Gardner timing error detector not depend on signal amplitude?

texts:
I have implemented a Gardner timing error detector with an interpolator and timing recovery in software. I noticed that the error term of the Gardner TED depends on the signal amplitude, which causes problems with the large variations of signal strength my receiver is expected to see.
With bursts (no timing sequence, only CW + sync word over 64 symbols), there are variations of 10 dB from the start of the transmission to after a few (15-20) symbols of frame data. For continuous signals (long CW + long timing sequence), this is less of an issue as the signal level tends to be relatively constant or varying rather slowly.
Is there a way I can make the Gardner error factor not depend on signal amplitude or should I instead focus on making sure the TED always sees a constant power level?


texts: Adjusting coefficients of a LP and HP filter for a new sampling rate

texts:
I am pretty new to DSP and try to undertand some concepts and the relations between them.

I currently learning about the Pan-Tompkins algorithm and trying to recreate it in Matlab (Pan–Tompkins algorithm). The sampling rate for this algorithm was 200Hz but I am working with data sampled in 512Hz.
The coefficients for the filters are:

LP:
b = [1 0 0 0 0 0 -2 0 0 0 0 0 1];
a = [1 -2 1];

HP:
b = [-1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 32 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1];
a = [1 1];

my question is: how do I adjust those coefficients for a new sampling rate?

p.s: I can't work with any of the Matlab's toolkits.


texts: Performance metrics for analysing the performance of various filters

texts:
We are working with speckle noise reduction in ultrasound images by using various despeckling filtering techniques. We compare the performance of each filtering technique by using Quality metrics such as mean square error,signal to noise ratio,peak signal to noise ratio,Edge preservation index,Speckle index.

We know that mean square error metric indicates how different the images being compared are. Therefore the lower its value is,the closer the estimated image is to the original image.
The peak signal to noise ratio must be high which indicates an improvement in speckle reduction.But when we calculate these metrics, mean square error is high and Peak signal to noise ratio is small. any one please say me why this happens so.

The original image that  we have used is shown below

Then we corrupted the original image with 3% noise

The formula that we used to calculate these metrics is given below

The result which we obtain while using median filter is given below


texts: Equalization knowing the Channel State Information

texts:
I am trying to simulate a simple communication system using channel estimation. I was able to estimate the channel response using the Least Squares channel estimator, and it is working properly (evaluated by the mean square error between the channel and its estimation). But once I have this knowledge, I need to equalize the signal at the receiver. I have tried to apply the zero-forcing method, but it's not providing a good performance.

My questions are:

What is the proper way of doing this equalization?
In a system where two transmitters send information at the same time and in the same frequency, I can estimate both channels at the receiver, but how can I equalize the signal once I know both channels?
Is there is any reference that anyone could indicate to read more about this topic and help answer these questions? It wold be enough for now.

Thanks in advance.