texts: Can we really use equal gain combining in the receiver in the real world?

texts:
Can we really use equal gain combining(EGC) in the receiver in the real world?
Because we have to know the phase of elements in the channel matrix first before we do the EGC,however,the channel should be randomly generated every second,so I wonder that if we can do EGC in the real world or not?

Or we can  do the EGC,but the result is just a little bit better of not doing EGC?


texts:
Say $a_x$, $a_y$, $b_x$, $b_y$ are signals with values in $\mathbb{R}$. Say they are independednt, normally distributed variables with given SDs $\sigma_{a_x}$, $\sigma_{a_y}$, $\sigma_{b_x}$, $\sigma_{b_y}$.

Say $f: \mathbb{R}^4\longrightarrow\mathbb{R}$ is a continuous function of these signals (i.e. $f(a_x, a_y, b_x, b_y)\in \mathbb{R}$ for all $(a_x, a_y, b_x, b_y)\in\mathbb{R}^4$.

General question

What is the SD of $f(a_x, a_y, b_x, b_y)$?
Example/Question

Say $a$ and $b$ are vectors (say in $\mathbb{R}^2$) and $\sphericalangle$ is the angle between $a$ and $b$, say defined as
$$\sphericalangle(a,b) := \mathrm{arcos}\left(\frac{ab}{\vert a\rvert \lvert b\rvert}\right)\quad\text{where}\quad a=(a_x, a_y)\quad\text{and}\quad b=(b_x, b_y).$$

Say $a_x, a_y, b_x, b_y$ are independent and normally distributed, have SDs $\sigma_{a_x}$, $\sigma_{a_y}$, $\sigma_{b_x}$, $\sigma_{b_y}$ respectively.

Then Gauss error proparation yields
$$\sigma_\sphericalangle:= \sqrt{ { \left(\frac{\delta \sphericalangle}{\delta a_x} (a,b) |\mu \ \cdot \sigma_{a_x}\right) }^2 + { \left(\frac{\delta \sphericalangle}{\delta a_y} (a,b) |\mu \ \cdot \sigma_{a_y}\right) }^2 + { \left(\frac{\delta \sphericalangle}{\delta b_x} (a,b) |\mu \ \cdot \sigma_{b_x}\right) }^2 + { \left(\frac{\delta \sphericalangle}{\delta b_y} (a,b) |\mu \ \cdot \sigma_{b_y}\right) }^2}$$

as the error of $\sphericalangle$ (where $..|\mu$ shall note evaluation at average $\mu$).

Question

Is that $\sigma_\sphericalangle$ the SD of $\sphericalangle$?


texts:
The problem I am having is related to sample rate conversion and more precise to sample rate reduction. I have been working on the paper Interpolation and Decimation of Digital Signals Tutorial Review in [1] and A digital signal processing approach to interpolation in [2] where the maths make my understanding difficult.

To start with,bellow is defined a new sequence $w'(n)$

$$
w'(n)=
\begin{cases}
w(n), &\text{if}\quad n=0,\pm M, \pm 2M, \ldots \\[2ex]
0, & \text{otherwise}
\end{cases}\tag{1}
$$

that is non zero only at multiples of $M$ of the old sampled sequence. Then it is stated that a convinient representation of $w'(n)$ is:

$$
w'(n)=w(n)\left(\frac{1}{M}\sum^{M-1}_{l\ =\ 0}e^\frac{j2\pi ln}{M}\right)\tag{2}
$$

The term in the brackets corresponds to a discrete Fourier series representation of a periodic impulse train with a period of $M$ samples. That in my understanding means that is a multiplication by a $\delta(n)$ in time domain if the term in the brackets refers to the inverse DFT of one $(1)$. If that is the case then shouldn't that be a convolution?

Then an expression for elaborating on the $z$-domain is derived:

$$
 y(m)=w'(Mm)=w(Mm)
$$

Then $\mathcal Z$-transform is utilized after:

\begin{align}
Y(z)&=\sum^{+\infty}_{m=-\infty}y(m)z^{-m}\\
&=\sum^{+\infty}_{m=-\infty}w'(Mm)z^{-m}\\
&=\sum^{+\infty}_{m=-\infty}w(m)\left(\frac{1}{M}\sum^{M-1}_{l\ =\ 0}e^\frac{j2\pi lm}{M}\right)z^{-\frac{m}{M}}\\
&=\frac{1}{M}\sum^{M-1}_{l\ =\ 0}\left(\sum^{+\infty}_{m=-\infty}e^\frac{j2\pi lm}{M}z^{-\frac{m}{M}}\right)
\end{align}

and we end up here.

How to conclude the result bellow?
What are the steps missing?
Where has $m$ gone?
And why do we have the exponential as argument of $W$ in equations $(3)$?

$$
Y(z)=\frac{1}{M}  \sum^{M-1}_{l\ =\ 0}  W\left(e^\frac{j2\pi l}{M}z^{\frac{1}{M}}\right)\tag{3}
$$

[1] R.E.Crochiere and L.R.Rabiner,"Interpolation  and Decimation of Digital Signals Tutorial Review",Proceedings of the IEEE,vol.69,pp 300-331,March 1981

[2] R.W.Schafer and L.R.Rabiner,"A digital signal processing approach to interpolation",Proceedings of the IEEE,vol 61,pp 692-702,June 1973


texts:
I am confused regarding the complexity of the Fast Fourier Transform (FFT). The Discrete Fourier Transform is:

$$\qquad X\left [ k \right ]=\sum_{n=0}^{N-1}x[n]W_{N}^{kn}\quad \text{where}\quad W_{N}=e^{\frac{-j2\pi}{N}}\tag{1}$$

The first step in deriving the FFT is to express $(1)$ as the sum of even and odd inputs

$$\qquad X[k]=\sum_{n\ \rm even}x[n]W_{N}^{nk}\ +\ \sum_{n\ \rm odd}x[n]W_{N}^{nk}\tag{2}$$

And after many other derivations the final equation $(3)$ is obtained as

$$\qquad X[k]=\sum_{r=0}^{N/2-1}x[2r]W_{N/2}^{rk}\ +W_{N}^{k}\ \sum_{r=0}^{N/2-1}x[2r+1]W_{N/2}^{rk}\tag{3}$$

Equation $(3)$ consists of $2$ length $N/2$ FFTs, which is more computationally efficient than performing one big length $N$ FFT.

My question is: isn't equation $(2)$ also $2$ length $N/2$ FFTs? In other words, why all the bother about deriving equation $(3)$? Isn't equation $(2)$ already more computationally efficient than equation $(1)$? Or perhaps I'm missing something?

Any help is greatly appreciated.


texts: Calculating power in frequency band using power spectral density

texts:
In a text I am using it states:

  The fourier transform $F(\omega)$ has both positive and negative
  frequencies and $F(\omega) = F^*(-\omega)$. The power spectral
  density $S(\omega)$ is similarly defined to have the positive and
  negative frequencies, and the power spectral density is similarly
  equally divided between the positive and negative frequencies. Thus
  power spectrum $S(\omega)$ is symmetric about $\omega = 0$. So the
  contribution to the total power that is contained in the band
  $[\omega_1, \omega_2]$ is $$P_{[\omega_1, \omega_2]} = \frac{1}{2 \pi} \int_{\omega_1}^{\omega_2} [S(\omega) + S(-\omega)]d \omega.$$

Question: If $S(\omega)$ is the power per frequency $\omega$, even though $\omega$ takes both positive and negative values, I still don't see why the power contained in the band $[\omega_1, \omega_2]$ is not simply given by $$\frac{1}{2 \pi}\int_{\omega_1}^{\omega_2} S(\omega) d \omega?$$ I am only beginning to learn about these things, hence there might be something trivial that I am missing. Thanks.


texts:
I have a Masters Degree in DSP but all I learnt there is theory. I would like to complete some projects at home which I can add in my resume to showcase my skills and to get jobs in industry. But I have not found anything after googling. So can someone here please guide me what should I do?
P.S.:

I am not interested in image processing.

I would like to start projects on Matlab/Octave or Python and then writing the same code in Embedded C for a DSP hardware.


texts:
I have two sets of MFCC coeficients which I get with fftExtract utility from two wav audio files.

How can I compare this sets (matrix 271x13 and 551x13) if they are similar ?

Probably I should use EMD distance ...  But Im sligthly lost in all this math and I dont
know how to interpret its results .

How can I make from this distance matrix if both audio signals are similar or not - (Yes/Not)

I used following code for EMD :

function D = distEmd( X, Y )

Xcdf = cumsum(X,2);
Ycdf = cumsum(Y,2);

m = size(X,1);  n = size(Y,1);
mOnes = ones(1,m); D = zeros(m,n);
for i=1:n
  ycdf = Ycdf(i,:);
  ycdfRep = ycdf( mOnes, : );
  D(:,i) = sum(abs(Xcdf - ycdfRep),2);
end

fftExtract

http://digitalmusics.dartmouth.edu/wiki/index.php/Music103_Tutorial


texts:
What is the algorithm for removing background noise from an audio clip of human speech? Say I have a 10 second audio clip of somebody talking, and there is background noise of them tapping on their keyboard. Now, I would like to remove the "keyboard tapping" signal while preserving the speech signal (always assuming a single speaker).

What are the necessary steps to achieve this? Is the best approach to use voice activity detection to try and isolate a representative signal of keyboard tapping, and then "remove" that pattern from the rest of the signal? What is required to "remove" the specific sound from the signal once identified?


texts: Realization of a filter based on its transfer function

texts:
How can we check whether the filter is realizable given its transfer function and What are the parameters the realization depends on?
Here is an example:

Show that a filter with transfer function
$H(ω)=\frac{2(10^5)}{ω^2 + 10^{10}} e^{-jωt_0 }$
is unrealizable. Can this filter be made approximately realizable by
choosing a sufficiently large t0.

--Modern digital and analog communication systems (3rd edition) - B.P. Lathi
Note: The question is of general nature. The example quoted above serves only as to make it clear, what is being asked. No one has to provide a solution for this problem.


texts: Z transform - Inverse System function - Why number of poles and zeros myst be equal?

texts:
I know that if a system is causal then the system function H(z) must have :
a) a ROC that spans from the exterior of the most distant pole and
b) the number of zeros must not be greater than the number of poles

I have found this exercise:

and the solution:

Why the number of poles and zeros must be equal?


texts:
I was wondering if there was a simple approach to smooth images, e.g using a Gaussian Kernel without introducing boundary effects which are also present when using Discrete Fourier Transforms?
For both approaches, the missing pixels at the left-hand side will be taken from the right-hand side and vice versa.

Right now I have a 5x5 Gaussian Kernel which I move iteratively over my MxN image. Up to now, in case I reach the last pixel of the image, the operator automatically uses  the first pixel and so on. Since my picture is not periodic I get artifacts even if I pad my picture on all sides with the intensity values of the corresponding adjacent column or row.

Do you gus have any suggestions?

Thanks!


texts: Often we need to estimate the time difference of arrival between two signals to find the location of a target. Many algorithms gives the time delay corresponding to a sample number or time delay is a function of the sampling frequency. To get away with this problem subsample time delay estimation procedure is used. One can see that if the sampling frequency is very high the subsample may not be necessary. Can anyone help to find when subsample is necessary and when it become irrelevant? I am trying to understand the other associated factor.


texts: Digital Lowpass Filters, basic questions

texts:
Forgive the newbie questions, I'm just starting out and wanted some basic information.

I have a signal coming from a digital gyro, most of the useful information is below 90hz but I am seeing significant noise spikes from mechanical motor vibration at 200-400hz (the spike moves depending on motor RPM).  This needs to be filtered in realtime to reject the noise with as little delay as possible.

When the sampling rate is increased, say from 1000hz up to 8000hz, does this increase the delay?
When sample rate increases, is the number of samples/cpu load per timestep for the same filter increased?
In the above situation, what is a good way to increase rejection of signals above 90hz while keeping the delay low?  ie: filter type, raise/lower sampling rate...etc.
An existing system I've seen uses an LPF followed by 2 bandstop/notch filters positioned over the noise spikes.  The designers say this 3 filter approach creates less latency than doing a single more aggressive LPF, does this sound right?


texts:
I am a little confused on CDMA and its interaction with the typical pulse shaping filter ex, raised root cosine.
Is the RRC filter bandwidth selected to be consistent with the CDMA spreading/update sequence? Or is it kept the same as prior to CDMA implementation so the number of pulses is consistent with the number of transmitted bits?


texts: How do bode plots work with unstable systems work?

texts: If I had a system with right-half s-plane poles, how would a frequency response work? Since a purely imaginary value for s, would cause the Laplace transform to diverge for such a system, what meaning would a bode plot have in such a situation?


texts: Recreating matlab second order filter in python

texts:
I want to create some custom filters in Matlab filterDesigner tool and then apply those filters to a large set of audio files. I want to use python to do the application of the filters to the audios.
I created a Butterworth filter in Matlab and exported it as an ASCII file. In there I get two datas. SOS Matrix and Scale Values. Then I read the file in python and got the SOS matrix and filtered my signal using the result as follows (using scipy.sosfilt)
from scipy import signal

x_sos = signal.sosfilt(sos, x)  # x is my signal

The result I get has the desired shape except for the amplitude being different (I applied the same filter in Matlab and compared the output from the python code). I understand that scale values needs to come into this equation. But how do I incorporate that? I don't understand this answer. This question has answers explaining the scale values but still I have no idea how I should use that in the python code

EDIT: (After getting the answer) To anyone who is interested in the final code, here it is. And there's an explanation again in this thread
from scipy import signal
import numpy as np

x_sos = signal.sosfilt(sos, x) * np.prod(scale)  # x is my signal


texts: Is there something as self-completing data transfer, like Sudoku?

texts:
I would like to know, since I am not that knowledgeable in this subject, if this idea is feasible or widely implemented:

If a data package is sent, I know parity-bytes are used to detect errors. Could we not use a similar concept to self-complete data packages?
Let's say a client requests an image to be sent from a server. The server first sents a function of each row and column, e.g. the sum of all pixel values in each row and column, so 20 values for a 10x10 image. The server then sends the values of exactly as many pixels as are needed to  let only one possible combination of all pixel values remain.

Would this require less or more data packages to be sent, in the best/average/worst case?


texts: Limitation of Wiener filter in image denoising

texts:
In image processing, up to which value of mean and variance can a Wiener filter be used to denoise the image?

Wiener filtering can be used when variance of the noisy image is low, but I don't know the specific values of variance and mean.

What are limitations of the Wiener filter?


texts: Understanding Spatial Frequency of Images (using matlab)

texts:
I wish to isolate a crosswalk from its surrounding, and I thought that as crosswalks have a repetitive pattern of thick white stripes, I might be able to isolate them by measuring their frequency range, and using a band pass filter with the right cut-off frequencies.

The problem is that I have no idea how to analyze the frequency of patterns in images. I know I should run fft2 on an isolated crosswalk, but I don't know what to do with the result.


texts: How to get the signal to noise ratio from noisy data set with sinusoid?

texts:
I'm new to DSP and trying to do some example problems for a class and am confused by this one.

Given a data set, and knowing it contains a signal that is a single sinusoid with random content noise, how can I calculate the signal to noise ratio? I also have the sampling rate.

I have been able to plot the autocorrelation, and can see the characteristic peak at 0 to indicate the noise. I'm just unsure if I can get the power of the signal and noise from this autocorrelation, or the power spectra density (taking the fft of the autocorrelation). Do I need to find the amplitude, frequency, and phase of the sinusoidal signal and subtract it from the original to do this? If so, how do I find those?

Any help is appreciated, I am new to DSP and learning to thanks!


texts: Combining Analog Video Signals for Diversity Receiver System

texts:
The system I want to build takes an analog video from two receivers (that are receiving from the same transmitter) and then combines the signals using equal gain combining (or maximal gain combining if a simple circuit exists to variably attenuator each signal).
My question is, how is this method practically implemented?
The best I've been able to find is that each analog signal is converted to digital and then has the gain adjusted before combining.

Is the combining of the signal usually done before or after the receiver?
Would the method of converting to digital and scaling gain introduce too much noise to the signal in practice?
What kind of ADC and variable amplifier would be needed in this application (NTSC composite video)?

Thanks for any help you can provide.


texts:
The following sequence is given which is supposed to be time-variant:
$$y[n] = \sum_{k=n_0}^n x[k]$$
I'm having difficulties proving the time-variance or finding a counterexample for it being time-invariant.
My idea (which proves it being time-invariant?) is:
\begin{align}y_1[n] &= T\{x_1[n]\}\\x_2[n] &= x_1[n-n_0]\end{align}
check, whether    $T\{x_1[n]\} = y_1[n-n_0]$
$$y_2[k] = \sum_{k=n_0}^n x_2[k] = \sum_{k=n_0}^n x_1[k - n_0] = \sum_{k=n_0-n_0}^{n-n_0} x_1[k] = y_1[n-n_0]$$
Where is my mistake? Working around the summation is giving me trouble.
EDIT:
Following the advice in the comments for $n = n_0$ the left equation equals to 1 while the right equation equals to 0, making both equations unequal and therefore the sequence time-variant?
Second sequence
This sequence is supposed to be time-invariant.
\begin{align}y[n] &= \sum_{k=n-n_0}^{n+n_0} x[k]\\
y_2[k] &= \sum_{k=n-n_0}^{n+n_0} x_2[k] = \sum_{k=n-n_0}^{n+n_0} x_1[k - n_0] = \sum_{k=n_0-2n_0}^{n} x_1[k] = y_1[n-n_0]\end{align}
I followed the same steps but I'm still not sure about for the case $n = n_0$ since when added in left and right equation
$$\sum_{k=0}^{2n} = \sum_{k=-n}^{n}$$
They seem unequal to me?


texts: Is match filter a solution of some differential equation?

texts: Matched filter is very common and termed as optimal filter for detection purpose. The question is : can we formulate a differential equation and show that the solution is matched filter.


texts: recovery plane from 4 point imaging using calibration camera

texts:
I have a camera and its K matrix (calibration matrix) also I have image of plane, I know the real points of the 4 corners and thier correspondence pixel. I know how to compute the H matrix if z=0 (H is homography matrix between Image and the real plane).
And Now I try to get the real point of the plane (3D point) with the rotation matrix and the transltion vector
I follow this paper :Calibrating an Overhead Video Camera by Raul Rojas in section 3 - 3.3.
My code is:

ImagePointsScreen=[16,8,1;505,55,1;505,248,1;44,301,1;];

screenImage=imread( 'screen.jpg');
RealPointsMirror=[0,0,1;9,0,1;9,6,1;0,6,1]; %Mirror
RealPointsScreen=[0,0,1;47.5,0,1;47.5,20,1;0,20,1];%Screen
imagesc(screenImage);
hold on
for i=1:4
    drawBubble(ImagePointsScreen(i,1),ImagePointsScreen(i,2),1,'g',int2str(i),'r')
end

Points3DScreen=Get3DpointSurface(RealPointsScreen,ImagePointsScreen,'Screen');

figure
hold on
plot3(Points3DScreen(:,1),Points3DScreen(:,2),Points3DScreen(:,3));
for i=1:4
    drawBubble(Points3DScreen(i,1),Points3DScreen(i,2),1,'g',int2str(i),'r')
end

function [ Points3D ] = Get3DpointSurface( RealPoints,ImagePoints,name)
M=zeros(8,9);

for i=1:4

M((i*2)-1,1:3)=-RealPoints(i,:);
M((i*2)-1,7:9)=RealPoints(i,:)*ImagePoints(i,1);
M(i*2,4:6)=-RealPoints(i,:);
M(i*2,7:9)=RealPoints(i,:)*ImagePoints(i,2);

end

[U S V] = svd(M);
X = V(:,end);
H(1,:)=X(1:3,1)';
H(2,:)=X(4:6,1)';
H(3,:)=X(7:9,1)';
K=[680.561906875074,0,360.536967117290;0,682.250270165388,249.568615725655;0,0,1;];

newRO=pinv(K)*H;
h1=newRO(1:3,1);
h2=newRO(1:3,2);

scaleFactor=(norm(h1)+norm(h2))/2;
newRO=newRO./scaleFactor;
r1=newRO(1:3,1);
r2=newRO(1:3,2);
r3=cross(r1,r2);
r3=r3/norm(r3);

R=[r1,r2,r3];

RInv=pinv(R);
O=-RInv*newRO(1:3,3);
M=K*[R,-R*O];
for i=1:4
   res=pinv(M)* [ImagePoints(i,1),ImagePoints(i,2),1]';
   res=res';
   res=res*(1/res(1,4));
   Points3D(i,:)=res';

end
Points3D(i+1,:)=Points3D(1,:);  %just add the first point to the end of the array for draw square

end

function drawBubble (xcenter, ycenter, radius, color, string,colorNum)
% draw circle with input xcenter, ycenter, radius, color and text
angles = linspace(0, 2*pi, 40);
bubx = xcenter + (radius * cos(angles));
buby = ycenter + (radius * sin(angles));
fill(bubx, buby, color, 'EdgeColor', color)
text(xcenter, ycenter, string, 'Color', colorNum, 'FontSize', 13, ...
       'FontWeight', 'normal', 'HorizontalAlignment', 'center')

end

My problem is that I don't get good results

1.The point 1 is at (0,0,0) and this is not the real location

2.the points are upside down

What I am doing worng?


texts:
I'm having 2 signals which have different sampling rates i.e., 1ms(A) and 4ms(B) and I've tried to upsample/downsample either of the signals based on the code snippet Resampling based on FFT which I believe is equivalent of what is descriped in How to resample audio using fft or dft

% FFTRESAMPLE Resample a real signal by the ratio p/q
function y = fftResample(x,p,q)

% --- take FFT of signal ---
f = fft(x);

% --- resize in the FFT domain ---
% add/remove the highest frequency components such that len(f) = len2
len1 = length(f);
len2 = round(len1*p/q);
lby2 = 1+len1/2;
if len2 < len1
  % remove some high frequency samples
  d = len1-len2;
  f = f(:);
  f(floor(lby2-(d-1)/2:lby2+(d-1)/2)) = [];
elseif len2 > len1
  % add some high frequency zero samples
  d = len2-len1;
  lby2 = floor(lby2);
  f = f(:);
  f = [f(1:lby2); zeros(d,1); f(lby2+1:end)];
end

% --- take the IFFT ---
% odd number of sample removal/addition may make the FFT of a real signal
% asymmetric and artificially introduce imaginary components - we take the
% real value of the IFFT to remove these, at the cost of not being able to
% reample complex signals
y = real(ifft(f));

i.e., reimplemented the fftResample in C++ using FFTW libraries.

The problem is that the output is distorted because of the zero-padding.
I'm relatively new to Singal processing and I would appreciate if users can help me understand:

if I'm using the right code snippet to upsample/downsample?
highlight the benefits of upsampling(A to 4ms) compared to
    downsampling(B to 1ms)?


texts: What are some good signal and image processing blogs worth following?

texts: Are there some good sites or blogs where I can keep myself updated on the latest news and papers about image and signal processing research, or I should just check out "classical" providers like IEEE Transactions, Elsevier, etc?


texts: Cross-correlation of maximum length sequence with noise

texts:
Background

I am using a maximum length sequence to find the impulse response of a linear system. I want to calculate the expected SNR. I know that my signal amplitude is scaled by the length of the sequence after cross-correlation, but I am not sure what happens to noise.

Idea

Since the sign of noise is random, the sign of an element in a maximum length sequence doesn't matter: [1,1,-1] may as well be [1,1,1]. Cross-correlation with a maximum length sequence then is like using a moving average filter except the coefficients are all 1's instead of 1/n where n is length of the filter. A moving average filter reduces noise power by the length of filter, and scaling by the length of the filter increases noise power by the length of the filter squared. Therefore the noise power increases by the length of the filter.

Questions

What is the cross-correlation of a maximum length sequence and thermal noise?

What is the cross-correlation of a maximum length sequence with a discrete sine (e.g. EMI from a clock signal)?


texts: How can I vectorize the computations for a first-order recursive filter?

texts:
I have a simple single pole low pass filter (for parameter smoothing) that can be explained by the following formula:

$$
y[n] = (1-a) y[n-1] + a x[n]
$$

The architecture that I'm using has access to single-instruction, multiple-data (SIMD) instructions that can perform multiple vectorized calculations in parallel. I would like to take advantage of this capability, but I'm not sure how to do that for a recursive filter like this. The problem is that every computation needs a previous result.


texts: What does the "smoothing optical flow" mean?

texts:
I have this algorithm for Correlation Based Pyramids of images:

• Create two Gaussian pyramids from the two input
images

• Compute optical flow using “5×5” regions on
smallest pyramid level.

• Smooth optical flow, and use as initial guess for
higher resolution.

• Continue with next level. Search close to guess
from higher resolution.

What is the third step asking me to do?


texts: Sinusoidal tone in Signal-to-Noise measurements

texts:
Mengapa kita menggunakan nada uji sinusoidal saat mengukur signal-to-noise pada saluran bicara daripada hanya sinyal ucapan itu sendiri ?
EDIT: Translation from Indonesian according to Google:
Why do we use a sinusoidal test tone when measuring the signal-to-noise of the speech channel rather than just the speech signal itself?


texts: why is my inbuilt function is not matching with the db5 CA and CD coefficients generated by my code

texts:
I am working on wavelet decomposition.As read by the theory i have implemented wavelet decomposition in MATLAB. CA and CD for wavelet decomposition of db1 matches with the inbuilt function but not with db5 CA and CD coefficients why?

My code-

f=10;
fs=200;
amp_x = 14;
amplitude for sinusodial 2
amp_y = 9;
Time vector
t=(0:1:500);
Create a sine wave
x = amp_x * (sin(2*pi*f/fs*t))+10*rand(1,501);

LO_D =[ 0.7071    0.7071];
HI_D =[-0.7071    0.7071];

C = conv(LO_D,x);
k = 2;
y= C(2:k:length(C));

D = conv(HI_D,x);
g=2;
y1= D(2:g:length(D));

[CA,CD] = dwt(x,'db1');

y and y1 matches with CA and CD coefficients of inbuilt function.

The same code if i change it to db5 im not getting matching results .please help me out

f=10;
fs=200;
amp_x = 14;
amplitude for sinusodial 2
amp_y = 9;
Time vector
t=(0:1:500);
Create a sine wave
x = amp_x * (sin(2*pi*f/fs*t))+10*rand(1,501);
%low pass and highpass filter coefficients of db5 using wfilters
LO_D =[  0.0033   -0.0126   -0.0062    0.0776   -0.0322   -0.2423    0.1384 0.7243    0.6038    0.1601];
HI_D =[ -0.1601    0.6038   -0.7243    0.1384    0.2423   -0.0322   -0.0776   -0.0062    0.0126    0.0033];

C = conv(LO_D,x);
k = 2;
y= C(2:k:length(C));

D = conv(HI_D,x);
g=2;
y1= D(2:g:length(D));

[CA,CD] = dwt(x,'db5');


texts:
A bit stream is mapped to symbols during baseband modulation, each symbol (say for QPSK) represents two bits. Each symbol has a corresponding I & Q value, which in QPSK, are used as the amplitude for a cosine and sine carrier wave, which are then summed for the QPSK passband signal. Where does the baud rate come into play? If it is the rate at which the bit stream is processed and modulating the carrier, how do I account for baud rate if I just have a bit vector converted into a complex vector of symbols, in say, Matlab?

When I have a receiver digitizing the passband signal I receive into I/Q samples, do these values correspond to the original I/Q values representing the symbol? If the sample rate 4x the baud rate, I should receive 4:1 received samples:transmitted symbols? I would then down sample into 1:1 Rx:Tx symbols, then map the symbols back to bits using the Tx constellation (assuming no rotation)?


texts: How is the energy conservation principle preserved in the case mentioned below?

texts:
Suppose there are 2 sinusoidal signals $\cos(\omega t)$ and $\cos(\omega t)$, which implies that their individual power is $1/2$ unit and total power is $1$ unit. When added, the signal $2\cos(\omega t)$, however, has power $2$ unit.

Obviously,

$$\rm Energy_{signal} = Power_{signal}\times Time_{period}$$

So, how does the principle of energy conservation hold in this case?


texts: Explaining a statement regarding aliasing in a document

texts:
I am reading this document and in the section "Aliasing and Anti-Aliasing", the author make the following statement:

  Consider a signal with two frequency components: $f_1 = 10\,\text{Hz}$ which is
  our desired signal and $f_2 = 20\,\text{Hz}$ which is noise. Let's say that we
  sample the signal at  $f_s = 30\,\text{Hz}$. The first frequency component,
  $f_1$, will generate the following frequency components at the output
  of the multiplier, $10\,\text{Hz}, 20\,\text{Hz}, 40\,\text{Hz}, 50\,\text{Hz}, 70\,\text{Hz}$ and so on. The
  second frequency component $f_2 = 20\,\text{Hz}$ will generate the following
  frequency components at the output of the multiplier, $20\,\text{Hz}, 10\,\text{Hz},
> 50\,\text{Hz}, 40\,\text{Hz}, 80\,\text{Hz}$ and so on.

It's probably obvious to most but could someone explain where those generated frequency components come from. Thanks.

The background is that I've started reading Richard Lyons' DSP text and I really like it so far except that now I hit the bandpass sampling section and I didn't follow it. So, I've been looking for other material. If anyone knows of some good material ( books or papers ) on bandpass sampling ( I've printed out the dsp.stackexchange material and will check that out next ), it's appreciated.


texts: Given a noise-free channel with limited bandwidth, why is the channel capacity infinite?

texts:
In the Shannon Hartley theorem, if N = 0 but B is limited, the theorem tells us that C is infinite.  Why is this intuitively true?

If bandwidth is limited, I can only post so many pulses per second over the channel.


texts: How does gaussian convution stacking work?

texts:
If I first convolve an image with a Gaussian filter with σ = 1, and then convolve
the output of that image with a Gaussian with σ = 2

Does this give the same effect as convolving the original image with σ = 3?


texts: Frequency spectrum of complex equivalent Baseband Signal - Frequency Content around carrier frequency

texts:
Basically my objective is to generate a complex equivalent baseband signal of an AM modulated wave in MATLAB. Then to view the difference in frequency spectrum of AM modulated RF signal vs complex equivalent baseband of (AM modulated RF signal).

I generated an equivalent complex baseband signal for an AM modulated signal and when viewing the spectra of complex baseband signal, I see frequency content centered at 0 but too many frequency content around carrier frequency as well.

Could anyone let me know why I am seeing frequency content around carrier frequency as well?

Here are the steps I followed:

Generated a carrier signal and message signal.
Did AM modulation.
Generated Inphase and Quadrature components.
Then represented complex baseband signal by $x_{bb}(t)=x_i(t)+jx_q(t)$
View the frequency spectra of AM modulated wave.
View the frequency spectra of Complex baseband signal.

When I compare the spectra, I see that the spectra of complex baseband signal is centered around 0 but has frequency components at the carrier frequency range (both positive and negative). Here are the screenshots:

Spectra - AM modulated signal

Spectra - Equivalent complex baseband signal

Here is my MATLAB code:

fs=15000;
t=0:1/fs:1;
fm=5;
m=sin(2*pi*t*fm);  %message signal
%plot(t,m,'b');
fx=1000;
x=cos(2*pi*t*fx); %carrier signal
xm=m.*x;  %AM
ax=subplot(2,1,1);
plot(ax,t,m,'b');
axm=subplot(2,1,2);
plot(axm,t,xm,'g');
hx=imag(hilbert(xm)); %hilbert transform

xit=(xm.*cos(2*pi*fx*t))+(hx.*cos(2*pi*fx*t));
xqt=(hx.*cos(2*pi*fx*t))-(xm.*sin(2*pi*fx*t));

xbb=xit+(1i*xqt); %equivalent complex baseband signal
figure(2)
plot(t,xbb);

%spectrum
nfft=1024;
F=fftshift(fft(xm,nfft));      %spectra of AM modulated wave
fvals=fs*(-nfft/2:nfft/2-1)/nfft;
figure(4);
stem(fvals,abs(F),'b');

%spectrum
nfft=1024;
F=fftshift(fft(xbb,nfft));  %spectra of equivalent complex baseband
fvals=fs*(-nfft/2:nfft/2-1)/nfft;
figure(5);
stem(fvals,abs(F),'b');


texts: Evaluating regular Fourier Transform from DFT

texts:
Suppose we know the DFT of a discrete limited sequence, some $X[k],\  k = 0, 1,\dots ,N-1$. How can we calculate the Fourier Transform of the same signal for a random frequency $\Omega$?

EDIT:

How about this, we can get the starting signal with inverse DFT, like: $$x[n] = \frac{1}{N}\sum\limits_{k=0}^{N-1} X[k]e^{j\frac{2\pi}{N}kn}$$
and then get the Fourier Transform normally with:
$$\begin{align}
X(e^{j\Omega}) &= \sum\limits_{n = -\infty}^{+\infty} x[n]e^{-j\Omega n}\\ &= \frac{1}{N}\sum\limits_{n = -\infty}^{+\infty} e^{-j\Omega n}\sum\limits_{k=0}^{N-1}X[k]e^{j\frac{2\pi}{N}kn} \end{align}$$

Is this correct?

EDIT2:
As @Matt L. correctly states, I want the Discrete-time Fourier transform rather than the normal Fourier Transform, of course. In my native-language textbook Fourier Transform means both FT for continuous signals and DTFT for discrete signals. I hoped that that was clear from the edit. :)


texts: FFT of a Square wave, what are frequency of each Harmonic here

texts:
After getting suggestion to first study FFT on square wave to understand FFT of discrete signal, I'm trying to understand the FFT of square wave.

I've coded a program, here is the details,
Frequency of Wave, 100KHz

f= 100e3;

Sampling Frequency, 1MHz

Fs = 1e6;
Ts = 1/Fs;   %Sampling Rate
% Ts = 0.000001

t = 1/f;       %Time period of 1 Oscillation = 1/f

% t = 0.00001
n = 0.000001:Ts:t;  %Generating Samples

% t/Ts = 10

i.e 10 Samples will be generated

% 1st  Sample at 0.000001s
% 2nd  Sample at 0.000002s
% 3rd  Sample at 0.000003s
% 4th  Sample at 0.000004s
% 5th  Sample at 0.000005s
% 6th  Sample at 0.000006s
% 7th  Sample at 0.000007s
% 8th  Sample at 0.000008s
% 9th  Sample at 0.000009s
% 10th Sample at 0.00001s

defining and plotting the square wave,
Since Amplitude is not defined it will be 1 by default.

x=square(2*pi*f*n);
subplot(2,1 ,1)
stem(n,x)

Output wave,

Now plotting its FFT

subplot(2,1 ,2)
F=fft(x)
plot(2*abs(F)./(t/Ts));

Output wave,

For the First point,

Fs / # of Samples

1000,000/10 = 100KHz <- First point = 100KHz

Number of  Harmonics should be = 10 but What are the Frequency of each Harmonic ?
also, why the amplitude has changed to 1.2 which is supposed to be 1


texts: Signal analysis: selection of oscillatory parts

texts:
I'm trying to analyze the following signal:

My aim is to select the high-amlitude oscillation parts limited by the red stars.
At the beginning, i have tried to set a threshold and keep only the upper parts. however,a part that the threshold depend on the signal itself, this method didn't work well and it allows the selection of the non desired parts.
does anyone have another idea?

Thanks in advance.


texts: How to segment phone call audio into silence/non silence?

texts:
My problem is that I don't know the energy of the background noise, so I can't just threshold the energy. The processing is done in real time, and I have about 500msec to decide.
Ideally, I'd want quiet consonants considered non-silence.


texts: Calculation of the correlation of two sinusoidals

texts:
First of all I apologize for the long question. I am working with a 3-axis accelerometer and a 3-axis gyroscope and I need to calculate the correlation between the different signals (acceleration in x,y,z and the quaternion x,y,z,w). In the end, I need to get 21 correlations, one for each pair. I was thinking of using the standard Pearson product moment correlation coefficient. Particularly the second formula in the "for a sample" section that uses the standard score, sample mean, and sample standard deviation.

Now, I am computing these correlations as features for a classifier. I am computing these features over a window size. I will now introduce some MATLAB code that will help illustrate what I am trying to do.

First we have a discrete sine and a cosine.

t = [0:1:600]; %time samples
f = 1/2; %input signal frequency
fs = 100; %sampling frequency
x = sin(2*pi*f/fs*t); %generate sine
figure(1);
stem(t,x,'r');
figure(2);
stem(t*1/fs*1000,x,'r');
hold on;
plot(t*1/fs*1000,x);
y = cos(2*pi*f/fs*t); %generate cosine
figure(3);
stem(t,y,'r');
figure(4);
stem(t*1/fs*1000,y,'r');
hold on;
plot(t*1/fs*1000,y);

Now I get Pearson's coefficient for a window of 50 measurements (notice this is equivalent to a quarter of a cycle) with the following code. This can easily be changed for windows of different sizes but I will leave at the moment in 50 for the sake of simplicity.

xMean = sum(x(1:50))/50;
yMean = sum(y(1:50))/50;
tmpX = x(1:50) - xMean;
tmpY = y(1:50) - yMean;
stdDevX = sqrt((sum(tmpX .* tmpX)) / 49);
stdDevY = sqrt((sum(tmpY .* tmpY)) / 49);
stdScoreX = tmpX ./ stdDevX;
stdScoreY = tmpY ./ stdDevY;
stdScores = stdScoreX .* stdScoreY;
r = sum(stdScores)/49

This gives an r(Pearson's coefficient) equivalent to -0.9183. This makes sense from the scatter plot of x vs y:

I get the scatter plot with this:

figure(5);
scatter(x(1:50),y(1:50));

Different windows have different Pearson's coefficients, either a value close to 1 or a value close to -1, which makes sense from the scatter plots of x vs y in the other windows. The second window has an r = 0.9183. The third window will have a value again close to -1. This windows are with no overlap whatsoever.

Now, if I change the window size to 100 samples, things change.

Same story with a window of size 200, but now it looks like a circle. Now it is visible (and calculations also confirmed it) that the correlation is 0 or a value very close to 0.

So my question after all this is… what is the correct way of calculating a correlation for this type of periodic signals? Or are both of these ways correct? Preferrably I would like to obtain a number close to -1, 0, or 1 (discrete) since these numbers will be input features in a classifier. I know that, for example, two sinusoidal functions that are in phase will have a correlation close to 1 (or 1, depending on noise but that is another issue). My main problem is the correlation in the example I showed. So, probably I should calculate over a window size that is the length of at least a period in order to get either a 0 or a 1, but I am not entirely sure.

Thanks!


texts: Frequency response of numerical differentiation by polynomial interpolation / finite difference

texts:
One can use polynomial interpolation (or finite difference) to do numerical differentiation. However, there seems to be a surprising lack of interest in obtaining frequency response of this numerical differentiation method. Is there any reference I can look up to, or can anyone provide what frequency response would be? I want to obtain frequency response of any filter order (any number of filter taps).

By higher order, I mean higher-order differences as in https://en.wikipedia.org/wiki/Finite_difference#Higher-order_differences or in case of polynomial interpolation, polynomial degree.


texts: Power Spectral Density of Digital Modulations

texts:
I want to simulate the power spectral densities of different modulations in MATLAB (without Simulink). Unfortunately I have only very basic understanding of signal processing, hence I struggle with this topic. Currently I have a very simple spectrum like this:
 There are different objects (A and B in the picture above) sending signals on different frequencies. I am not simulating any bits or calculating a sine, I am just assuming that there is something going on on a certain frequency with a certain power level.

Now, I want to simulate different digital modulations, e.g. PSK, BPSK, GMSK, FSK. I understand that the spectrum should look something like this:

There is not a single center frequency, but also side bands. A and B shall illustrate different modulations. My question now is, how can I calculate these different modulations? In text books I may find figures displaying the 'power spectral density' of ASK, but often times there is no formula how to calculate this. I guess the spectral density depends on the actual bit sequence, but I just need the envelope.

In one text book ("Digital Modulation Techniques" by F. Xiong) I found this formula for ASK:
  $$\Phi_s(f) = {A^2T \over 4} \mathrm{sinc}^2\big[T(f - f_c)] + {A^2 \over 4} \delta(f - f_c)\big]$$
where $T$ is bit duration, $A$ is carrier amplitude, $f_c$ is carrier frequency.
At first I found this formula promising, because it does not require a bit sequence, but I do not understand the last addend.

Are there any formulas like this for other modulation techniques?
If not, how could I calculate these graphs?


texts: I am studying windowed Fourier ridges for my project. I know that it computes the instantaneous frequencies and amplitudes from the spectrogram. But I don't really know how it relates to signal processing, can anyone explain this concept to me?


texts: Identifying steps in walking accelerometer

texts:
I have some recorded walking data from an android phone and want to determine when someones foot hits the ground using the acceleration data. I'm making an assumption here that the person will "bob" vertically as they walk (i.e. when their foot pushes off their body will rise, and as a foot is approaching the ground their body will fall)

Given that, I have taken the $X$-axis acceleration data (which corresponds to the persons vertical axis), removed gravity, and filtered out the noise giving me something like this:

My question is: where in this signal does a foot hit the ground?

I have 3 options:

Would it be all of the valleys in the signal?

This would be correct if I were measuring position information but acceleration is a bit different. I think the marked points are actually the moments of maximum deceleration (or, acceleration towards the ground), so at these times the body is still in motion and not stationary

Would it be the points where acceleration is 0?

0 acceleration would mean the person/device is moving at a constant speed, which I believe would be the case upon foot landing because I'm assuming the speed is 0 and therefore not accelerating. But this would also be true on the "rising" part of the walk cycle - as a person pushes off the ground they will accelerate in +X, but once they reach the peak of the rise they would no longer be accelerating and you would probably see a 0 again

Would it be the 0 acceleration points, but only after a valley?

A valley would, I think, represent point of max acceleration towards the ground (i.e. theyre on the way down), so a 0 acceleration value after a valley would presumably be once theyve stopped accelerating downwards and the foot has landed?

I'm thinking the answer is #3 but wanted to see if my logic is correct, or whether my interpretation of acceleration graphs is completely wrong

EDIT: this data is recorded using a chest harness carrying the phone, so the phone is not on any of the limbs. The main assumption being made is that the torso will move up and down with each step. Also, in the phones coordinate system positive X is pointing towards the persons head

EDIT2: Using a low pass butterworth filter, instead of a bandpass, to eliminate any high frequency noise. The blue vertical lines are where I visually observe the foot hitting the ground so I can use it as a guide (even though it will be off by a few milliseconds):

And here is the filtered signal (red) over the raw signal (black):


texts: Matched Filters for attenuation detection?

texts:
I have been reading about matched filters and I'm interested in possible uses for them,

I have seen they are used very often in radar, in such a way that you send out a signal, and then you read the incoming signals to look for the reflection of the pulse you sent out using matched filters.

I think I'm right in saying that by the time difference from when you sent out the signal to when you found it again, you can detect how far the object is from you.

What I am interested to know is, can you detect the attenuation of the signal from matched filters and not just the time delay?

Hence giving you some information on the material properties of the reflecting object


texts:
In wavelet packet image compression, different types of entropy methods can be used, like Shannon and log-energy.

Shannon entropy uses this equation $\mathrm{ent}= -\sum (x^2 \times \log(x^2))$, where log-energy equation is $\mathrm{ent}=\sum(\log x^2)$, what this equation means?
which technique can achieve better results for images (in terms of PSNR)? and how can I determine the appropriate technique?


texts: How to Combine Low and High Frequencies of Two Images in MATLAB?

texts:
I would like to combine two images A and B in the following way:

1) I want to take a Fourier transform of both of them

2) For image A I want to apply a weighted filter, which gives more emphasis for low frequencies

3) For image B I want to apply a weighted filter, which gives more emphasis for high frequencies

4) I want to combine these frequencies and take the inverse Fourier transform

Can someone give me any guidelines where I should start (which functions etc.) to do this in Matlab? =) I'm kinda learning about Fourier transform and I want to play around with images. I was hoping if someone could give an example of code how this could be done etc.

Please note that I'm new at this stuff and I'm not very familiar with all the terminology yet. I'm an amateur trying to learn about Fourier transform by doing an experiment with it :)

Thank you for any help! =)

P.S. I would appreciate if someone could give me a minimal code snippet showing me what I need to do =)


texts:
I have seen some literature where the covariance is updated first, like $(P_k)^{-1} = (P_k^-)^{-1} + H^T R^{-1} H$, where $P^-$ is the a priori estimate of the state covariance $P$. Then, the updated covariance is used to calculate the Kalman gain $K$.

And I have seen other literature where the Kalman gain is first calculated using the a priori estimate of $P$, then $P$ is updated as $P_k^- - K_k H P_k^-$.

Could anybody provide any unification on this? Are the two methods mathematically equivalent?


texts: Convolution equivalent to matrix multiplication?

texts: Is it possible to write the full convolution between the image and the filter as a matrix multiplication operation? If so, can someone give a simple example of how that works?


texts:
Following this paper (see the Matlab example at the end), I am trying to make a least-squares algorithm in Octave, but for type II (I know about firls). For type I FIRs, it works flawlessly, but not type II. This is my attempt (in Octave):

M = floor(N/2);
oddN = 1 - mod(N, 2);
n = [1:2*M + (~oddN)] - 0.5*(~oddN);
m = [0:M];
if(oddN)
    q = [wp+K*(1-ws), wp*sinc(wp.*n) - K*ws*sinc(ws.*n)];
else
    q = [wp*sinc(wp.*n) - K*ws*sinc(ws.*n)];
end
Q = (toeplitz(q(m+1)) + hankel(q(m+1), q(m+M+1)));
b = wp*sinc(wp.*(m+0.5*(~oddN)))';
a = Q\b;
if(oddN)
    h = [a(M+1:-1:2); 2*a(1); a(2:M+1)];
else
    h = [a(M+1:-1:2); a(1); a(1); a(2:M+1)];
end

Could someone please point out the mistakes?

I found one of the mistakes: q derives from the integral whose second part is simplified as -K*ws*sinc(ws), but the 0.5 offset in sampling does not discard the K*sinc(k+0.5) term. That was just a bad omission from my part. So I added it:

n = [~oddN : 2*M] + (~oddN)/2;
q = wp*sinc(wp.*n) + K*(sinc(n) - ws*sinc(ws.*n));

The results changed, but they are still wrong. Even so, the a(1) terms in the middle seem to need multiplying by a constant. I used the hammer a bit and ~2.25 seems to be a reasonable value. It doesn't change the fact that the results are wrong (or that the hammer is not a solution).


texts: How can I find the delay between two signals?

texts:
I am using NARX in Matlab.
Is there any method to find the delay between input and output signals?
My aim is to decrease memory length by finding the delay.
TIA.


texts: Separate two sinusoids with very close frequencies (5 Hz difference) in a 44.1 Khz audio signal, with good time-domain resolution

texts:
Here is an audio signal which is a mix of (WAV file available here):

a 2000 Hz sinusoid, beginning at 1.00 sec,  without fade-in/fade-out

a 2005 Hz sinusoid, beginning at 1.031 sec, with a slow fade-out at the end

background noise

(In reality it's even more complex: the sinusoids can also vary in amplitude...)
How to separate the signal into two signals (the two sinusoids), but also with a good temporal resolution?
Attempt 1
I zero-padded the signal to the next power of 2 (final length: 524288), did a real-FFT.
The size of the real-FFT vector h is 262145 frequency bins to cover the frequency range [0, 22050hz], so each bin has a width 0.084 Hz. Pretty good news, we can distinguish the two sinusoids with this!

Now we can isolate the two sinusoids with:
h1 = h.copy()   # the real-FFT
h1[:23750] = 0  # zeroing bins outside [23750, 23800]
h1[23800:] = 0
x1 = irfft(h1)  # inverse real-FFT

h2 = h.copy()
h2[:23810] = 0  # zeroing bins outside [23810, 23860]
h2[23860:] = 0
x2 = irfft(h1)  # inverse real-FFT

It works, but the time-domain resolution is very bad (it's normal because the frequency resolution is very high!) : the sine is very poorly localized in time-domain. Instead of a fast attack in the separated sinusoid, we have a slow fade-in...
Of course, zeroing bins outside of the frequency range of interest (bins [23750:23800]) is not optimal, I should have used a (non-rectangular) window in the frequency domain instead:
h1[:23750] = 0
h1[:23800] = 0
h1[23750:23800] *= window

But even with such a windowing I doubt I can avoid a slow time-domain resolution after the separation.
Attempt 2
Use a STFT instead of a global FFT of the signal. This helps to localize, but ... in order to have a good frequency resolution to be able to separate the two sinusoids, we have to take a big FFTSIZE, such as 16384. Then each of the
8193 frequency bins (real-FFT) will have a 2.7 Hz width!
Not enough to distinguish or separate the two sinusoids that have only 5 Hz of difference... So this approach will fail.

I know this is probably an example of the time-frequency trade off / uncertainty principle, but in this precise case, is there something else we can do to improve the separation?
Remark: I already read Separating waves of very close wavelengths but it did not help for this case.


texts: Intuitions on Kumaresan-Tufts algorithm for exponential fit

texts:
I am analyzing a transient signal presumably consisting of superposed exponentials. Such a case is indicated for the Prony analysis, but my data aren't noiseless enough, so I have turned to the Kumaresan-Tufts (KT) algorithm.

After reading the original article (Estimating the Parameters of Exponentially Damped Sinusoids and Pole-Zero Modeling in Noise, 1982) and a bit of googling I made use of the Matlab package Complex Exponential Analysis and more or less things work

My concern is now an intuition on the process - or, better, on its input parameters (cause FFT-like thinking is of course out of question):

What should I expect of increasing or decreasing of model order?
How can I assess the number of modes decinig for signal reconstruction and how can I pick them from the output parameters (dampings, frequencies, complex amplitudes)?
Is there any recommended signal treatment before the KT method is applied? Such as detrending the data for FFT.
Does time-inverting of the signal providing any help? It would mean damped exponentials instead of growing.


texts:
I am doing the project in audio fingerprinting, in that after I convert the signal into mono, I have to apply the hanning window, the size of the window is 370ms, and I have to overlap the signal with a factor of P=28/32.

How can i do that?

What is meant by 28/32?

How can I overlap this signal with this overlapping factor?


texts: I am working with phoneme clustering. I have been using $k$-means clustering have done $k$-means on 2 different data: spectrogram vector and MFCCs. I would like to know how I can evaluate the result of clustering quality? I mean I need to know standard ways to validate the clustering. It would great if someone explains their experience.


texts: Whats the optimal window function to use for analyzing real-time data samples?

texts: Say you wanted to run a X point FFT on the last X audio samples that were played. The problem being, using a normal hann window function would place emphasis on the "middle" of the audio sample. However, the most current audio samples played (the audio samples toward the end) are of the greater importance for algorithms such as real-time beat detection. The question is, which window function best emphasizes the frequencies at the end of the waveform?


texts: I was referring to this link https://cdn.selinc.com/assets/Literature/Publications/Technical%20Papers/6059_HowMicroprocessor_Web.pdf?v=20180606-230156 . I am not very clear about the derivation of the cosine filter and how and where it comes from?


texts: Calculating the PDF of a waveform from its samples

texts:
A while ago I was trying different ways to draw digital waveforms, and one of the things I tried was, instead of the standard silhouette of the amplitude envelope, to display it more like an oscilloscope.  This is what a sine and square wave look like on a scope:

The naïve way to do this is:

Divide up the audio file into one chunk per horizontal pixel in the output image
Calculate the histogram of sample amplitudes for each chunk
Plot the histogram by brightness as a column of pixels

It produces something like this:

This works fine if there are a lot of samples per chunk and the signal's frequency is unrelated to the sampling frequency, but not otherwise.   If the signal frequency is an exact submultiple of the sampling frequency, for instance, the samples will always occur at exactly the same amplitudes in each cycle and the histogram will just be a few points, even though the actual reconstructed signal exists between these points.  This sine pulse should be as smooth as the above left, but it isn't because it's exactly 1 kHz and the samples always occur around the same points:

I tried upsampling to increase the number of points, but it doesn't solve the issue, just helps smooth things out in some cases.

So what I'd really like is a way to calculate the true PDF (probability vs amplitude) of the continuous reconstructed signal from its digital samples (amplitude vs time).  I don't know what algorithm to use for this.  In general, the PDF of a function is the derivative of its inverse function.

PDF of sin(x): $\frac{d}{dx} \arcsin x = \frac{1}{\sqrt{1-x^2}}$

But I don't know how to calculate this for waves where the inverse is a multi-valued function, or how to do it fast.  Break it up into branches and calculate the inverse of each, take the derivatives, and sum them all together?  But that's pretty complicated and there's probably a simpler way.

This "PDF of interpolated data" is also applicable to an attempt I made to do kernel density estimation of a GPS track.  It should have been ring shaped, but because it was only looking at the samples and not considering the interpolated points between the samples, the KDE looked more like a hump than a ring.  If the samples are all we know, then this is the best we can do.  But the samples are not all we know.  We also know that there's a path between the samples.  For GPS, there's no perfect Nyquist reconstruction like there is for bandlimited audio, but the basic idea still applies, with some guesswork in the interpolation function.


texts: Recover signal from its power spectral density

texts: I was wondering whether we can recover a signal $x(t)$ from its power spectral density using signal processing techniques.


texts: How to solve a time delay from an FIR filter?

texts:
I have designed a FIR high-pass filter $\{h_n\}$ with length 257 to filter a signal $\{f_n\}$ with length 512. To convolve signal with filter in time,
I performed an FFT with length 768, padding zeroes at tail of $\{h_n\}$ and $\{f_n\}$, and then I multiplied their frequency responses $\{H_n\}$ and $\{F_n\}$ . To convert it back to time I performed an IFFT producing a filtered signal $\{\tilde{f}_n\}$.
However, $\{\tilde{f}_n\}$ has a time shift comparing to the origin signal, and their amplitudes also have a deviation.
I know there is an 'Overlap-Add Method' to solve this problem, but my calculation resources are limited.
Is there any method to solve such problem using minimal times of FFT/IFFT?


texts: How to calculate a delay between two signals in frequency domain?

texts:
Calculation of a delay between two signals is generally done by finding the maximum in the cross-correlation of the two signals (in time domain):
$$
\tau_{delay} = \arg\max_\limits{t}(f(t) \star g(t))
$$
A cross-correlation can be easily computed in the frequency domain by pointwise multiplication of the Fourier transform of one signal with the complex conjugate of the Fourier transform of the other signal:
$$
f(t) \star g(t) = f^*(-t)*g(t)=\mathcal{F}\{F\}^*\cdot \mathcal{F}\{G\}
$$

Now I wonder how I can compute the cross-correlation peak in the frequency domain, thus without the need of an inverse transform. I can imagine two easy cases:

The case for two sinusoidal signals seems trivial, as the (smallest) delay then agrees with the phase (divided by the angular frequency) of the Fourier component (of this pointwise multiplication) for the respective frequency.
This can be generalized to two identic signals with a fixed delay, as the (unwrapped) phase response equals then the angular frequency multiplied by the delay, for all Fourier components. To avoid ambiguities, let's assume that the delay is smaller than half the period of the highest frequency that is considered.

However, when the two signals are not exactly identic, I suppose some kind of "weighed sum/integration" should be taken of the phase response, where the weighting factor scales with the relative contribution of the respective frequency, i.e., with the magnitude of the respective Fourier component. However, I am not able to find the exact formula rigorously.

Some "internal brainstorming" that might lead to the solution: if we consider one of the two functions ($f$) as a "basis" function, and project the other signal ($g$) on this function, the result of this projection (inner product) should contain the delay we are looking for. I guess this is only possible if the "basis" is first constructed as the analytic representation of $f$, such that phase differences can be taken into account. I also think that this inner product should be conserved after Fourier transform.


texts:
I was wondering why the energy of a signal is defined as

$$E=\int_{-\infty}^{\infty}|x(t)|^2dt$$

I mean, I know that it is a definition, but is there a physical reason for it being defined like that?


texts: Differential Pulse Code Modulation (DPCM)

texts:
In the following figure we see an encoder and a decoder of a Delta Modulation System, which is a simple DPCM System with a 1-bit quantizer:

We assume that signal $x(n)$ is known, that $\hat{x}(-1)=0$ and that there is not noise in the telecommunication channel.

I want to define the predicted signal $\tilde{x}(n)$, the difference signal $d(n)$, the quantized difference signal $\hat{d}(n)$, the encoded signal $c(n)$ and the reconstructed signal $\hat{x}(n)$

I have written the following equations:
$$d(n) = x(n) - \tilde{x}(n) \; \; \; \; (1) $$
$$\tilde{x}(n) = 0,8 \cdot \hat{x}(n-1) \; \; (2)$$
$$\hat{x}(n) = \tilde{x}(n) + \hat{d}(n) \; \; \;\;(3)$$
In addition, duo to the feedback, the error $e(n)$ for the reconstruction (between the initial signal on the transmitter and the reconstructed signal on the receiver) is equal to the error of the quantizer. So, I also get:
$$\hat{d}(n) = d(n) + e(n) \; \; \; \; (4) $$
$$\hat{x}(n) = x(n) + e(n) \; \; \; \; (5) $$

After solving the system of the equations, it seems to me that the number of the uknowns exceeds the number of the linear independent equations. If it so, I need one more equation since not all the equations above are linear independent.

So my question is what am I missing? Thanks in advance!


texts: Overview of regression methods between the data

texts:
I am trying to find different approaches (classic and modern) in order to find a regression function between data (in time domain) with high-frequency noise.
Consider the result of this code (in Matlab):

fs = 10000;
t = -1:1/fs:1;
target = tripuls(t,50e-2);

noise = target + randn(1, length(target));

plot(noise,'o')
hold on
plot(target,'r')

I would like to find methods that, starting from noise find a function similar to target (for convenience, I attach here an image).

A classic solution can be to implement a low-pass filter. Did you know other (and maybe better) way to reach the expected result?

Thanks


texts: Why are there comb-like hills in a periodogram?

texts:
I am playing with the periodogram of MATLAB. I created a simple script to observe how it behaves:

rng(1);  %# initialize the random number generator

Fs = 1000;  %# Sampling frequency
duration = 0.1; %# seconds

A = 1; %# Sinusoid amplitude
f = 150; %# Sinusoid frequency
eps = 0.01;

t = 0:1/Fs:duration;
x = A * sin(2*pi*f*t) + eps * randn(size(t));

periodogram(x,[],1024,Fs);

I have no problem with the code and can write my own periodogram function using the algorithms given in the documentation but I wonder the theoretical reason behind the comb-like hills which are not 150 Hz. What do I get those instead of getting a single spike over 150 Hz? Is there anything special in the distances of the peaks of these hills?


texts: Updating FFT algorithm accordingly when upsampling/resampling

texts:
I have some experimental data from a shaker (vibration exciter) running at sine $10\textrm{ Hz}$ for test purposes (recorded using an ADXL345 accelerometer with $f_s=512\textrm{ Hz}$. I would like to upsample these data for varies ratios, e.g. a factor of $\times 10$.
How do I update the FFT in MATLAB accordingly? As of now, upsampling appears to shift results in the $x$-direction as shown in the figure below (black is the original signal, blue is the resampled signal with cutoff at $1/2\times f_s$):

Peak at approximately 80 Hz (black) appears at approximately 40 Hz (blue).

But I would like to increase the resolution within the existing data without shifting results in the $x$-direction?

MWE:

close all; clear all; clc

fs=512; % sampling frequency

load data.mat

t=X{1}; % time
y=X{4}; % accelelation

y=y./1024; % calibration of data

%% resampling

%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%
ups=10; % upsampling rate
%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%

dns=1; % downsampling rate
fu = fs*ups; % upsampling frequency
% tu = t(1):1/fu:t(end)-1/fu;
ytu=resample(y,ups,dns);

%% FFT and plot

figure;

x=y;
[N,m] = size(x);

freq = 0:fs/length(x):fs/2; %frequency array for FFT
xdft = fft(x); %Compute FFT
xdft = 1/length(x).*xdft; %Normalize
xdft(2:end-1) = 2*xdft(2:end-1);
semilogy(freq,abs(xdft(1:floor(N/2)+1)),'k-')
hold on

x=ytu;
[N,m] = size(x);

freq = 0:fu/length(x):fu/2; %frequency array for FFT
xdft = fft(x); %Compute FFT
xdft = 1/length(x).*xdft; %Normalize
xdft(2:end-1) = 2*xdft(2:end-1);
semilogy(freq,abs(xdft(1:floor(N/2)+1)),'b.')

xlabel('Frequency (Hz)');
ylabel('Accel (g)');
grid on;

Here is the link to data.mat.


texts: Is spectral leakage due to windowing 'different' for the DTFT and DFT?

texts:
I'm currently trying to improve my foundational understanding of spectral leakage due to windowing. The more I have been reading about it, the more I'm starting to think that 'windowing spectral leakage' can mean two completely different things depending on whether you're talking about the DTFT or DFT. I'm hoping to get some confirmation/clarification on what I discuss.

With regards to the DTFT, this transform is essentially the Fourier transform equivalent for a discrete signal. This transform assumes that the signal is periodic with that period being infinity. With regards to windowing, spectral leakage is induced due to the fact that you are multiplying a discrete signal $x[n]$, by some window function $w[n]$. The corresponding convolution in the frequency domain causes the spectrum of $X[f]$ to be smeared or 'leaked'.
Now moving onto the DFT, the idea of spectral leakage induced by windowing seems to be a completely different concept. Unlike the DTFT, the DFT assumes the signal is periodic (with that period being less than infinity). In this context, it is often stated that spectral leakage will occur if the collected signal does not contain an integer multiple number of points to the digital frequencies contained (i.e. if you have a $3\textrm{ Hz}$ signal with sample rate $10$ samples/sec, and collect $18$ samples, you will have spectral leakage). Thus if you recorded a sinusoid at interval where it's not able to 'finish' a complete cycle, the next period will restart and there will be a discontinuity in the signal. While this makes sense that you would be inducing 'false' frequencies or leaking frequencies, it seems to be a completely different idea from the convolving concept of the DTFT.

My primary question:

Is 'window spectral leakage' completely dependent on context (DTFT vs DFT)?
Furthermore, is it safe to say that the spectral leakage induced by each case separate from one another? As in, spectral leakage from the DTFT is induced by frequency convolution, and spectral leakage from the DFT induced by the assumption of periodicity? Or is there something that ties these two things together?


texts: Frequency modulation implementation without trig functions

texts:
https://dsp.stackexchange.com/a/1087 This answer suggests generating a sine wave using complex phasor multiplication. At each step the phasor is multiplied by another complex number (omega) to rotate it by the desired step. Then the imaginary part of the phasor is a sine wave. Can this be extended further to implement frequency modulation without having to use a trig function as part of updating the angular frequency? I've played with it a bit, but haven't figured out a clever way to remove the trig function when updating the frequency.

My initial idea: the phasor would be updated like so:

phasor = phasor * (omega * modulator)

The modulator would sweep back and forth between exp(j*freq_mod_max) and exp(j*freq_mod_min). I haven't figured out how the modulator would be implemented.

Side note: why would you want to do this? I did some quick and dirty performance tests and using this phasor multiplication method is quicker for generating a sine wave (about 3-4x quicker). Though calling the sine function directly is still of course quick enough so this is a somewhat academic curiosity at this point.


texts: How are filters for polyphase decimation like?

texts: I know in naive decimation the signal should filtered by a low-pass filter with cutoff frequency equals to Nyquist rate.I don't know how the filter banks like for polyphase decimation. Are them still low-pass filters? and what are their cutoff frequency.


texts: Can convolution of one signal with different signals give the same answer?

texts:
Let us consider $x_1(t)$, $x_2(t)$, $x_3(t)$, all the same within some some duration 0 to $T$ but all different outside this interval.  Now let us multiply each of these signals with $w(t)$, a window function - nonzero from 0 to $T$ but 0 outside this interval.  So the multiplication of this $w(t)$ with each of $x_i(t)$ will give the same signal.  This should be same as Inverse Fourier Transform of $W(\omega)$ convolved with $X_i(\omega)$.  Does this mean convolution of the same signal with different signals can give the same result?  Any comments?  Is there a pitfall in my interpretation?

Also I see mathematically equations neatly entered in questions and answers on this site?  Where should I start to learn on how to write the equations?


texts: What are some deep results in sampling theory?

texts:
I have taken a basic course or two in signal processing and I've seen the Nyquist sampling theorem - an interesting and surprising result giving sufficient conditions for reconstruction of a reasonably regular signal in terms of spectral support.

I was told there are many much deeper results in this field than this canonical theorem. Unfortunately, I find the literature rather inaccessible. Hence, I'm asking here for examples of deep results in sampling theory.

I am most interested in results which rely on other branches of mathematics, namely topology and algebraic geometry.

I understand the question is vague but also think that it is clear enough to give answers to. I don't mind making this a community wiki since it has no single correct answer.


texts: Do we really need two cameras for 3D vision?

texts:
I heard that we can see the world 3D or we can detect distance because we have 2 eyes. but I think we can do it with only one eye! (test it!)

So why we can do it only with one eye? is it possible to have stereo vision for a robot only with one camera?


texts: Figuring out how an RF pulse affects the MRI spins

texts:
So, I understand how the slice selection gradients work in MRI. So, the frequency offset introduced by the slice selection gradient at a location $z$ relative to the MRI isocenter is given by:

$$
\Delta f = \gamma z G_z
$$

Putting some numbers to it, for a gradient of 5 mT/m and a slice thickness of 3 mm, we introduce an offset of $639.7$ Hz per slice.

Now, say I have an RF pulse which has a bandwidth of 2000 Hz and a center frequency of 10 kHz. Given the location of the spin along the slice direction, I want to know whether this spin will be excited by the RF pulse or not.

Is it possible to figure that out from this information?

I have a feeling it should be possible from the relation above but I am not sure how the bandwidth parameter will come into play here.

Is it a simple matter of checking that if the resonant frequency experienced by the spin falls in the RF center frequency $\pm$ bandwidth / 2.0?


texts: Implementation of wideband beamformer for planar array

texts:
I would like to obtain a good reference (or references) on the implementation of a wideband beamformer for a small planar (rectangular) array comprised of 4 rows by 6 columns, for 24 elements in total.  This is a passive acoustic array with air-coupled microphones; the array does not transmit.  Each element in the array has the same spacing distance.  Given discrete input signals from each element, I would like to compute the beamformer/beamsteering output for a beam at a particular azimuth and elevation.  I need to be able to compute multiple beams to electronically steer the array. Is it possible to compute multiple beams all at the same time?

I suppose that I might need to use Fourier-domain techniques involving the FFT.  Are these techniques justified for my small array?  I've found an IEEE conference paper (http://www.curtistech.co.uk/papers/wideband.pdf or http://dx.doi.org/10.1049/ic:19980121) that may deal with this technique, but how would I apply this to a 2D planar rectangular array?  I am interested in the mathematics operating on discrete signals rather than  descriptions of the process.

Since this is an acoustic array, the elements are evenly spaced on a 4 rows by 6 columns grid.  The spacing between the elements is less than $\lambda/2$, where $\lambda$ is the minimum wavelength of the sound wave.  The wavelength range is between 3.31 cm and 331 cm.  This is a receiving array.  The target area should be in the far field of the array, but I am also wondering if it is possible to resolve signals in the near-field as well.


texts: Converting floating point IIR to fixed point arithmetic, given its difference equation

texts:
Is there some general way to perform the task given in the question? I have a filter given as follows:

Tdouble MoogVCF::run(double input, double fc, double res)
{
  double f = fc * 1.16;
  double fb = res * (1.0 - 0.15 * f * f);
  input -= out4 * fb;
  input *= 0.35013 * (f*f)*(f*f);
  out1 = input + 0.3 * in1 + (1 - f) * out1; // Pole 1
  in1  = input;
  out2 = out1 + 0.3 * in2 + (1 - f) * out2;  // Pole 2
  in2  = out1;
  out3 = out2 + 0.3 * in3 + (1 - f) * out3;  // Pole 3
  in3  = out2;
  out4 = out3 + 0.3 * in4 + (1 - f) * out4;  // Pole 4
  in4  = out3;
  return out4;
}

It's a 4th order IIR low pass filter which works in a feedback loop that generates a nice resonance. It sounds great, but my only problem is that my current CPU is too slow to handle it. It has no FPU and floating point operations take approx. 10x more cycles. I was thinking about converting it into fixed point arithmetic but after long time spent over a piece of paper and calculator I just gave up. Nothing seems to work. Do you know of any method with which I could convert the coefficients and add some bitshifts to make it a fixed point arithmetic filter?

I have 16bit DAC so the samples I process are 16 bit signed but the CPU is 32 bit so the filter can operate on 32bit integers.


texts: Matlab: Different 1-D FFT results with and without zero-padded

texts:
Recently I am implementing FFT with C code. So I use Matlab for algorithm verification. As Radix-2 FFT was chosen, I've learned that zero-padding technique could help to reach 2^N limitation without messing up the result. However, hesitation has shown.

Let's say a set of acceleration data with the size of 150 (x) comes across. If FFT function is called directly like:

X = fft(x);
abs_X = abs(X);
figure
plot(abs_X(2:75)); % discard first data (DC)

The result looks reasonable:

With this result I can easily indicate the highest amplitude and get the main frequency with frequency step. However, after zero-padding the data to the size of 256 and then apply fft, a confusing result is obtained.

x = [x; zeros(106,1)]; % pad 0 to 256
X = fft(x);
abs_X = abs(X);
figure
plot(abs_X(2:128)); % discard first data (DC)

Now I cannot identify the main frequency since the decreasing pattern is shown.

From here I have learned that the fft function in Matlab will automatically adjust the data size with zero-padding in order to perform FFT. So theoretically two operations demonstrated above should perform same results. I know that zero-padding will increase the visual resolution (not data resolution). But the difference between the beginning of the figures, where main frequency can't be determined, is confusing me.

Can anyone help to clarify the situation? Or do I misunderstand something important about the statement "Zero padding is going to perform a same result."? Any comment is appreciated. Thank you for browsing my question!


texts: How to upscale an image without affecting pixel aspect ratio?

texts:
For example, if we upscale an image from 100 x 100 to 400 x 300 then we see the pixel stretch in the image which I want to eliminate.

The easiest way I could think is to crop the image to maintain the pixel aspect ratio i.e. upscaling to 300 x 300 or 400 x 400

Edit (Clarification): I want upscale an image without causing pixels to stretch, upsampling algorithms like nearest neighbor, bilinear interpolation, etc cause pixel stretching, especially when there is a change in aspect ratio e.g. from 1:1 to 16:9


texts: Calculating power and energy of a signal in matlab

texts:
I am a little new to matlab and would like to get some help on how to write a program or a function that can calculate the energy/power of a signal depending on whether the signal is an energy or a power signal.
Here is the program the I tried:

function [ ] = power_or_energy( s,t)
%UNTITLED Summary of this function goes here
%   s is the input signal
%   t is the independent variable

s2=s.^2;
energy=int(s2,t,-inf,inf);

if energy~=inf && energy~=-inf
    disp ('the signal is an energy signal with energy equal to ' );
    disp(energy);
end

power=limit((int(s2,t,-t/2,t/2))/t,t,inf);

if power~=inf && power~=-inf
    disp('the signal is a power signal with power equal to ');
    disp(power);
end

But the problem here is that when i want to check it for a function like

t=0:0.01:10;
s=sin(t);
power_or_energy(s,t);

then i get the following error:
"Undefined function 'int' for input arguments of
type 'double'."

So how to i actually calculate the energy for such signals?

Also the above function works fine with the below input:

syms t;
s=sin(t);
power_or_energy(s,t);

Please help me in writing a function that can be used to find the energy of any input signal.


texts: How to detect corners in a binary images with OpenGL?

texts:
I have binary 160x120 images such as:

I would like to detect corners of those white blobs. They are previously closed by mathematical morphology so there shouldn't have any inner corners. In this specific case, I would want 16 corners, like:

My first attempt was using some OpenCV functions like goodFeaturesToTrack or FAST but they are particularly slow (plus FAST is very unstable). My idea would be to do such a computation on the GPU, as my source image comes from it. I looked for ideas on the web on how to write such shaders (I'm using OpenGL ES 2.0), but found nothing concrete. Any idea how I could start such an algorithm?


texts: Is it possible from 1 order cascading low pass filters to make a butterworth filter

texts:
Today I have heard, that is "impossible" to make low pass butterworth filter from cascading 1 order low pass filters(RC).

I am trying to discover is that true, but I can find little to none information on their difference.

Can you please explain to me, is that true, and if it is why. Thanks?


texts: Is there a closed form expression for main-lobe width increase given a window?

texts:
We know that when we window a signal, we increase the main-lobe width. Let 'main-lobe-width' here be the null-to-null bandwidth of the main lobe. Let us further more say that the main-lobe width of a square window is '1'.

What I would like to know, is given a particular symmetric window, is there a closed form solution for how much the main-lobe will increase by, relative to that of a square window?

I can look up this percentage increase for various windows just fine by the way. I am asking if there is a closed form solution for what the percentage increase is, given an arbitrary symmetric window.


texts: Source separation on multiple observations of the same signal

texts:
I possess an ensemble of signal observations:

$x_i[n]=s[n]*g_i[n]$, $i=1,2,....,N$ where $N$ is a very large number compared to individual signal lengths (signal lengths are identical). Here, $s[n]$ is a time localized signal (e.g. a positive deflection) while $g_i[n]'s$ are arbitrary signals.

I have two intentions:

Recovering the original signal, $s[n]$.
Recovering 'some information' for each $g_i[n]$. I am specifically interested in how spread these signals are. So, rather than exact $g_i[n]'s$, some statistics about their temporal centrality would be more than enough.

I am eager to hear about possible techniques for achieving these. I should note that, since these will be a part of  an academical work, I prefer neat and structured techniques. Lastly, I prefer computationally inexpensive ways since I will perform this operations multiple times.


texts:
Consider a circular antenna array (UCA) with $N$ antenna elements as part of the CW bistatic radar system.
For each element we receive a signal of length $L$: $x_n[k]$, $k=1...L$, $n = 1...N$.
Оne of the antenna elements is directed towards the transmitter.
Thus, we can find cross ambiguity functions (CAF) as cross-correlation functions between the reference channel and the remaining $N-1$ channels.
Further, in each CAF, we detect targets and determine the ranges to them.
My question is how to define directions on a target.
Standard explanations include beamforming. The way I understood this method: the samples from each of the channels, except for the reference for same time point, are used to form the directional diagram, i.e. from the $N \times L$ of the signal samples, $L$ diagrams are formed. Each directional diagram includes lobes directed to all targets as well as side lobes.
How to find out which direction corresponds to each target with its own range and how to exclude side lobes if their amplitude can be higher than the amplitude of the lobe directed at a target with a weak reflected signal?


texts: Estimate Delay of a Known Signal Delayed by Sub Sample Resolution

texts:
Given a known signal $ x \left( t \right) $ and its delayed version $ y \left(t, \tau \right) = x \left( t - \tau \right) $.
Both are sampled by Sampling Frequency $ {F}_{s} $ to generate the signals $ x \left[ n \right] $ and $ y \left[ n \right] $.

Assume they are sampled on a time support much larger than their non zero time (Namely, the whole part of the signal where it is not vanished is sampled).

Given those 2 sampled signals, using Cross Correlation (Matched Filter) the Delay Parameter $ \tau $ can be estimated.

Yet, its estimation is limited to the grid defined by the sampling Frequency which means if $ \tau < {T}_{s} = \frac{1}{{F}_{s}} $ the estimation isn't accurate (Actually, the part which is the modulo of the time interval can't be estimated).

Is there an effective method to estimate the delay in Sub Sampling resolution?

Thank You.


texts: What is the physical interpretation of the absolute value of a fourier transformed signal, $\left| F(t)\right|$?

texts:
If one has some oscillating voltage signal, for example:
$$V(t) = V_{max}\cos(2 \pi \nu_{0}t) e^{-\gamma t}$$
and you take the Fourier transform of this in the usual way to get:
$$\hat{V}(\nu) = V_{max}\int^{+\infty}_{0} \cos(2 \pi \nu_{0}t) e^{-\gamma t} \ dt$$
Now I often see in literature the absolute value squared taken of the result, so $\left|\hat{V}(\nu) \right|^{2}$ the physical interpretation of which it is the energy per unit frequency.

My question is, if $\left|\hat{V}(\nu) \right|^{2}$  is energy per unit frequency, then what is the interpretation of $\left|\hat{V}(\nu) \right|$?


texts: Adding White Gaussian noise to a voice signal

texts:
I'm trying to add White Gaussian Noise to an audio file. However; the energy of the noise should be 1/10th of that of the signal. My first attempt is as following:

[y,Fs] = audioread('drum.wav');
%sound(y,Fs);
noisy_sig = addnoise(y,10);
sound(noisy_sig,Fs);

function noisy_sig = addnoise( sig , SNRdb )
   sig_power = norm(sig,2) / length(sig);
   % noise power is equal to sigma^2
   sigma2 = sig_power / 10^(SNRdb/10) ;
   noisy_sig = sig + sqrt(sigma2)*randn(size(sig));
end

When I listen the resulting wav, there is no difference between the original and noisy ones. Am I doing something wrong? Any help would be appreciated.


texts: How do I design a speed-efficient high-pass filter to replace slow polynomial filter?

texts:
I have the recorded membrane potential of a neuron, and I want to get a more even baseline by removing certain slow changes from the signal (see figure). The best way to achieve this thus far has been to make a polynomial interpolation of the data, and subtract it from the raw data signal. I divide the data into windows of $20000$ samples, and then apply the filter to each window separately.

However, beside not being a very neat solution, it is by far too slow to be incorporated in the analysis program I am constructing (where the raw data signal is of order of magnitude $10^8$ samples). I have experimented with MATLAB's designfilt, but so far all the filters I have created with this function distort the small variations with $y$ amplitude around $1$ in the figure, and they have to be kept relatively intact.

I know this is not a question with a single simple answer, but any ideas about how to speed up the execution of my polynomial subtraction filter are most welcome.

I have also included MATLAB code below to generate the displayed figure. If there is a way to upload files I could upload some of the raw data values too.

x = load_x(); % Too large to be uploaded as a text file (even a 20000   snippet will clog the post

polynomial_nbr = 18;
width_indx = 20000;

y = polynomial_filt(x, polynomial_nbr, width_indx);

clf
t = 1:length(x);
plot(t, x, 'b', t, y, 'r')
legend('Pre-filtered', 'Post-filtered');

function y = polynomial_filt(x, polynomial_nbr, width_indx)

y = nan(size(x));

for i=1:width_indx:length(x)

  indx = (i-1 + (1:width_indx))';
  indx(indx>length(x)) = [];

  [p, ~, mu] = polyfit(indx, x(indx), polynomial_nbr);

  y(indx) = x(indx) - polyval(p, indx, [], mu);

end

end


texts: How to process the .wav file to get meaningful information like the number of inhales/exhales

texts:
The plot of the data from the wav file.

Link to the sound file - https://drive.google.com/open?id=0BxAPlY66uX7qWWtDZUw3WkpjSGs&authuser=0

Also the FFT plot is like this

I am new to the field of digital signal processing.

I have a .wav file in which I have recorded the breathing sound in a closed environment. I am trying to figure out what kind of technique I can apply to process the .wav file to get meaningful information like the number of inhales/exhales.

I have read the file in python using Scipy library and plotted the data. Also I have applied Fast Fourier Transformation to plot the frequency amplitudes. At this point I am feeling lost, how could I get meaningful information like the number of exhales/inhales?


texts: What is wrong with RMS and Maximum delay spread Matlab code?

texts:
pls, help check this matlab code. The code is to determine the RMS and Maximum delay spread. After running the code below, I was not getting desirable results, the Max. delay spread is giving the same value in all the data except where i have more than two spikes. I also got some negative max. delay spread in some cases, I guess I did do one or two things right. Pls, can someone help to check the work-ability of this Matlab code below? Thanks in anticipation.
this is the matlab code:
clc

clear

format short

ts= 1e-06;  % ts=1

T = readtable('Rx8.csv');

kk=table2array(T);

z=length(kk);

allData=zeros(z,2);

num_5=0;

sum_p_i=0;

Mean_delay=0;

Mean_t2bar=0;

for ii=2:z
% if(ii)

s(ii,1)=kk(ii,1)+1i*kk(ii,2);

pdp=abs(s(ii,1))^2; % Values for PDP

allData(ii,1)=pdp; % Store PDP in matrix
%

tau=(ii-1)*ts; % delay_n

allData(ii,2) = tau; % Store delay_n in matrix
%
num_5=num_5+(pdp*tau); % Equation (2.18)

sum_p_i=sum_p_i+pdp
%
t2bar=(pdp*tau^2)/pdp; %nth_delay_square_bar

end

tbar=sqrt(num_5/ sum_p_i);

for jj=2:z

tau_rms= sqrt(( allData(jj,1)* (allData(jj,2)-tbar)^2)/sum_p_i);
end

[aaa, bbb]=max(allData(:,2))

[aaa, bbb]=max(allData(:,1))

[B,I] = sort(allData(:,1),1,'descend')

delay1=allData(I(2),2)

delay2=allData(I(1),2)

% delay2=allData(I(3),2)

maxdelay=delay2-delay1

tbar

tau_rms

plot(allData(:,2), allData(:,1))

% RMS_delay=sqrt(Mean_t2bar-Mean_delay^2)

%
% t_pdp=1:length(pdp)

% plot(s(:,2),pdp)
%
% [a, b]=max(pdp)

% [a2, b2]=min(pdp)
% plot

% Amp1=sum(abs(s.^2))

% Amp2= rms(s)^2

% Amp3=norm(s)


texts: multiplication in frequency domain equals convolution in time domain mismatch

texts:
We know that a convolution in the time domain equals a multiplication in the frequency domain.
  As per the article below,

https://www.dspguide.com/ch9/3.htm

In order to multiply one frequency signal by another, (in polar form) the magnitude components are multiplied by one another and the phase components are added.
 To prove this, i created two sinusoidal signals in MATLAB,

Freq1 = 1000;
Freq3 = 3000;
Fs = 16000;
T = 1/Fs;
Nos = (0:128-1)*T;
Amp = 1.0;
Signal1 = Amp*sin(2*pi*Freq1*Nos);
Signal3 = Amp*sin(2*pi*Freq3*Nos);

In the time domain, i convolved these two signals,

Signal4 = conv(Signal1,Signal3);

Taking a 32 point FFT of the two signals & and a 32 point FFT of the convolved signal,

NFFT = 32;
freqdata1 = fft(Signal1,NFFT);
freqdata2 = fft(Signal3,NFFT);
freqdata4 = fft(Signal4,NFFT);

in the frequency domain, i multipled the magnitude components of the two individual signals and added the phase component of the two signals. I compared this Magnitude and phase value with the Convolved signal's phase and magnitude value. I expected the values, [Newmag' NewPhase'] & [Mag3' Phase3']  to be similar since the a convolution in time domain equals a multiplication in the frequency domain. But they are not. What am i missing here ? What have i done wrong ?

for ii = 2:((length(freqdata1)/2)+1)

   sig1_cc = real(freqdata1(1,ii));
   sig1_dd = imag(freqdata1(1,ii));
   Mag1(ii-1) = sqrt((sig1_cc^2)+(sig1_dd^2));
   Phase1(ii-1) = atan(sig1_dd/sig1_cc);

   sig2_cc = real(freqdata2(1,ii));
   sig2_dd = imag(freqdata2(1,ii));
   Mag2(ii-1) = sqrt((sig2_cc^2)+(sig2_dd^2));
   Phase2(ii-1) = atan(sig2_dd/sig2_cc);

   sig3_cc = real(freqdata4(1,ii));
   sig3_dd = imag(freqdata4(1,ii));
   Mag3(ii-1) = sqrt((sig3_cc^2)+(sig3_dd^2));
   Phase3(ii-1) = atan(sig3_dd/sig3_cc);

   Newmag(ii-1) = Mag1(ii-1)*Mag2(ii-1);
   NewPhase(ii-1) = Phase1(ii-1) + Phase2(ii-1);

end

[Newmag' NewPhase']
[Mag3' Phase3']


texts: How can I detect the presence of echo in an audio?

texts:
The presence of echo during a phone call can be extremely irritating. Echo consists typically of the speakers sound being picked up by the microphone with a certain delay (10-200ms). One of the causes of this are low quality headsets (no noise cancellation) or a bad setups (the microphone & speaker are placed too close together).

My question is following: Is it possible to detect echo in an audio recording of only one call party? (This is a recording of one side of the call but due to the echo, you can hear both sides; the second party with a lower amplitude obviously) I do not need to cancel the echo, I just want to detect its presence to a certain degree and accordingly recommend a change in the setup.

I appreciate any suggestions and hints, thank you in advance :)


texts: DC and Nyquist values are higher than all other bins

texts:
I am using a mean centered signal and applying fft to it, however the dc and nyquist values are a lot higher than all other values. Wondering if I am doing anything wrong or need to add anything else. The first image is the signal and the other is the fft output.

load Dog_1_interictal_segment_0001.mat

%interictal_segement is a struct with a 16 x 10000 matrix named data
data=interictal_segment_1.data;

%take first two seconds of the first eeg channel
data=data(1,1:10000);
freq=5000;
N=freq*2;
hz=linspace(0,freq/2,N);

averaged_mean=data-mean(data);
channel_frequency=abs(fft(averaged_mean)/N).^2;

figure(1)
plot(averaged_mean)

figure(2)
plot(hz,channel_frequency,'ko-','markerfacecolor','m','linew',1)
xlabel('Frequency (norm.)'), ylabel('Power');


texts: Since DMRS is also being used as a reference signal in PUCCH and PUSCH channel, so why an extra SRS signalling is required.


texts:
I'm trying to segment a specific part in a gray-scale image and I want to achieve this by applying the OSTU algorithm, but after using this function:

im = imread('ready.jpg');
oim =otsu(im,10);
figure,imshow(oim, [1 10]);

I could show the output image, but I do not know how can I save the result to be able to use it later as a mask.


texts:
I am trying to implement the BPSK-63 modulator using liuid dsp library.
Here is my code:
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <getopt.h>
#include <math.h>
#include <complex.h>
#include <liquid/liquid.h>

#define SAMPLERATE 11025

int main(){
    modem mod   = modem_create(LIQUID_MODEM_BPSK);

    nco_crcf q = nco_crcf_create(LIQUID_NCO);
    nco_crcf_set_phase(q,0);
    nco_crcf_set_frequency(q, 2*M_PI*(2000.0/(float)SAMPLERATE));

    const char *message = "1010111100101111011001010100101110100"; //RYry

    std::complex<float> x=_Complex_I*0;
    unsigned int symbol=0;
    unsigned int outs=0;
    int ctr=0;
    float ous[1];
    int offset=0;

    while(1){
        if(ctr++ == (SAMPLERATE/63)) {
            ctr=0;
            char sym = message[offset]-'0';
            offset++;
            if(offset >= strlen(message))
                offset = 0;
            symbol = sym & 0x01;

        }

        modem_modulate(mod, symbol, &x);
        nco_crcf_mix_up(q, x, &x);
        nco_crcf_step(q);

        ous[0] = x.real();
        fwrite((void *)ous,1,sizeof(float),stdout); // | aplay -f FLOAT_LE -r 11025 -c 1
    }
}

It looks like it works, however it sounds a bit different (has some clicks) to working modulator from FLDIGI program.
However when I feed the audio into the fldigi, I can see a carrier in the right spot, however it demodulates a series of "ot e" like:
ot e
ot e
ot e

instead of RYryRYryRYryRYry.
So, it sees five characters instead of four and detects a space and a newline when there is none.
Also, i see some harmonics, but i don't think this is a problem, since fldigi provides an adjustable filter on it's waterfall.
What i'm missing?


texts: Neural networks in system identification - What type of activation functions?

texts:
I made a free software for all operative systems, even Android. It's called Deeplearning2C. It can train a neural network and generate C code and MATLAB-code. C-code for embedded systems and MATLAB-code for simulation.

I tried a example when I trained a neural network with linear data with constraints = nonlinear model in other words. It works very well. Here is the result when I use identity as activation functions.

I identification the model from the features (inputs)

$$-y(t-1), -y(t-2), ... , -y(t-n), u(t), u(t-1), ... , u(t-m)$$

And the label (output) $y(t)$. This is the regular transfer function or ordinary differential equation identification method.

But the issue I have is that I have a lack of experience of system identification when it comes to neural networks. I'm used to state space models (Subspace identification) and regular parameter estimation (Recursive Least Square). I wrote a library for that too.

My questions for you are that:

What type of activation functions should I use?
How many layers is "good" enough?
How many neurons should I have?
How do I know...how to know what to select? I'm seeking practical experience for identification with neural networks.

Here is my software Deeplearning2C. It's using Deeplearning4J as the core, but the user interface is made by me.

Also here is the MATLAB/GNU Octave library for linear system identification with recursive least square and subspace identification. I wrote a C-library too for recursive least square identification as well.


texts:
I'm working on a frequency tracking problem with noise, where the amplitude of the noise is orders of magnitude higher than the amplitude of my signal (~1000x). Some details:

The signal is roughly sinusoidal, less than 1 Hz with non-stationary frequency and amplitude
The noise is caused by motion, and is (in many cases) not periodic
I have a 3-channel accelerometer reference for the noise. My target signal has some correlation with motion, so adaptive filtering may have issues

Empirical Mode Decomposition and even simple windowed FFT work well in the no-motion case, but fall apart as soon as motion kicks in - the broadband effects of the noise mask the frequency-domain peak of the signal. I tried RLS adaptive filtering, but I believe the aperiodicity of the noise prevents the filter from converging to anything useful. Would something like Extended Kalman Filtering be able to handle this case? What other options are there when your noise overpowers your signal?

Update

Phase-locked Loops:

According to this article, PLLs...

  ... can demodulate FM signals with very high accuracy and reliability, or it can detect signals buried in noise, but it can't do both in a single configuration, because the two tasks require very different setups and assumptions. In its role as an FM detector, a PLL doesn't reject noise very efficiently, and as a weak-signal detector, it can't decode FM very efficiently. The reason should be obvious — to detect FM modulation of a given bandwidth, the PLL's feedback loop low-pass filter (Figure 2, green) must be opened up enough to allow the modulation's bandwidth to pass unimpeded, but this causes the PLL to become more susceptible to noise.

Kalman Filtering:

I've been testing out an Unscented Kalman Filter (UKF) for the problem. I found examples of people using Extended Kalman filters, which work ok for the simplest problem of tracking a sine wave because it's easy to derive the Jacobian of a sinusoid. I opted to try the UKF because of the ease of incorporating further non-linear constraints / future problems where the Jacobian isn't easily accessible.

I set up the filter using a state variable of $\begin{bmatrix} A & \omega & \phi & b\end{bmatrix}^T$, where the measured signal associated with the state is $z = A \sin (\omega - \phi) + b$.

The process model assumes that $A$, $\omega$ and $b$ remain constant, and the phase advances according to the current frequency: $\phi[t+1] = \phi[t] + \omega[t] \times dt$ ($\phi$ is constrained to $[-\pi, \pi)$).

The problem I find with this model is that solutions are not unique. For instance, in the case pictured below, the filter has converged by sending the frequency to 0, and manipulating the phase variable to influence the signal:

In each subfigure, the black line represents the estimate from the filter. The top subfigure is the fit of the model against noisy measurements (blue dots). The remaining subfigures show each variable of the state vector - colored dots are the true value that we would like the filter to converge to.

Even if I could restrict the argument of the sin to one process variable, it would always be possible to set the frequency to 0 and use the amplitude or baseline shift to fit the signal.

My sense is there are clever ways of designing the process model to constrain the behavior of the filter. Can anyone suggest a better way to set up the problem?

For reference, I have attached a python 3 notebook with the code I am using.