How is signal to noise ratio actually measured by receiver equipment? This sounds like quite a basic question but it surprised me, how is SNR actually measured?
You have the incoming signal:

It seems like the SNR is just the visual comparison of the peak signal 'amplitude' to the noise 'amplitude'. In this case the signal 'amplitude' clearly includes noise as well, because it is summed together i.e. it is the sum of the noise 'amplitude' at that frequency and the signal 'amplitude' at that frequency because 'amplitude' in this case referring to volt-seconds under the curve of a complex sinusoid at that frequency times the signal (which contains the noise sinusoid and signal sinusoid at that frequency, which creates a greater amplitude sinusoid at that frequency and therefore there is more area (volt-seconds) under the curve when that sinusoid is multiplied by a complex sinusoid of the same frequency.
I'm assuming the received signal power would be calculated by squaring the $V_{RMS}$ of the samples taken and I assume noise is also sampled at a separate occasion to get the $N_{RMS}$ (or so called noise variance, but real noise doesn't have a 0 mean i.e. no 0 frequency component), or it uses a frequency spectrum analyser or something to get the voltseconds at unused frequencies and then integrate it with respect to frequency i.e. $s^{-1}$ to get the average voltage and then square it to get the average power, the problem with that being that you can't measure the noise on the frequencies that are used, so I would think that it is acquired separately, like in the wifi layer 1 preamble and training sequences.
But the signal comes in as the signal+noise, so surely it's actually (signal (${V_{RMS}}^2$ + noise ${V_{RMS}}^2$) / noise ${V_{RMS}}^2$) that is being measured, unless noise is deliberately subtracted from the signal i.e. received signal (${V_{RMS}}^2$ - noise ${V_{RMS}}^2$) / noise ${V_{RMS}}^2$ to get the SNR? I'm assuming that the traditional $P_{\text{signal}}$ means the signal without the noise and not the received signal as it is seen.

The SNR is simply the mean-square of a demodulated symbol divided by the variance of the signal, or in dB:
$$\DeclareMathOperator{\SNR}{SNR}\SNR =10\log_{10}(\mu^2/\sigma^2)$$
$$=20\log_{10}(\mu/\sigma)$$
A typical metric for this in receiver equipment is the Error Vector Magnitude (EVM) where an "error vector" is the Euclidean distance from the actual sample at the optimum timing location in each symbol to the actual symbol location in a reference waveform (as the distance to closest decision boundary, just prior to decision). EVM is the normalized rms of the error vectors over multiple samples, where the waveform and the decision boundaries are scaled to either the rms amplitude or a maximum amplitude in the constellation, depending on which standard is used. We can therefore use one sample for every symbol to compute a sufficient number of errors to derive a statistical quantity.

See this post for more details on computing EVM from which the above graphic was copied: How to calculate EVM in %age of an Equalized Constellation in 16QAM?, but to describe it simply we statistically estimate with enough samples the standard deviation of this normalized error. The EVM is the average of the noise magnitudes as a percentage of the normalized constellation which represents the signal. (The average of the magnitudes is not the same as the magnitude of the average ($\mu$) so an additional small adjustment factor for this is needed in translating EVM to SNR that depends on the signal constellation used and the normalization that was used for computing EVM).
The relationship between EVM and SNR is therefore:
$$\SNR =  20\log_{10}(1/\text{EVM}) - K$$
Where $K$ is an adjustment factor depending on the peak-to-average of the constellation and how the constellation was normalized to the peak for computing EVM.
This would be affected by the portion of the noise as in the OP's figure that is within the bandwidth of the receiver after demodulation just prior to symbol decision (which is the noise bandwidth we are concerned with: the waveform bandwidth).


Why is the second-order-section cascade approach giving "worse" results? It seems that as a general rule, for digital IIR filters of order more than 2, one is encouraged to use a cascade of second order sections, for reasons connected to the poor stability of high-order filters. See, for example, my question What is the largest "safe" order for the digital Butterworth filter of a given signal? and the answers given there, or the "Limitations" section of https://uk.mathworks.com/help/signal/ref/butter.html.

Now suppose I have a digital signal $x(t)=\sin(2\pi t)$, $t \in \{0,0.08,0.16,\ldots,99.92,100\}$.

Suppose I want to filter this signal in the passband $[f_\mathrm{min},f_\mathrm{max}]=[0.5,3]$ using a zero-phase Butterworth bandpass filter, with filter order $60$; I believe the transfer function of the filter should be
$$ H(2\pi if) \ = \ \frac{f^{30}}{f^{30}+\left( 0.4f^2 - 0.6 \right)^{30}} $$
where 0.4=1/(3-0.5) and 0.6=1/((1/0.5)-(1/3)).

Presumably, the correct result for the filtered signal will be almost identical to the original signal.

In MATLAB, I can perform the filter by a cascade of second order sections, using the code

n=15;
Wn=[0.5 3]/6.25;
[z,p,k]=butter(n,Wn);
[sos,g] = zp2sos(z,p,k);
d=filtfilt(sos,g,x);

or I can perform the filter directly from the transfer function coefficients, using the code

n=15;
Wn=[0.5 3]/6.25;
[b,a]=butter(n,Wn);
d=filtfilt(b,a,x);

The results for the two approaches are as follows:

The upper graph gives the result of the SOS approach, and the lower graph gives the result of the direct approach.

It seems that for quite a significant proportion of the duration of the signal (around the start and end), the SOS method is giving considerably inaccurate results.

    Why is the direct approach giving (what I presume to be) a much more accurate result overall than the SOS approach? Should some caveats be added to the general advice that one should always use SOS?

I don't think that has anything to do with SOS vs transfer function representation, but it appears to be an artifact how filtfilt()  handles initial conditions. The documentation is annoyingly vague about this and simply says it

  filtfilt minimizes start-up and ending transients by matching initial
  conditions.

If you implement filtfilt() manually you get the expected result: both results are basically correct, but SOS has less  numerical noise. See code snippet below.

In general, using filtfilt() is tricky: ideally you have an input signal that has lots of zeros at the beginning and the end, where "lots" means more than the length of the impulse response. Otherwise the results are highly dependent on the initial states of the filters and how you manage them.

Numerically, SOS is always better. Take a look at your "a" coefficients: the ratio of the largest to smallest coefficient is in the order of 1e10. Adding numbers with that different in magnitude creates numerically issues, even if you use double precision. For the SOS representation, that ratios is only about 14, that's about 9 orders of magnitude better!

%% signal
t = (0:.08:100)';
x = sin(2*pi*t);

%% original
n=15;
Wn=[0.5 3]/6.25;
[z,p,k]=butter(n,Wn);
[sos,g] = zp2sos(z,p,k);
d1=filtfilt(sos,g,x);
[b,a]=butter(n,Wn);
d2=filtfilt(b,a,x);

%% manual implementation of filtfilt()
y = filter(b,a,x);
y1 = filter(b,a,flipud(y));
y = g*sosfilt(sos,x);
y2 = g*sosfilt(sos,flipud(y));
plot([y2 y2-y1]); legend('result','difference');


Group delay (MATLAB $\tt grpdelay$) I have been trying to figure out how it is I have an IIR filter with large negative group delay.  It's my understanding that if a filter is causal and stable, the group delay can't be negative at all frequencies.  Is this correct?

The essence of my problem can be seen in the following MATLAB code.  This first filter gives four poles inside the unit circle and four zeros at zero.  The group delay is both positive and negative with a mean of zero.  I believe the contribution of each pole and zero adds and negates one sample of average group delay, respectively.

(1a) grpdelay(1,[1 0 0 0 .1],100000,'whole')

Now, if I create a semi-stable filter, the poles on the unit circle contribute 1/2 sample, so the average group delay is -2.  Actually it's -2 everywhere.  If the poles are outside the unit circle, there's no contribution; so, the average is -4 for the four zeros at zero.

(1b) grpdelay(1,[1 0 0 0 1],100000,'whole')
(1c) grpdelay(1,[1 0 0 0 2],100000,'whole')

OK so far.  Now if I just change the last coefficient to be complex valued, I believe this just rotates the poles, but now there is no contribution for poles inside the circle, poles on the circle subtract 1/2 sample instead of adding, and poles outside the circle subtract 1 sample.  And thus, you get negative group delay everywhere, even when the poles are inside the circle.

(2a) grpdelay(1,[1 0 0 0 .1j],100000,'whole')
(2b) grpdelay(1,[1 0 0 0 1j],100000,'whole')
(2c) grpdelay(1,[1 0 0 0 2j],100000,'whole')

I briefly glanced at the code for grpdelay but did not make a serious attempt to understand it.  Also, if you flip this around and use an FIR filter, changing the last coefficient to complex does not change the average group delay contribution by each zero.

Example (1a):

Example (2a):

This is a bug, either in your code, or in your version of the Matlab function grpdelay. As you correctly pointed out, making the last denominator coefficient imaginary just rotates the pole-zero diagram. This means that for the given example, the zeros stay where they were (in the origin), and the pole angles are rotated by $\pi/8$; or course, the pole radii remain unchanged.

For the group delay function this means that it only gets shifted (along the frequency axis), but its maximum and minimum values (and its overall shape) must remain the same.

Using Octave I get the expected result:

g1 = grpdelay(1,[1 0 0 0 .1],1000,'whole');
g2 = grpdelay(1,[1 0 0 0 .1j],1000,'whole');
plot([g1(:),g2(:)]), axis([0,1000,-.4,.5]), grid on


What is the relation between kernel functions, kernels used in convolution and null spaces of a matrix? I have recently started learning about machine learning and have come across kernels and null spaces. I understand that null space is the set of all vectors that satisfy the equation A.v = 0 (Where A is a matrix). I have been taught that null space is a set of vectors that are squished to 0 when transformation matrix A is applied. Then I came across SVM where kernel functions are used. I read that null spaces are also called kernels of a matrix. My questions are as follows.

Are both kernel functions and null spaces same? If they are related, how are they related?
What is the relation between kernel functions used in SVM and null
    space of matrix? If yes, how? Is this derived from the null space of
    transformations used in SVM?
What is the reasoning behind kernels used in kernel convolution? For
        example: Is the Gaussian kernel used a representation of the
        transformation? How are they related to null spaces? If they are related, how are they related?

Could you please answer these questions? I am very confused. Thank you in advance.

A kernel in the context of digital signal processing refers to the impulse response $h[n]$ of a filter. Particularly for finite length FIR kernels, filtering is carried out by the convolution operator; hence the kernel of the filter is the kernel of the convolution sum. This kernel maps an input vector $x[n]$ into an output vector $y[n]$ i.e.,
$$ y[n] = h[n] \star x[n] = \sum_{k=-\infty}^{\infty} h[k]x[n-k] ~,$$the sum truncated to a finite length for practical FIR convolutions.

As you have described, a null space of a linear mapping is the set of vectors $x$ whose image is the null vector in the range space of the mapping given by the matrix A. As the name implies it is a set of vectors and not a kernel function (sequence). The only kernel in this context is the matrix A which maps input vectors to output vectors, just as in the context of a convolution operator. $$y = A x$$


image processing find liquid level (only with images processing, no AI) I want some ideas to get the level of a water/liquor in a cup. i tried many ideas (contours, edges,...) but nothings is working. So any ideas that i can explore.

Constrain the environment. Best is to make sure the glass is always the same (position and shape), take a picture of the glass empty and use that as reference. Otherwise, you'll need some algorithm to first detect the position of the glass (if the background is constrained, you can use that).
If the liquid is isolated, I'd recommend applying a sobelY filter and something to reduce noise. Take a vertical line out of the graph and identify each step (top glass, water level, bottom glass).


Intuitive interpretation of Laplace transform So I am getting to grasps with Fourier transforms. Intuitively now I definately understand what it does and will soon follow some classes on the mathematics (so the actual subject). But then I go on reading about the laplace transform and there I kind of lose it. What is the moment of a signal? Why is the fourier transform a special case of the laplace transform? How can I come to grips with the Laplace transform?

Ive looked at these sources before I asked this question:

What is meant by a system's "impulse response" and "frequency response?"

How to distinguish between the different frequency domains?

Amplitude vs Frequency Response

Why is the Fourier transform so important?

http://en.wikipedia.org/wiki/Laplace_transform

Why is the fourier transform a special case of the laplace transform?

The Laplace transform produces a 2D surface of complex values, while the Fourier transform produces a 1D line of complex values.  The Fourier transform is what you get when you slice the Laplace transform along the jω axis.  For instance, a simple lowpass filter $H(s)=\frac{1}{s+1}$ has a single pole in the S plane to the left of the origin:

Viewed from the side, the magnitude of this Laplace transform forms a surface, with the pole acting like a tent pole that raises the amplitude to infinity at that point (and an implied zero at infinity that drops the amplitude to zero the farther away from the origin you get in any direction):

If you now take the value of the surface along the jω axis only, that's the Fourier transform.  It's the red curve in the image above, which you can see forms a lowpass filter.  If you moved the pole farther away from the origin, the tent would move in the same direction, and the slice along the jω axis would drop, both reducing the gain (which we compensate for by adding overall gain) and increasing the cutoff frequency.  I've been meaning to make some animations of stuff like this...

http://www.maximintegrated.com/en/app-notes/index.mvp/id/733

https://dsp.stackexchange.com/a/9579/29


Calculate Collision force using Audio Signal I am having an audio sample hitting an aluminium object (a small piece of hard ball 200g) at varying speed.

I am trying to calculate the amount of force the object is hitting with using the audio signal to find the amount of damage its causing to the aluminium

I am able to find any direct or indirect method to calculate that. I am not sure if its easy do it.

Any leads/directions?

Answer:
Whenever the impact happens, that portion of audio has to be detected( if you want to automate it) or else you can identify those sections manually.

To automate it you need to identify characteristics of the audio during the impact, for ex: the frequency contents, or time spread. A good approach to combine these two characteristics is to either do a wavelet decomposition (morlet wavelet) or do a short time Fourier transform.

Once you have a method to identify these regions then you would need a reference (may be using an already available audio meter) to map these signals to audio (in DB), this would depend on how you quantize the analog audio, Sampling rate, FFT size etc. From this identified reference you can then map other audio impacts to suitable audio decibels.


Detecting phones in a word audio I am working on speech recognition for Hindi language. My approach is to map speech to text based on phones. For this, first split the speech into the words and then map each word to text by detecting phones. Till now I am able to split the speech audio file into separated words audio files (I used energy levels for this).

Number of different phones is 1023. For each phone, I have extracted MFCC coefficients for every 10ms frame. (So different phones have different number of frames and for each frame I store 38 coefficients. This serves as templates for each phone.) Also for each word in speech, I extract MFCC coefficients for every 10ms frame.

But I am unable to find efficient algorithm to map word into sequence of phones. I have been thinking of two approaches for this:

Naive bruteforce approach in which I simply try to match every possible group of frames in word to every available template (using DTW) and then select the closest match. But this is very inefficient method.
Split word into phones. Then for each phone select most appropriate template from 1023 available templates (using DTW). But in this approach, I am not getting how to split word into separate phones.

Is there any another approach to solving this problem? Or can it be done efficiently using second approach?

Though acoustic speech features (MFCCs, etc) & lexical words are correlated in daily language, traditionally in speech recognition software, acoustic and linguistic information are assumed independent.

You can see this in popular open-source speech recognition packages like Sphinx3/4 (Java/C using HMMs vs DTW) or Kaldi (C++/neural nets), one must train both an acoustic model & a language model. For the English language, there are a few off-the-shelf models available, but in most domain-specific cases, it is best to train a new one or adapt an existing model. I'm not sure about the state of acoustic/language models for the Hindi language, but here is a paper you may have seen.

I recommend looking for existing Hindi models to use within Sphinx or Kaldi.


What are the various definitions for $\rm SNR$, and their associated methods for measuring it? The definition of $\rm SNR$ seems to be somewhat of a tower of babel in industry. What definitions of $\rm SNR$ are there (feel free to site application), and how exactly can it be measured for that applications?

My specific questions on $\rm SNR$ are:

How may we measure $\rm SNR$ for a communications system, if we have not yet been able to attain optimal bit sampling timing, and all we have to work with are all the signals in the receiver up to and including the envelope of the I and Q channels? See this post for context.
Once we attain optimal bit sampling, and softbits are attained, how best to measure $\rm SNR$ (or $E_bN_0$)? One way I use is:
$$
10\log_{10}\left[ \frac{\textrm{mean}\left\{\lvert s_n \rvert^2\right\}}{\textrm{var}\textrm\{\lvert s_n \rvert\textrm\}}\right],
$$

however I understand this is not suitable for low $\rm SNR$ cases. What other ways exist?

It's not practical to come up with a comprehensive list of definitions of signal-to-noise ratio, because different measures are relevant in different applications. Here are a few common measures that I've come across and/or used in the past (you'll find that they are bent toward communications applications):

"Pure" SNR: I call this "pure" because it is a literal interpreation of the term "signal to noise ratio." This measure is just the ratio of the power in the signal of interest to the total noise power over the signal's bandwidth. Any noise present outside of the signal bandwidth can be filtered out without removing any signal power, so it can be ignored. This definition isn't often used for digital communications; there are more useful quantities that I reference below. SNR is more relevant for analog modulation.

$\frac{E_b}{N_0}$: This is the most common measure of signal to noise ratio used for digital communications. Referred to in the comments above by Dilip as the "BEND" ratio, Eb/N0 is the ratio of the energy per (information) bit ($E_b$) to the (assumed to be white) noise power spectral density ($N_0$) over the signal bandwidth. This is often used as a figure of merit for digital communications systems because it is normalized by the system's symbol rate and modulation format (since it is the energy per bit). This allows very different system implementations to be compared for efficiency on a level playing field, in an apples-to-apples fashion.
In some other specific applications, you might also see $\frac{E_s}{N_0}$ (energy per symbol instead of energy per bit in the numerator) or $\frac{E_c}{N_0}$. I've seen two different meanings for the latter: first, in coded communications systems, you might see this, where $E_c$ stands for the energy per coded bit, in contrast to $E_b$, which is the energy per information bit. Another meaning: in direct-sequence spread spectrum systems, you might see $E_c$ refer to the energy per chip instead.

$\frac{C}{N_0}$: This quantity is often termed "carrier to noise ratio." One place that I've seen this used as a common measure of signal quality is in GPS-related literature. This is very similar to the "pure" definition of SNR, except that the (again, assumed to be white) noise power spectral density value $N_0$ is used in lieu of the total noise power over the signal bandwidth. Since bandwidth is a somewhat squishy term that has many different definitions (there is no one particular measure of what a signal's "bandwidth" is that is accepted to be correct), this ratio is somewhat more precise and less ambiguous than the "pure" SNR metric.


Filtered Signal vs File Compression Paradox 1. Original Situation

I have an original signal as a column data matrix n channels data x:mxn (single), with m=120019 the numer of samples and n=15 the number of channels.

Also, i have the filtered signal as a filtered column data matrix x:mxn (single).

The original data is mainly random, centered at zero, from sensor pickups.

Under MATLAB, am using save with no options, butter as highpass filter, and single for casting after filtering.

save essentially apply a GZIP level-3 compression over a binary HDF5 format, hence we could assume the filesize is a good estimator of the information contents, i.e. maximum for a random signal, and close to zero for a constant signal.

Saving the original signal creates a 2MB file,
Saving the filtered signal creates a 5MB file (?!).

2. Question

How it is possible the filtered signal has a bigger size, considering the filtered signal has less information, removed by the filter?

3. Simple Example

A simple example:

n=120019; m=15;t=(0:n-1)';
x=single(randn(n,m));
[b,a]=butter(2,10/200,'high');
 xf=filter(b,a,x);
save('x','x'); save('xf','xf');

creates 6MB files, both for for the original and filtered signal, which is bigger than the previous values due to using pure random data.

In a sense, indicating that the filtered signal is more random than the filtered signal (?!).

4. Evaluative Example

Consider the following:

A filter created from a random signal $x_r$ from gaussian noise $\sim N(0,1)$, and a constant signal $x_c$ equal to $1$.
Disregard the data type, i.e. let's use only double,
Disregard the data sizes, i.e. let's use one column data vector of 1MB, $n=125000$, $m=1$.
Lets consider the $a$ parameter as the Randomness Index for testing: $x=\alpha x_r+(1-\alpha)x_c$, meaning $\alpha=1$ is fully random and $\alpha=0$ fully constant.
Consider a highpass butterworth filter with $w_n=0.5$.

The following code:

%% Data
n=125000;m=1;
t=(0:n-1)';
[hb,ha]=butter(2,0.5,'high');
d=100;
a=logspace(-6,0,d);
xr=randn(n,m);xc=ones(n,m);
b=zeros(d,2);
for i=1:d
    x=a(i)*xr+(1-a(i))*xc;
    xf=filter(hb,ha,x);
    save('x1.mat','x'); save('x2.mat','xf');
    b1=dir('x1.mat'); b2=dir('x2.mat');
    b(i,1)=b1.bytes/1024;
    b(i,2)=b2.bytes/1024;
    i
end
%% Plot
semilogx(a,b);
title('Data Size for Filtered Signals');
legend({'original','filtered'},'location','southeast');
xlabel('Random Index \alpha');
ylabel('FIle Size [kB]');
grid on;

With the following chart as result:

This simulation reproduces the condition of the filtered signal always having a notorious bigger size than the original signal, which contradicts the fact that a filtered signal has less information, removed by the filter.

+1 on very interesting and insightful experiment.

Some thoughts:

It's not true that filtered signal has less information. It depends on your input signal, filter type,  and cut-off frequency.
When you high-pass the noisy signal, you're removing the slowly changing components. That makes your signal composed of 'more frequently changing random numbers', thus more random. Of course, that depends if your input signal contains high frequencies or not. Your input is noise, so contains every high frequency. But if your input is a more ordered signal, it will lose much of its energy after a certain HP cut-off frequency, the output becomes near zero, less random, less size.  I think if you increase the cut-off frequency of your HP filter pretty high, after a certain point, the file size will decrease.
One other experiment would be to pass the signal through a LP filter with a low cutoff frequency and see the difference.
Based on the same theory in 1., you're high-passing your signal, essentially removing the DC part, xc, and leaving it with noise xr.


Matching two audio inputs (signals) this is for an application in which passwords are stored in the form of audio inputs .
In order to provide user access , the password(in audio form) will be matched with already stored password(audio only) in the database . I wanted to know how to measure the extent of similarity in two audio patterns . the audio patterns can be a clap sequence , knock sequence , whistle etc

Look up acoustic fingerprint or audio fingerprinting and dynamic-time warping for some basic technologies used to compare audio files.


What real object could be modeled by this transfer function? What object from real world could be modeled by this transfer function?
What could parameters b, p1 and p2 stand for?

You transfer function models the behavior of a system. That may be a mechanical system, an electric system, etc. In your case you have a second order low pass filter.

If you calculate the step response of your system, you will notice that the system cannot follow the step. It will either grow fast, with a high slope and then start to decrease its slope when approaching the step, to never ever reach it (at least not in a finite lapse).

On the other hand, the system could grow faster to the step, reach it, then have an overshoot and start descending, reach again the step but continue descending, grow again, etc. This oscillation goes mathematically forever, or at least until its amplitude is negligible.

All this will depend on the position of the system poles.

Lets consider that the system models the behavior of, lets say, an elevator trying to reach the first floor. Here the mass of the elevator, the spring coefficients, etc. all can be used to model the system, and the step (input) is the wish to reach the first floor. The actual response of the system will be the step response, which will model the actual movement of the elevator reaching the first floor.


Signal reconstructed using IFFT: amplitude scale larger than input I am trying to test reconstruction of signal from IFFT.
Below is the python program. I break the signal into overlapping windows and do FFT after multiplication with window function. Then I apply IFFT to the FFT output and then overlap-add to reconstruct the original signal. My problem is the reconstructed signal's amplitude is almost 10 times the input. Why is this so ? How can I scale it ?

import numpy as np
import pandas as pd
import matplotlib.pyplot as plt

winsize=512
shift=16

aircutrawdata=pd.read_csv('input.csv',sep=',',skiprows=8000,skipfooter=3500).values
#%%
vibaxisdata=aircutrawdata[:,8]
data_len=len(vibaxisdata)
window_func=np.hanning(winsize)
fftdata=np.zeros((0,int(winsize/2)))
rec_frame=np.zeros(len(aircutrawdata))
for frame in range(0, data_len, shift):
# =============================================================================
#     if frame>0:
#          break
# =============================================================================
    kiri=vibaxisdata[frame:frame+winsize]
    if len(kiri) != winsize:
        break
    windata = window_func * kiri
    fftframe=np.fft.fft(windata,n=winsize)
    ifft_frame=np.real(np.fft.ifft(fftframe,n=winsize))
    rec_frame[frame:frame+winsize]+=ifft_frame

fig,(ax1,ax2)=plt.subplots(nrows=2,ncols=1)
ax1.plot(vibaxisdata,color='b')
ax2.plot(rec_frame,color='r')
plt.show()

The wave on the top is the original input and the wave at the bottom is the reconstructed signal. You can see that the signal is rather faithfully reconstructed except the amplitude. It is almost ten times the original. Is that how it is supposed to be ?

It's because you are probably adding back overlapping portions together, which adds same signal $N$ times where it seems $N \approx 10$ in your case. To test your algorithm use a zero overlapping case and see that IFFT produces the signal without a gain factor.


Find fourier series of my signal T = 3*(1/25);
fs = 1000;
t = 0:1/fs:T-1/fs;
x = sawtooth(2*pi*25*t,1/2);
x = (x + 1)*0.5;
subplot(3,1,1)
plot(t,x)
grid on
xlabel('Time t (sec)')
ylabel('Amplitude')

I want to take DFT(Discrete Fourier Transform) of my signal using FFT and FFTShift and find frequency spectrum. However,the graphics I found did not seem right to me. Is there anyone who can help? If there are errors in what I do, I want to find and fix it.

Following code for second graph:

y = fft(x);
n = length(x);
f = (0:n-1)*(fs/n);
power = abs(y).^2/n;
subplot(3,1,2)
plot(f,power)
xlabel('Frequency')
ylabel('Power')

First, your period is $T = 1/25s=0.04s$ which means you do not need to take integral of 3 periods for getting fourier series. A plot of the DFT is given below when I took only first $T=0.04s$ is given by

As you can see, all even multiples of $25\text{Hz}$ is zero, while odd multiples of $25\text{Hz}$ is decreasing with increasing frequency. From the plot $|X[5]|^2 = 0.002914$, while $|X[7]|^2=0.0008386$. Even multiples of $25\text{Hz}$ have very low values ($\lt 10^{-12}$). Taking DFT is not the correct method, so read on.

The above values are expected because fourier series of your even symmetric signal with non-zero DC has all $b_k =0$. For computing $a_k$, your signal $x(t)=2t/T$ for $0 \le t \le T/2$ and $x(t)=-2t/T-2$ for $T/2 \le t \le T$. I am not going into details of how to compute fourier series coefficients but I will give hints
$$
a_0=\frac{2}{T}\int_0^T x(t)dt\\
a_k=\frac{2}{T}\int_0^T x(t)\cos(2\pi kt/T)dt \\
=\frac{2}{T}\int_0^{T/2} 2t/T\cos(2\pi kt/T)dt + \frac{2}{T}\int_{T/2}^T (-2t/T-2)\cos(2\pi kt/T)dt\\
=\frac{2((-1)^k-1)}{\pi^2k^2}\\
a_k=\frac{-4}{\pi^2k^2}\text{  for k odd}\\
a_k=0\text{ for k even}
$$
So as you can see, the even harmonics are all zeros, while odd harmonics decrease by a factor of $1/k^2$. So
$$a_5=\frac{-4}{\pi^2 25}\\
a_7=\frac{-4}{\pi^2 49}
$$


How to do smoothing without affecting phase I have signals captured from two channels of a measurement system. There is phase and magnitude differences between these two signals.

I want to apply nonlinear smoothing to the magnitude, but the phases should remain the same after the process.

This is what I have done both channels:

Perform FFT
calculate magnitude
calculate phase
apply nonlinear smoothing to the magnitude
get real and imaginary part from smoothed magnitude and phase
perform inverse FFT.

However the resultant signal is not as desired. It appears as a two sided time signal or has some additional signals in the end.

Is there a solution to this problem?

Any smoothing in the frequency domain is similar to convolution filtering.  And convolution filtering (with anything other an impulse or zero) will spread the width of a signal.

If you want a real result, you may need to check to make sure the input to your IFFT is conjugate symmetric.

If you don't want any phase shift, you could try copying the phase from before your filter into the result after your filtering and before the IFFT.


Understanding the Conditions for Recovering a Discrete Time Signal Through Sampling So, we had a very brief introduction to the Nyquist-Shannon sampling theorem (the discrete time version).

While discussing this, we have seen that multiplying a discrete time signal $x[n]$ by an impulse train $p[n]$, we get a sampled signal $x_p[n]$, and if certain conditions are fulfilled, then it is possible to recover the original signal.

In the frequency domain, what it means (informally) is that we are taking $X(e^{j\omega})$ and shifting it by integer multiples of $2\pi/N$ (i.e the sampling frequency) and "adding" all of this to the plot of $X(e^{j\omega})$ giving us the plot of $X_p(e^{j\omega})$ (in other words we have replicas of $X(e^{j\omega})$)

After looking on the internet and doing a couple of drawings, I've seen that one important condition to recover the initial signal $x[n]$ is that $\omega_s > 2\cdot \omega_m$, where $\omega_s$ is the sampling frequency and $\omega_m$ is the frequency of the original signal (i.e the frequency of $x[n]$).

However, the teacher, in his lecture only mentioned the following condition:

"If the original signal satisfies the following condition: $$X(e^{j\omega}) = 0, \text{ for } \frac{\pi}{N} < |\omega| \leq \pi$$ then the sampled signal $x_p[n]$ contains all the information about the original signal and we can easily recover $x[n]$ from $x_p[n]$".

This makes a little sense to me, as I understand that if $X(e^{j\omega}) \neq 0 $ for all $\omega$, then we do, indeed, have overlapping (i.e aliasing) and it is not possible to recover the signal. But I do not understand why $X(e^{j\omega})$ has to be $0$ specifically for $\frac{\pi}{N} < |\omega| \leq \pi$.

Also, how come the teacher did not speak of the other condition (namely that $\omega_s > 2\cdot \omega_m$) ? Has he forgotten it, or is there something I am missing ?

Lets have a higher level of the idea of signal reconstruction from samples.

When you try to reconstruct something from partial information it is important to know what you know about the result required.

For instance, if you know the signal is a Line if $ \mathbb{R}^{2} $ then all you need is 2 samples, any samples of it.
You can even generalize that to any polynomial with having enough points samples from it.

You can even find more interesting rules (For instance, can you think about the case of having single Sine Function with Zero Phase ad Unit Amplitude, how many point will you need?).

So when dealing with reconstruction we need assumption of the model and then derive the requirements from the samples.

If we assume the signal has limited bandwidth, namely from a given frequency and above it vanishes in the Fourier domain then the Nyquist Shannon Sampling Theorem states a simple rule, sample the signal in a uniform manner with sampling rate which is larger then twice the bandwidth of the signal.
Once you do that, you have enough information to reconstruct it.

So in your case, saying that the signal vanishes from certain frequency and above is just stating it obeys the assumption of the model.


PC-to-PC data communication over Sound FSK / PSK / DSSS For fun and work I've to implement a java based application that can Encode text or data to sound and transmit using PC's speaker and on the receiving side is a microphone with decoder software.

I was thinking of using FSK (syn or async or suggest) for the implementation:

Text (data) ---modulate---> 10khz-20khz audio carrier wave ---> Over the air transmission ---> microphone -> demodulate ---> Display text or data.

My main considerations are:

~200bps bandwidth if not more
resilient to noise upto certain level
preferably 16khz - 20khz carrier wave with 44.1khz sampling (there's less noise in this area and its not really audible to older people and mics, speakers are cheap, readily available in laptops)
Not too complex coding logic.

What Modulation would work best with least coding efforts? Any recommendations on libraries/samples for achieving FSK/BFSK/PSK or DSSS using java?

This is actually a really tough problem because of the channel characteristics. Most computer speakers have fairly limited bandwidth, have significant non-linearities and the room acoustics are often time variant.

Life becomes A LOT easier if you can just run a cable from the headphone output of one PC into the line input of the other.


What's the Difference Between LMS and Gradient Descent Adaptation? I found algorithms that seems the same to me, but they are described with different names (in field of adaptive filtering).

For example:

LMS - least-mean-squares seems to be GD - stochastic gradient descent
Often the stochastic gradient descent is called just gradient descent what seems to be something different (but still similar) according to wikipedia.
Extension of GD (or LMS) with normalization seems to be always called NLMS, why not NGD?

Why are the used names different? Are the algorithms really the same thing or am I wrong?

The LMS algorithm is based on the idea of gradient descent to search for the optimal (minimum error) condition, with a cost function equal to the mean squared error at the filter output. However, it doesn't actually calculate the gradient directly, as that would require knowing:

$$
E(\mathbf{x}[n]e^*[n])
$$

where $\mathbf{x}[n]$ is the input vector and $e[n]$ is the filter output error at time instant $n$. The value of this expectation would typically be something that isn't known ahead of time, so it's instead approximated (in most cases) by the instantaneous value $\mathbf{x}[n]e^*[n]$ instead. This (described in more detail at Wikipedia) is the key feature of the LMS filter structure.

Gradient descent just refers to the method used to hunt for the minimum-cost solution; it doesn't force the use of any particular cost function. An LMS filter is a specialization of gradient descent that uses a mean-squared error cost function and the aforementioned approximation for the gradient at each time step.

Stochastic gradient descent would refer to techniques that don't directly compute the gradient, instead using an approximation to the gradient that has similar stochastic properties (e.g. the same expected value). LMS would be an example of an algorithm that uses a SGD approach, using the approximation I described above.


Interpolation techniques for extremely limited hardware I am trying to playback samples on a retro-computer. The hardware is extremely limited, so I am forced to use software PWM as my main method of playback. The limitations of the platform mean that I cannot use a particularly high PWM frequency and the resulting sample resolution is unlikely to exceed 5 or 6 bits. Let us say just for the sake of definiteness that my PWM frequency is about 20 kHz.

Now, my hardware is also quite limited in terms of memory. Samples recorded at 20 kHz are not really practical. As a way of compromise, I am trying to use samples recorded at 5 or 10 kHz, i.e. at frequencies that are factor of 2 or 4 lower than my PWM frequency. I need to interpolate, but sophisticated processing is basically impossible. I tried simple "upsampling" (i.e. inserting zeros between existing samples) and it produces quite a lot of ringing artefacts. I also tried "zero-order hold" and "first-order hold", as defined in What are the relative merits of various upsampling schemes? Frankly, I cannot easily choose between the latter two; both ways are better than simply inserting zeros, but the ringing in both cases is still substantial.

Clearly, more sophisticated filtering is required. My problem is that software PWM is expensive and my computer is weak. All computations will have to take place in 8-bit registers. There is no hardware multiplication/division (except for multiplications/divisions by 2 done by shifting). A FIR filter can probably be implemented by using look-up tables, but I guess no more that 3-5 filter taps are plausible.

So, the question is, what  filtering methods are likely to give me the most bang for a buck the least ringing for my very limited platform? Should I try to experiment with FIR filters tweaked to get 8-bit friendly numbers or should I do something entirely different? What filtering strategies do you think are likely to work best in the limited computing environment?

In terms of my homework, I did read through the post cited above and several other posts related to it. I've done some reading on Google. I understand how to construct moving averages, including polynomials (="higher order holds"), but the discussions seem to suggest that they are not the best option. FFT is impossible. I am currently playing with various 3-5 tap FIR filters, which seem to work, but are approximately maximum that I can do computationally. I tried to implement CIC filter (and failed to get a decent result). I have not looked at IIR filters yet. So, basically, I am asking this question in the hope that you could help me to focus on something a bit more specific, maybe a family of filters that are more likely to be useful in my case, maybe even specific filters.

For retro hardware, the low pass filter after the PWM is probably the inertia of the speaker cone and/or your ear's internals. So first make sure the ringing isn't beyond your requirements at the full PWM sample rate, or due to quantization (use some tiny 20kHz sample rate waveform loops for testing).

For upsampling, I might try a zero or first order hold followed by an IIR filter (1st order or a biquad low-pass) chosen so that it only has a tiny number of "1"'s in the binary coefficients, allowing shift-and-add recursive filtering with deterministic latency.  If you can find coefficients with a tiny enough number of "1"'s in total, a cascade of biquad IIR's is also possible.

On the data preparation side, you can pre-noise-filter your samples on a more powerful system, instead of just depending on the retro-CPU's quantization.  You an also simulate your retro system's output filter, and tune the noise filtering to move some of the noise outside the passband.


Reconstructing input from the output of a known discrete-time IIR filter I have known discrete-time IIR filter.  I pass a signal x(n) through the filter and obtain the output y(n). How would I design a filter to obtain x(n) from y(n)?
In general, is it possible to design such a filter? If not then under what conditions is it impossible to retrieve the input x(n)?

In general, is it possible to design such a filter?

If by "in general" you mean "for any arbitrary IIR filter", then the answer is no.  In fact, for any real filter, retrieving the actual $x(k)$ in a useful way is problematic, at best.
The basic problem is that a filter filters information out, and leaves you with a residue.  At some point, the initial information is simply lost, and you cannot retrieve it.
About the only time that an inverse filter will work is if you're dealing with some natural phenomenon that shapes the spectrum of $x(k)$ but does not significantly remove any information that you need.

If not then under what conditions is it impossible to retrieve the input $x(k)$?

Mathematically, if you have a filter $H(\omega)$, and it is zero for one or more values of $\omega$, then you simply cannot reconstruct any components of $x(n)$ at those frequencies.
Moreover, if the filter defined by $y(k) = h\left ( x(k) \right)$ has pure delay or unstable zeros (i.e., if it is not minimum phase), then $x(k)$ cannot be reconstructed.  There are techniques for coming up with a best guess of the signal, however.
Practically, the filtering problem is never $y(k) = h\left ( x(k) \right)$.  It is always $y(k) = h\left ( x(k) \right) + n(k)$, where $n(k)$ is some noise process.  So where the filter's response comes close to zero in the frequency domain, you cannot practically reconstruct $x(k)$.

I have known discrete-time IIR filter. I pass a signal x(n) through the filter and obtain the output y(n). How would I design a filter to obtain x(n) from y(n)?

In increasing order of both difficulty and possibility of success:

Cheat, and just measure x(n).
Calculate the transfer function, $H(\omega)$, for your filter.  If the filter is minimum phase, then $\dfrac{1}{H(\omega)}$ will be realizable and stable.  Implement that.  Your output will be $\hat X(\omega) = H(\omega) \left(Y(\omega) N(\omega) \right) $ -- which means that at any frequency where the original filter really attenuates the signal, your estimate of $x(k)$ will be noisy.
Study up, and implement a Kalman filter to estimate $\hat x(k) \simeq x(k)$ from $y(k)$.  If you know the filter structure exactly, this is probably the best way to get an estimate $\hat x(k)$ that is optimal in the least-squares sense.  Assuming your Kalman filter is stable, and the turn-on transient is acceptable, you can derive a steady-state Kalman filter from it.

A corollary: there's probably some variant of the Wiener filter that would do this, too -- but it's been a very long time since I've needed to understand those critters, so I can't say for sure.  Certainly, a steady-state Kalman and a properly constructed Wiener filter -- if it exists -- would basically be the same thing, since they're solving the same problem.
If mean-squared error is not the thing you need to minimize, or if you don't know the filter structure exactly, then some other optimal state estimation method will be your friend.  I'd start with H-infinity filters, but I'm not sure that's where I'd end up.
Yet another wrinkle in this is that -- especially if the system isn't minimum-phase -- if you can get by with a delayed estimate (i.e., you're estimating $x(k - k_0)$ for some delay $k_0 > 0$) then it would be better than any estimate of $x(k)$, and this, too, would drop out of a steady-state Kalman or Wiener filter design exercise.

If you take that last course, note that especially if you're dealing with a non-minimum-phase filter you will be able to get much better estimates of delayed versions of $x(k)$.  I.e., if there's some delay $\kappa$ that is perfectly acceptable, then you can make a Kalman estimator that finds $\hat x(k - \kappa) \simeq x(k - \kappa)$ with less (possibly much less) error than any filter could find $\hat x(k) \simeq x(k)$


Filtering a signal using Autoregressive (AR) filter and finding the coeff of AR filter using Yule Walker equation in MATLAB I have a random signal x of 1000 samples and I've to generate y1 by filtering x using an autoregressive (AR) filter of order n (let say 4). I've tried the following code to implement the autoregressive definition (As in Autoregressive model only the past values of the model output, and the present value of the model input are used)

x = rand(1,1000);
a=0.1
M=4;
for i = 1:length(x)
sum=0;
y1(i)=0;
for j = 1:M+1
    if (i-j>0)
     sum = sum + a*y1(i-j);
    end
end
y1(i)=sum +a*x(i);
end

The problem is that, I need to multiply the previous output y1 and current input x with a=0.1, otherwise it gives an unstable filter which leads to infinity. Please clarify which mistake I'm making while filtering x using AR filter or I'm doing it the wrong way?

Furthermore after doing this step, I need to find y2 by adding some amount of additive white gaussian noise (AWGN) in y1

y2 = awgn(y1,25);

and then using x and y2, I have to find the coefficients of my autoregressive filter using Yule Walker equation for order 4. I've tried this

r = autocorr(y2,4)';
R = toeplitz(r);

w = inv(R)*r;   % as Yule Walker eq: w=R^-1*r

But it's not giving correct coefficients (or filter taps).
Please somebody help me if there's something I'm missing in understanding the question fully or making some mistake somewhere in the code?

Before I go into the specific issues, I'd like to note that there is a minor confusion with respect to M not representing the AR process order. More specifically an AR process of order M=4 would use the M previous values of y1 rather than M+1 values (in the loop with index j). So, with M=4 your implementation is really that of an AR process of order $p = 5$.

Coefficients and stability

  The problem is that I need to multiply the previous output y1 and current input x with a=0.1, otherwise it gives an unstable filter which leads to infinity. Please clarify which mistake I'm making while filtering x using AR filter or I'm doing it the wrong way?

The stability of the AR process can be determined from the location of the poles of the transfer function. To find those poles, we start off with the recurrence equation of the AR(p) process:

$$
  y[n] = a x[n] + a \sum_{j=1}^{p} y[n-j]
$$

From which the transfer function in the $z$-plane is:

$$
\begin{align}
  H(z) &= \frac{a}{1-a\sum_{j=1}^{p} z^{-j}} \\
       &= \frac{az^{p}}{z^{p}-a\sum_{j=0}^{p-1} z^j}
\end{align}
$$
Hence the poles are given as the roots of the $z^{p}-a\sum_{j=0}^{p-1} z^j$ polynomial. In your specific case of $p = 5$, the pole location are given by the following plot:

As you may notice, there is a pole near $z=2$ (more precisely at approximately 1.96595) when $a=1$ which indicates that the system is not stable. On the other hand, when $a=0.1$, all the poles are located inside the unit circle and the system is stable. Beware that even with $a=0.1$, starting with $p=11$, the system also has a pole outside the unit-circle and starts to be unstable.

In other words, there are no specific mistaking in your filtering code (ie. it correctly implements an unstable filter when a=1).

Yule-Walker equation

  But its not giving correct coefficients (or filter taps). Please somebody help me if there's something I'm missing in understanding the question fully or making some mistake somewhere in the code?

There are indeed a few issues with that part of the code. Perhaps the most significant issue is that the solution

$$
\begin{align}
    \vec{a} &= \mathbf R^{-1} \vec r \\
      &= \begin{bmatrix}
           r_0 & r_1 &\cdots & r_{p-1} \\
           r_1 & r_0 &\cdots & r_{p-2} \\
           \vdots & \vdots & \ddots & \vdots \\
           r_{p-1} & r_{p-2} &\cdots & r_0 \\
         \end{bmatrix}^{-1}
         \begin{bmatrix}
           r_1 \\
           r_2 \\
           \vdots \\
           r_p
         \end{bmatrix} \tag{1}
\end{align}
$$
where
$$
\begin{align}
  r_j &= \frac{c_j}{c_0} \\
      &= \frac{E\left\{y_n y_{n-j}\right\}}{E\left\{y^2_n\right\}}
\end{align}
$$

applies to an AR process of the form
$$
  y[n] = x[n] + \sum_{j=1}^p a_j y[n-j]
$$
where $x[n]$ is a zero-mean random variable. Unfortunately, the line

x = rand(1,1000);

Generates a uniformly distributed random variable in the range [0,1] with a mean of $\frac{1}{2}$. An example of random variable with zero-mean could be a uniformly distributed random variable in the range [-1,1]:

x = 2*rand(1,1000)-1;

For an input which does not have a mean of zero, or equivalently for an AR process of the form:
$$
  y[n] = K + x'[n] + \sum_{j=1}^p a_j y[n-j]
$$
where $x'[n] = x[n] - K$ is zero-mean, the equations can be updated to:
$$
\begin{align}
  \begin{bmatrix}
\mu \\
c_1 \\
c_2 \\
\vdots \\
c_p
  \end{bmatrix} &=
\begin{bmatrix}
  1   & \mu & \mu &\cdots & \mu \\
  \mu & c_0 & c_1 &\cdots & c_{p-1} \\
  \mu & c_1 & c_0 &\cdots & c_{p-2} \\
  \vdots & \vdots & \vdots & \ddots & \vdots \\
  \mu & c_{p-1} & c_{p-2} &\cdots & c_0 \\
\end{bmatrix}
\begin{bmatrix}
  K \\
  a_1 \\
  a_2 \\
  \vdots \\
  a_p
\end{bmatrix} \tag{2}
\end{align}
$$
where $\mu = E\left\{y_n\right\}$.

Implementation for zero-mean AR process

Finally there are a few small issues in the implementation of the equation. For one thing, the vector returned by autocorr(y2,p) contains p+1 values: from 0 to p lags. You need the first p values for the $\mathbf R$ matrix (leaving out lag p), and the last p values for the $\vec r$ vector (leaving out the zero lag value since $\vec r$ starts at $r_1$, not $r_0$).

So you can construct the Yule Walker equation with:

R = toeplitz(r(1:end-1)); % r(0), ..., r(p-1)
w = inv(R)*r(2:end);      % as Yule Walker eq: w=R^-1*r

So putting it all together, you should get something like this:

x = 2*rand(1,1000)-1; % uniformly distributed over [-1,1]
a=0.1
p=5;
M=p-1;
for i = 1:length(x)
  sum=0;
  y1(i)=0;
  for j = 1:M+1
    if (i-j>0)
      sum = sum + a*y1(i-j);
    end
  end
  y1(i)=sum +a*x(i);
end

y2 = awgn(y1,25);

r = autocorr(y2,p)';
% Or if you don't have autocorr, replace with the following 2 lines:
% r = xcorr(y2,p)';
% r = r(p+1:end);           % keep only non-negative lags

R = toeplitz(r(1:end-1)); % r(0), ..., r(p-1)

w = inv(R)*r(2:end);      % as Yule Walker eq: w=R^-1*r

Implementation for non-zero-mean AR process

The solution of equation (2) can be constructed in a similar fashion as for equation (1). The main thing to note is that in this case the $\mathbf R$ matrix elements $c_0, \dots, c_{p-1}$ are not normalized by $c_0$ (which would have yielded $r_0, \dots, r_{p-1}$). In this case you may replace autocorr with xcorr, which returns a sum of product without normalization (you would then need some adjustment since is xcorr is double sided, meaning it contains 2*p-1 values, including the values for positive, zero and negative lags). This would give you the following:

x = rand(1,1000);
...

mu = mean(y2);
r = xcorr(y2,p)';
r = r(p+1:end);           % keep only non-negative lags
R = toeplitz(r(1:end-1)); % r(0), ..., r(p-1)

A  = [1 mu*ones(1,p);mu*ones(p,1) R];
b  = [mu;r(2:end)];
w  = inv(A)*b;

% separate the bias term K from the AR process weights:
K  = w(1);
w  = w(2:end);


How to find transfer function looking bode plot I have some Bode plot (only in amplitude) as these:

Now, I must to find their transfer functions, in particular the exact value of their zeros and poles.

It is clear that their are a low pass filter and a high pass filter.

I'm looking for a method for extract the numerical value of their transfer function.

Thank you and bye,
Giacomo

locate the 3db point. (The point where the line before and after the curve would meet.) Take note of the magnitude and frequency f.

Find the angular frequency at this point. w = 2 * pi * f

If the graph is going up then it's a zero. H(s) = w + s

If it's going down then it's a pole. H(s) = 1 / ( w + s )

Once you've got that. check the slope. If it's 20dB/decade then you're set. If it is 40, you'll need to square it, 60 and you'll cube it etc.

Finally you'll need to find the gain constant to make it correct at low frequency (or DC)

So you're answers will look like this:

H(s) = Gain * ( w + s )
H(s) = Gain / ( w + s )


In edge detection, what value is being compared to the threshold value? When we say we set a threshold to the edge finding, what values are we exactly thresholding on? I know a threshold value to limit only certain edges to show up, but what value of the edges is compared with the threshold?

In both the sobel and canny edge detection, they are comparing with a threshold value but I am never sure what values of the edges are being compared.

Its simple. The gradient matrix is computed while performing the Canny/Sobel edge detection. No the gradient varies according to the edge. When you digitize this gradient function, you get a value which tells the magnitude of the gradient at that point.

When you threshold. These magnitudes are compared to the threshold value. If you want weak edges, you keep low threshold values. In such a case all minor edges will show up as well.
If you want only strong 'clear' edges, then you keep higher threshold values.


What is "find-grained" image? What does exactly and precisely the term "fine-grained" mean in image segmentation?
For example, undeniably, the picture below is a fine-grained image:

But is it correct to call below picture as a fine-grained image?

or

I think it means a picture with details. In my idea, by accepting this definition, every picture out there can be called "fine-grained".
Could anybody explain as precisely as possible what exactly a "fine-grained" image is?

Thanks in advance.

Just looking at the Google results I find:

FGComp2013 which states:

  This challenge will target fine-grained classification, i.e. classification among categories which are are both visually and semantically very similar. This is a very difficult regime which is even challenging for humans without careful training, and is critical for establishing a more detailed understanding of the visual world.

The Stanford Dogs Dataset has its primary reference which states:

  First, being a fine-grained categorization problem,
  there is little inter-class variation. For example the basset
  hound and bloodhound share very similar facial characteristics
  but differ significantly in their color, while the
  Japanese spaniel and papillion share very similar color but
  greatly differ in their facial characteristics.

So the upshot is that "fine-grained" means that the classes that are trying to be distinguished between are visually very similar.

In terms of a single image, it means that it is hard to see the differences across the entire image.  In that case, your top image is fine-grained. The dog and bird images are not.  Though if you zoomed into the hair of the dog or the feathers of the bird, I expect they would be.


Understanding an alternative representation for the difference equation of a FIR filter? I recently solved a problem which required that I compute the convolution of an input signal with the impulse response of a filter. I got the resulting discrete signal which is verified to be the correct answer:

  $$2\delta[n] + 6\delta[n - 1] + 3\delta[n - 2] + 3\delta[n - 3] + \delta[n - 4] + 4\delta[n - 5] + 2\delta[n - 7]$$

This is apparently equivalent to:

  $$2\delta[n] + 6\delta[n - 1] + 3u[n - 2] - 2\delta[n - 4] + \delta[n - 5] - 3u[n - 6] + 2\delta[n - 7]$$

Where $u$ is the unit-step.

I understand that the unit-step is zero when $u[n - k]$ has ($k > n$), and one when ($k \leq n$). But beyond that I can't understand why this representation is identical. How can I to transform my output to match this format, or to verify they are the same?

Given that
$$
u[n-k_1] = \sum_{k=k_1}^\infty \delta[n-k]
$$

you can take the difference between two steps with different shifts ($k_1 < k_2$) and obtain:
$$
u[n-k_1] - u[n-k_2]
  = \sum_{k=k_1}^{k_2-1} \delta[n-k]
$$

Then the difference

$$
3u[n-2] - 3u[n-6]
  = 3\left(\delta[n-2] + \delta[n-3] + \delta[n-4] + \delta[n-5]\right)
$$

which appears in your second equation can be replaced to give:
$$
\begin{align}
  y[n]
    &= 2\delta[n] + 6\delta[n - 1] + 3u[n - 2] - 2\delta[n - 4] + \delta[n - 5] - 3u[n - 6] + 2\delta[n - 7] \\
  &= 2\delta[n] + 6\delta[n - 1] - 2\delta[n - 4] + \delta[n - 5] + 2\delta[n - 7]
     + 3u[n-2]-3u[n-6] \\
  &= 2\delta[n] + 6\delta[n - 1] - 2\delta[n - 4] + \delta[n - 5] + 2\delta[n - 7] + 3\left(\delta[n-2] + \delta[n-3] + \delta[n-4] + \delta[n-5]\right) \\
  &= 2\delta[n] + 6\delta[n - 1] + 3\delta[n-2] +3\delta[n-3]+ \delta[n - 4] + 4\delta[n - 5] + 2\delta[n - 7]
\end{align}
$$

which corresponds to your first expression.


2D Haar Transform Computational Complexity I want to learn the computational complexity of the 2D Discrete Haar (Wavelet) Transform (DHWT). The number of operations (divisions, summations etc.) is the main focus. I think there are many different ideas on that. There are even papers claiming that we do not need any multiplications for it. I am confused. Any help? Thanks.

Edit: I should have added the orthogonality requirement for the Haar transform.

This is just to address the lack of multiplications:  why would there be?

The coefficients of the Haar wavelet (in 1D or 2D) are always +1 or -1.  There is no point in doing a multiplication to find the value or its negative.

Check out these notes, particularly this page:


How to perform downconversion using a DSP I read that to perform downconversion, first, I need to generate a sine wave with the same phase and frequency as the transmitter's and then multiply the signal received by the generated sine wave. I don't know how to do that, I think I have a way to generate the synchronized sine wave, but  I don't know how to perform the multiplication.

I read that I can use a mixer or a nonlinear circuit to perform the multiplication, but I have a DSP and it doesn't have a frequency mixer available. Can I do it digitally ? The carrier has 125kHz and the modulating signal has 20kHz. I can use a sampling rate up to 1MHz and the DSP can reach up to 70 MIPS.

I'm sorry if my explanation is incomplete, I am really studyng the problem, so if you need more information you can ask and I will try to make the question more complete.

With my best regards,

Daniel.

Yes, you can do the operation digitally; this is a very common operation in modern software-defined radio systems. The implementation is straightforward: you just implement the same steps that you would see in a hardware implementation. First, multiply by a sinusoid at the carrier frequency, then apply a lowpass filter to remove any images.


Implication of Cramer Rao Lower bound and relation with SNR : non-technical overview I need some clarifications on  CRLB which has been bugging me for quite some time.

Clarification: When computing the Fisher Information(FI), the mathematical expression contains $\theta$ and not $\hat{\theta}$ wnere $\theta$ is the true known parameter that we wish to estimate and $\hat{\theta}$ is its estimate. This means that FI is calculated based on the true, known value of the estimate. So, when implementing the expression of FI for any model, we plug in the true known value of theta and then compare it to the variance (MSE) of the estimate. Am I correct?
Confusion: The Mean Square error (MSE) between the true and the estimated volume may decrease with increasing Signal to Noise Ratio. CRLB is a lower bound for MSE. In that case, would CRLB decrease with increasing SNR?

This means that FI is calculated based on the true, known value of the estimate.

No! The Fisher estimation is based on the expectation of the likelihood function. The likelihood is just a normal probability – just that the role of the "original" variable and observation have been reversed.

I don't understand your example – spheres don't have a CRLB,  the problem of estimating a certain parameter has.

Maybe you're confusing the "set of all estimators for $\theta$" (which has a CRLB) with "a specific estimator" (which simply has a variance) or a with "a specific estimate" (which has an error).
That is, the Fisher info is a function of the "original" $\theta$.

  So, Can Fisher Information be high and does it vary with the noise?

So imagine a binary source: there's only two possible $\theta$, and they are sent through a memoryless AWGN channel – that is, you observe $x$, and it's either caused by $\theta=0$ or $\theta=1$, and figuring that out is your estimator's job (which gives you $\hat\theta$).

The likeliness function $f(x;\theta)$ hence gives you an idea of how likely an outcome $x$ is explained by the original effect $\theta$ was.

Obviously, if there's very little noise, that function will have a very "sharp" shape – there's little chance a 0 sent will be interpreted as 1 in the receiver.

Increasing the noise power makes it possible that this actually happens – the likelihood changes, and it becomes much more probable that the other $\theta$ value is the reason for what you're seeing.


Segmentation of half transparent material, e.g. glass I'm totally stuck on an issue regarding the segmentation of glassy objects.
I need to get the object as precise as possible.
My approaches were different.
At first I tried to remove the background, so that only some sharp contours are left.
But that only works for objects which have sharp edges / gradients.
Otherwise the object itself is also removed.
I posted two different images.

I tried to remove the background via morphological operations, like grayscale dilatation and a divison on it.
but it didnt help much.
after it, I tried a k-means with k=3 for getting the modified background separated from the gray and black values of the glass.
That wasn't successfull in some cases, but not overall/in average.
I also tried to make a canny edge detection with an overall blured filter, but that lead to weaker results in form of open contours, a lot of noise, etc. pp.

Canny with automatic threshold results:

testimg = imread('http://i.imgur.com/huQVt.png');
imshow(testimg)
imedges = edge(testimg,'canny');
imshow(imedges);

Same goes for the second image.

As you can see, there is a lot of noise inside and outside and doubled edges from the glas border. Even there are gaps in the edges.

So, I need your advices for getting a general approach for dealing with this problem of half-transparent materials, not for just these two images.

1) Other ideas for removing the background without damaging the object?

2) Other segmentation methods for getting the object separated from the background?

If it's possible, then with Matlab, IPT or statistical toolbox hints. Any other hints are also welcome!

Thank you for your answer in advance.
Sincerely

Why not just use a simple 2D FFT (guassian) high pass filer?

I did this real quick using MATLAB

Shard #1 using high pass FFT:

The same thing is done on #2.

Shard #2 using high pass FFT:

As you can see, the background and glass area is wiped out, and only the edges are traced.  I did not spend any time on it, but you can threshold the HP filtered output to have more crisp edges, or push the HP cut off higher.

Is this more the results you are looking to get?


What Noise Will Change the Frequencies Components of a Signal? Normally, when talking about adding noise to simulation data, the choice are AWGN, pink noise and so on. However, in practice, the noise in the experiments changes the locations of frequency peaks as well. How can I simulate this kind of noise?

I'm not sure I understand your question.

But usually the noise model is additive.
Namely, due to the linearity of the DFT, at each bin you will have the summation of the signal contribution and the noise contribution at this specific frequency bin.

Now, let's say we have 2 bins.
At the first the signal is 7 at the second is 6.
The noise is 3 and 5 respectively, now the summation will shift the DFT peak though the signal isn't changed.


Current harmonics via FFT -- do peaks matter? Typically, one would perform an FFT on a current waveform and measure it's magnitude at multiples of the line frequency (e.g. 50Hz, 60Hz) to get the harmonics, and subsequently the Total Harmonic Distortion.

If I know that my line frequency should be 50Hz but the FFT peak is not exactly at 50Hz (say, it is 55Hz), then what is my fundamental?  50Hz or 55Hz?
How about the harmonics?  Should I then measure at multiples of 50Hz or 55Hz?
If say my fundamental peak is exactly 50Hz, but my second FFT peak is at 110Hz, where should I measure the 2nd harmonic?  100Hz or 110Hz?  In other words, should I seek out the peaks?
Lastly, THD is usually defined to take the sum at multiples of the fundamental.  What if my FFT peaks are not exactly at multiples of the fundamental?

Update:
- I was trying to measure THD for my power supply with line frequency of 50Hz. If measure the peak of the FFT, it's 50.35Hz. If I record that as the fundamental, the fractional points add up in the subsequent harmonics (e.g. for the 29th harmonic, the value becomes 1460.26 vs 1450 if 50Hz was used as the fundamental).  From the FFT plot, clearly 1450 is the peak, while 1460 is not.
- Should I look for peaks at each harmonic?

You can interpolate the FFT bin results to get better estimates of the actual peaks in the spectrum.  Use a high quality interpolator, such as a windowed Sinc kernel.


Echo hiding technique I'm trying to figure out how to use DSP techniques in order to both add artificial echoes to a sound file, but more importantly to be able to detect these echoes. I've been reading about it here:
http://www.snotmonkey.com/work/school/405/methods.html#echo

Can someone point me to possibly some example code using Apple's Accelerate framework for being able to detect hidden echoes in sound coming from the microphone?  Also, for detecting the echoes, it is important that the application doesn't need to have knowledge of the original sound, meaning it should be detectable in any sound that's being played that has these echoes.

Thank you for any information.

So far I have not experimented with audio steganography, but anyway, maybe these ideas are of some use for you.

Adding echoes to a signal is pretty easy as this process just adds a delayed version of the signal back to itself. You probably want to use a gain factor $f_\text{gain}$ for the additional signal that is smaller than one, $f_\text{gain} \leq 1$, to avoid having an echo that is louder than the original sound. In DSP terms this process corresponds to FIR filtering with an impulse response that consists of a delta impulse at time zero ($\delta[n]$) and a time-shifted and weighted delta impulse:
$$
\mathbf{b} = \delta[n] + f_\text{gain}\cdot \delta[n-L_\text{delay}],
$$
with $L_\text{delay}$ the delay time in samples
and
$$
\mathbf{a} = 1
$$
because it's an FIR filter. To get multiple echoes simply use more than one delayed delta impulses.

For detecting hidden echoes in the signal, maybe it's worth a try whether you can achieve satisfactory performance with an adaptive filtering algorithm. For example, the (normalized) least mean-squares ([N]LMS) algorithm comes to my mind. I'm sure you can find implementations in a number of programming languages to download from the internet. Using the signals both as input signal and desired signal to the adaptive filtering algorithm would lead to a self-cancelling behaviour. Maybe, from analyzing the filter coefficients, you can find out about whether there are echoes in the signals and, if so, what the delay times are.

If you experiment with the Cepstrum method described in the web page that the link in your question points, take care to have your signal blocks long enough to capture the echoes.


Spectral Analysis vs. Spectral Line Analysis We usually talk about "spectral analysis" but some resources (this paper or this doc) talk about "spectral line analysis".

Does this make sense to you, i.e. are the 2 fields actually different or the 2 names refer to the same thing?

Spectrum analysis is more general: it involves looking at the entire spectrum of a given signal.

Spectral line analysis assumes that the spectrum contains several peaks (lines) of interest at specific frequencies.  The aim then is to find the precise frequency, magnitude, and phase of those peaks (lines).


Convolution sum vs auto-correlation problem I had an exam one of these days and one of the questions was:

"Knowing that an auto-correlation estimator of a sinal x[n] could be defined by:

$$
R_{xx}[k] = \sum_{n=-\infty}^{+\infty}h[n].x[n+k] , k=0,1,...,n-1
$$
and that the non-periodic convolution:
$$
y[n]=x[n]*h[n] = \sum_{k=-\infty}^{+\infty}h[k].x[n-k]
$$
show how to calculate the autocorrelation with a convolution operation"

And here is my reasoning:

let's say that:
$$
h[n] = x[-n]
$$
Therefore, we would have some $y[n]$ for the following convolution:
$$
y[n] = x[n]*h[n]=x[n]*x[-n] = \sum_{k=-\infty}^{+\infty}x[k].x[k-n]
$$
So it would be nice to invert the signal $y[n]$, finding then a $y[-n]$. However, we will get the same result as before since convolution is a commutative operation:
$$
y[-n]=x[-n]*h[-n]=x[-n]*x[n]=x[n]*x[-n]=\sum_{k=-\infty}^{+\infty}x[-k].x[-k+n]
$$
Therefore, $R_{xx}[n]=y[n]=y[-n]$, since the convolution sum is symmetric with respect to zero.

Question:

Is my reasoning correct by formal mathematics?
The professor gave me a zero in this answer because he says that my assumption $R_{xx}[n]=y[n]=y[-n]$ is wrong because the correct calculation of the autocorrelation has a complex conjugate in one of the terms (so it might not be symmetric with respect to zero). I agree with it even though he has not placed this formalism in the autocorrelation definition. So, following the exam's definitions, is my reasoning correct?

Also, can I say that this will always be true whenever the signal is purely real?
And by consequence of the question above: can some of these be demonstrated for complex signals?

There are two things to consider. First, looking only at the exam question you can say that the convolution of $h[n]$ with a sequence $y[n]$ is given by

$$(h*y)[k]=\sum_nh[n]y[k-n]$$

Now define $y[n]=x[-n]$:

$$(h*y)[k]=\sum_nh[n]x[n-k]$$

Consequently, we have

$$\sum_nh[n]x[n+k]=(h*y)[-k],\quad y[n]=x[-n]$$

So you need to flip the signal $x[n]$, compute the convolution, and flip the result. This is true, given the sum in your exam question. No need to worry about complex conjugates because we just used the definition of the autocorrelation sum in the question.

The other issue is the definition of the autocorrelation for complex sequences, and how it can be implemented by convolution. The usual definition is

$$R_{xx}[k]=\sum_nx^{*}[n]x[n+k]$$

The convolution of $x^*[n]$ with a sequence $h[n]$ is

$$(x^**h)[k]=\sum_nx^*[n]h[k-n]$$

For $h[n]=x[-n]$ we have

$$(x^**h)[k]=\sum_nx^*[n]x[n-k]=R_{xx}[-k]$$

So for complex signals, you need to compute the convolution sum of $x^*[n]$ with $x[-n]$, and flip the result to arrive at the autocorrelation.

Note that $R_{xx}[k]=R_{xx}^*[-k]$, so instead of flipping the result, you can also use its complex conjugate.


Sine Wave spectrum showing multiple frequencies - Audacity I am generating a sine wave using audacity
with 2000 Hz, 30s duration, amplitude 0.5

48 kHz is the sampling rate of the recording device.

I just play and record the generated sound, this is the spectrogram that I get. (Top is original while bottom is recorded) Frequencies other than 2kHz are also showing up in the spectrogram with significant power density. (Using 2048 point hanning window for FFT).

This is the time domain representation for the generated and recorded signal, this is shown to ensure clipping is not happening.

Why does the spectrogram show these higher frequencies and how could I remove them.

I could not find answer to my query on other similar posts in DSP stack exchange thus am creating this post.

Following could be the reason for extra frequencies -

Non Ideal characteristics of speaker used to play the generated sine wave.
Non Ideal characteristics of microphones used to record the played sine wave.

Also one needs to make sure that given sampling rate(48kHz) is supported by the speaker/sound card and the recording microphone otherwise extra frequency components will be introduced.


Impulse response if initial conditions are given In most books, impulse response of LTI systems are calculated by assuming initial rest condition,  but how to calculate response of an impulse input if there are some initial conditions present ?
My approach - first calculate zero input response (only using initial conditions that are given) and then calculate zero state response (taking impulse as input) assuming zero initial condition and then add both of them and we got total response.
Problem
But I think this method will give wrong answer because due to impulse input there must be some changes occur in initial conditions, so I wanted to know

Due to change in initial conditions, which equation will change, zero input or zero state or both?

What is best approach to get solution of problem like this?

Is by using Laplace transform this problem can be solved
without any change in initial condition?

The reason you're having trouble with this is because the impulse response of a system tells you something considerably different from a system's behavior when it is starting with non-zero initial values -- so no one bothers to treat this in the literature.

Problem But I think this method will give wrong answer because due to impulse input there must be some changes occur in initial conditions, so I wanted to know

A way to reconcile the ambiguity you see (and, typically, something you need to do if you're getting down and dirty with the details around the Dirac delta functional) is to make a distinction between $t = 0$, $t = 0^+$, and $t = 0^-$, where $0^+ = 0 + \epsilon$ and $0^- = 0 - \epsilon$ where $\epsilon$ is that famous infinitesimal from calculus that's greater than zero but arbitrarily small.
Then recall that the Dirac delta functional is zero for $t < 0^-$, and zero for $t > 0^+$, and -- strictly speaking -- undefined for $t = 0$

Due to change in initial conditions, which equation will change, zero input or zero state or both?
What is best approach to get solution of problem like this?

With a bit of hedging about with what you mean by initial conditions, your original idea -- that you can just add the impulse response and the system behavior with initial conditions -- is correct.
You lead yourself to this by recalling that if $\delta(t)$ is undefined at $t = 0$, then so are its effects.  So to be entirely picky and technical, the impulse response is only defined for $t > 0^+$, not before.
Then you hedge your initial conditions statement by specifying the initial conditions at $t = 0^-$.  Without that $\delta(t)$ in there, the state of the system at $t = 0^+$ would equal its state at $t = 0^-$, so (to an engineering level of rigor) this works.
With these details specified, you can just add the two responses.

Is by using Laplace transform this problem can be solved without any change in initial condition?

You can use the Laplace transform as an aid to solving this, but at least the way that I was taught to insert initial conditions into a Laplace-domain solution was to use Delta functionals and their derivatives.  You are, essentially, doing the same thing I did before, it's just that instead of taking the initial conditions at $t = 0^-$, you're finding the correct weights for $\delta(t)$, $\delta^2(t)$, etc., so that the states of your system match your initial conditions at $t = 0^+$.
Then to add in the impulse response, you just add in an extra one times $\delta(t)$.  The effect will be the same as your proposed solution.


Length of DFT defines whether I see the harmonic or the dominant frequency I have a number of signals that are periodic. I use an fft transformation to obtain the dominant frequency of each signal. In order to increase the frequency resolution I zero pad the signal before the fft. However I noticed that when the frequency is low, the peak of the fft corresponds to the first harmonic. A sample of my signal looks like this:

My code is the following:
N = 4000; %signal length
minf= 0.5; % Hz
maxf= 20; %Hz
resf = 0.01; % Hz resolution
fs = 2000; % Hz sampling frequency

NFFT = fs/resf;
Freq = fft(mysignal,NFFT);
Freq = Freq((minf/resf)+1:(maxf/resf)+1,:);

In this case the fft will perform well and detect the 1.12Hz frequency. However I noticed that when my signal's frequency is below 1Hz, as dominant frequency I will get the 1st Harmonic as the pick of the fourier domain plot and not the dominant frequency. In that case I can change the line  NFFT = fs/reso_freq; to NFFT = fs/reso_freq/2; which will actually give me the correct frequency, but it will not work for when my signal is above 1Hz. Can someone explain to me why this is the case and if there is a universal way of finding the correct frequency without knowing what to expect beforehand.
Thank you

Zero padding does not increase the frequency resolution— it just interpolates more samples in frequency but does not give you any more information. The frequency resolution in Hz of an non-windowed (rectangular window) Fourier Transform is 1/T where T is the length of your time domain waveform in seconds. You must have a longer data sample in order to get higher resolution or see lower and lower frequencies. An extreme case is that you need an infinite amount of time in order to observe a true DC signal.

Note that zero padding gives you more samples of the Discrete Time Fourier Transform (The DTFT is the result if the time axis was allowed to go from minus infinity to plus infinity; by zero padding you approximate that). The Discrete Fourier Transform in contrast has a time axis this only goes from 0 to N-1.


Intersymbol Interference -- Root raised cosine filter redundant if doing impulse response equalization? I'm implementing DBPSK (differential binary phase-shift keying) over audio (44.1kHz sample-rate) using software.  Eventually I'll use multiple audio frequencies at once, so I'm using FFTs and inverse FFTs.

As expected, as the window size (symbol length) gets smaller, ISI becomes more and more of a problem and the BER (bit error rate) increases.

It seems like using RRC (root-raised cosine) filters is an established way to reduce ISI.

If as part of my data transmission I send out a special known signal that can be detected, then presumably I can determine the impulse response of the channel and apply an inverse transformation to re-construct (ignoring noise) the original signal without multi-path distortion.  Why is RRC useful if this transformation can be done?

There are three reasons to avoid ISI through pulse shaping rather than correcting it via channel correction methods.

As Bruce pointed out, it is a simpler solution and requires fewer computations.
Channel correction usually amplifies the received noise.
Received noise will corrupt the channel model, making the ISI removal imperfect to some degree.


Finding minimum in the middle of a histogram I'm building an image processing pipeline and one step in the process requires thresholding the image based on saturation.

Example code using OpenCV-Python:

H, S, V = cv2.split(img)
mask = cv2.threshold(S, low, 255, cv2.THRESH_BINARY)

I know how to set the low value manually, if I look at the histogram:

(the red arrow indicates the value I would pick)

My question is: how can I programatically find that value?

Some non-rigorous descriptions of the problem:

find the lower boundary of the right-most "hump"
find the minimum value that's not at the edges

To me, this looks like you want something like

start from right. Add up bins' values, as long as

the added bin is larger than the last or
the value added by that bin divided by the current sum is still above a fixed threshold

Adjust the threshold to fit your needs.

Often, things like smoothing the histogram (effectively: applying a low pass filter to it) help. Or you do something like finding a function fit, something like a polynomial of fixed degree (let's say 5) that minimizes some error function (for example, sum square absolute error in all bins), and then analytically find local extrema.


Changing limit of integral make time-variant to time-Invariant system $S_1:x(t)\longrightarrow y(t) =  \int_{-\infty}^{3t}{x(\tau)\,d\tau}$

$x(t-t_0)\longrightarrow y_1(t) =  \int_{-\infty}^{3t}{x(\tau-t_0)\,d\tau}$

$=  \int_{-\infty}^{3t-t_0}{x(z)dz} \neq y(t-t_0) $    (Time Variant)

My teacher showed my the Above relation, shows that the system is Time Variant but if i change its upper limit to $3t-3t_0$ i.e $\int_{-\infty}^{3t-3t_0}{x(z)dz}$ then the system becomes Time-Invariant, my teacher told me this and now i am confused can someone tells me why upper limit $3t-3t_0$ this makes it time-Invariant

The system-1 whose input/output equation is
$$ y(t) = \mathcal{T}\{x(t)\} = \int_{-\infty}^{t} x(\tau) d\tau $$ is time-invariant.

Whereas the system-2
$$ y(t) = \mathcal{T}\{x(t)\} = \int_{-\infty}^{3t} x(\tau) d\tau $$ is time-varying.

On the other hand the system-3
$$ y(t) = \mathcal{T}\{x(t)\} = \int_{-\infty}^{3t-3t_0} x(\tau) d\tau $$ is still time-varying for any (fixed) value of $t_0$.

I presume that your teacher have tried to underline the fact that what makes the system-2 time-varying is the fact that $y(t-t_0) \neq \mathcal{T}\{x(t-t_0)\}$ which becomes as the equation shows:

$$y(t-t_0) = \int_{-\infty}^{3(t-t_0)} x(\tau) d\tau  = \int_{-\infty}^{3t-3t_0} x(\tau) d\tau  $$

and
$$\mathcal{T}\{x(t-t_0)\} = \int_{-\infty}^{3t} x(\tau-t_0) d\tau = \int_{-\infty}^{3t-t_0} x(\tau) d\tau $$

These two would be equal if th upper limit of the second equation would be $3t-3 t_0$. That's what your teacher told. Since they are not equal therefore the system-2 is time-varying.

Note that artifically making the limit as $3t-3t_0$ is quite meaningless, as the fixed value of $t_0$ has nothing to do with the variable amount of shifts that occurs in $x(t-d)$ or $y(t-d)$.


What does $[-\pi \; \pi]$ band-limited signal mean? I know that a band-limited signal is a signal which have a limited set of frequencies.

What is then the  relation between band limiting and $[-\pi \; \pi]$  since $-\pi$ and $\pi$ are related to angles and frequencies are in Hertz ?

Now first of all concept of bandlimited signal exists in the context of analog signals not digital ones. Digital signals are always (by their definition and existance) bandlimited.

Second, frequency in analog signals is expressed either as $f$ in Hertz (cycle per second, repetation per second) or as $\Omega$ in radians per second, where they are related as $$\Omega = 2 \pi f$$ So for example a frequency of $f=1$ kHz would equivalently make a radian-frequency of $\Omega=2 \pi 1k \approx 6283$ rad/sec.

When a bandlimited analog signal $x_a(t)$ is sampled by a period of $T_s$ (and digitized afterwards in practice), the resulting discrete-time sequence $x[n]$ will have a discrete Fourier transform expression of
$$ H(e^{j\omega}) = \frac{1}{T_s} \sum_k X_a( \frac{\omega + 2\pi k}{T_s})$$

such that it will be frequency normalized to the discrete-time frequency range of $-\pi \leq \omega < \pi$. Where $\omega$ is the discrete-time frequency in radians per sample. Where $X_a$ is the continuous-time Fourier transform (frequency content) of the analog signal $x_a(t)$.

for $k=0$ the base period of $H(e^{j \omega})$ becomes $$ H(e^{j \omega}) = \frac{1}{T_s} X_a( \frac{\omega}{T_s}) $$ where the frequency normalization becomes more apparent as setting $ \Omega = \frac{\omega}{T_s} $ yields the those analog radian-frequencies for $X_a(\Omega)$ when $\omega$ of $H(e^{j \omega})$ is selected from $\omega \in [-\pi,\pi)$


Low frequency and high frequency noise of ecg Is baseline wander in ecg is low frequecny noise in ecg and power line interference is high frequency noise in ecg? Is both the artifacts present in EMG signal?

The frequency content of baseline wander is usually in the range below 0.5 Hz, while power line noise is somewhere around 50/60 Hz accompanied by some harmonics.


Uncorrelated and zero mean image noise Consider a noisy image g(x,y) formed by the addition of noise n(x,y)
  to an original image f(x,y); that is

g(x,y) = f(x,y) + n(x,y)

  where the assumption is that at every pair of coordinates (x,y) the noise is
  uncorrelated and has zero average value.

What does it mean for noise to be uncorrelated and have zero average value. I know the term uncorrelated for random variables, and mean of probability distribution function, but I can't understand these concept for images.

How can noise have zero expected value at every pair of coordinates (x,y), or be uncorrelated ?

Let us consider two random variables $X$ and $Y$. You are aware of the following terminologies with respect to (wrt) the random variables:

(Zero Mean) : $\mathbb{E}(X) = 0$ and $\mathbb{E}(Y) = 0$,

(Uncorrelated) : $\mathbb{E}(X Y) = \mathbb{E}(X) \mathbb{E}(Y)$. It is interesting to figure out how it is different from independence.

Now, an image specified, say $g(x,y) = f(x,y) + n(x,y)$ is a discrete random random process. The notion of zero mean and uncorrelated are not very different here, as for any two points $(x,y) = (x_1,y_1)$ and $(x,y) = (x_2,y_2)$, the image intensity at these pixel locations are random variables.


One more inverse z transform ... a bit more complicated Thank you for all the help teaching me. I am looking at one more inverse z transform and not understanding what to do with it.

$F(t) = m \cdot a(t)$ (i.e. Force = mass $\times$ acceleration)

$F(s) = m \cdot d(s) \cdot s^2 $ (i.e. Where d(s) is displacement and initial position is 0)

Using $s = \frac{1-z^{-1}}T$:

$F(z) = m \cdot d(z) \cdot \frac{1-z^{-1}}{T}^2 $

$F(z) = m \cdot d(z) \cdot \frac{z^{-2} - 2z^{-1} + 1}{T^2} $

I have displacement calculated per sample already from another equation, so I can just substitute that in for $d(z)$. But what on earth do I do with that to get $f[n]$?

Thanks and sorry for all the questions. I unfortunately have zero education on this and it's hard to find simple to understand resources explaining how to do these things.

Okay it's:

$F[n]=m \cdot \frac{d[n−2]−2d[n−1]+d[n]}{T^2}$


Binary and Quadrature Phase Shift Keying I was learning Modulation Techniques and was on the topic Phase Shift Keying (PSK). As far as I have understood,

In Binary PSK we encode a single bit on a one signal element and we differentiate between two bits (or two signal elements) by their different phase.  And that it's favorable to keep the phase different maximum between the two different signal elements, so we use a phase difference of 180 degrees.

In Quadrature PSK, we want to transmit two bits per signal element. Since there are four different possible combinations of two bit numbers, we want four distinguishable signal element to distinguish between the four combinations. So we use signal elements with four different phases and with a difference between them of 90 degree (which is maximum possible).

Untill here, I seem to understand it. Correct me if anything from above is wrong.

I read from a source that in QPSK, the four different phases are 45°, 135°, 225°, and 315°. The excerpt from the said text is as under

  It makes sense to seek maximum separation between the four phase options, so that the receiver has less difficulty distinguishing one state from another. We have 360° of phase to work with and four phase states, and thus the separation should be 360°/4 = 90°. So our four QPSK phase shifts are 45°, 135°, 225°, and 315°

My first question is, aren't phases relative? What are these phase shift of 45°, 135°, 225°, and 315° relative to? My common understanding makes me wonder wouldn't just four sine waves, one not shifted at all, second sifted 90° relative to the first, third shifted 180° relative to the first, fourth shifted 270° relative to the first make more sense and be appropriate to represent different data elements?

The same text states that

  The term “quadrature modulation” refers to modulation that is based on the summation of two signals that are in quadrature. In other words, it is I/Q-signal-based modulation. QPSK is an example of I/Q-based modulation.

My second question is, how exactly is four different phases achieved by modulating the amplitude of just two- I & Q signals? At first, it seems to me that I/Q signaling can only be used to shift a signal 90° (i.e., 45° in each direction).

Also, Whether Binary PSK is also done using I/Q-based signals?

There's another question I and Q components and the difference between QPSK and 4QAM which (maybe) addresses my problem, but it's too complicated and mathemaical for a beginner like me.

(1) The phases are indeed relative. The numbers 45, 135, 225 and 315 are relative to the transmitter's oscillator. The receiver, of course, has a different oscillator, so it will need to do some procesing to estimate the phases.

Related to this subject, there is also differential modulation, where 1s and 0s are transmitting by changing the phase relative to the previous symbol.

(2) This identity is useful to understand I/Q modulation:
$$A_I \cos(2\pi f t) + A_Q \sin(2\pi f_c t) = S_{IQ} \cos(2\pi f_c t - \theta_{IQ}) $$
where
$$ S_{IQ} = \sqrt{A_I^2 + A_Q^2}$$ and $$\theta_{IQ} = \tan^{-1}\frac{A_Q}{A_I}.$$
So, as you can see, you can achieve any phase $\theta_{IQ}$ by adding the I and Q signals with appropriate amplitudes.

(3) Binary PSK is just an I signal. Since the two phases in BPSK are zero and 180, the Q portion of the signal is always zero.


Explain the formula of convolution and convolution theorem As stated in the title, I have two questions

Why convolution is defined as $(f*g)(x) = \int_{-\infty}^\infty f(t) g(x - t) dt$ instead of just $\int_{-\infty}^\infty f(t) g(x + t) dt$ ? Why we need to flip $g$ ?
Is there any intuitive explanation of convolution theorem (Fourier transform) ?

Because the front part of the input hits the filter first.

The graphical diagrams shown in most text books with an exponentially decaying signal multiplied by a step function show this directly.

While Causal response is not required for valid convolution, it is essential that the filter output responds and evolves from the head of the input.

Convolution also commutes, so the head of the impulse function also hits the signal first.

Correlation overlays the pair
 of functions which why it is typically used in measurements associated with similarly measures.


generate a repeated square pulse train in python Generate a pulse train with -2dBFS pulses of 20 sample duration every 8 seconds. I have referred here to generate a pulse train but not able to repeat it. How can I repeat this audio signal every 8 second and between aforementioned signal there has to be silence?

from scipy import signal
import matplotlib.pyplot as plt
import numpy as np

sampling_rate = 16000
time_seconds = 1
number_of_pulses = 20
pulse_frequency = 1/(time_seconds/(number_of_pulses/2))
t = np.linspace(0, time_seconds, sampling_rate, endpoint=False)
signal_dbfs = 0.8 #-1dBFS
plt.plot(t, signal.square(2 * np.pi * pulse_frequency * t) * signal_dbfs)
plt.ylim(-2, 2)
plt.show()

Actually, as mention by @applesoup your question is not clear but as far as I understood you can repeat the signal using

numpy.matlib.repmat(a, m, n)
Repeat a 0-D to a 2-D array or matrix MxN times.
Parameters:
    a : array_like
    The array or matrix to be repeated.
    m, n : int
    The number of times a is repeated along the first and second axes.
Returns:
   out : ndarray
   The result of repeating a.

For further better understanding, you can look this link


Random signals as power signals Why are random signals considered as power signals (i.e. signals with infinite energy and finite average power)?

Does this make any sense? What does it mean for random signals to have infinite energy even though we know that real-life signals (usually with inherent randomness in them) have finite energy!

Note that the condition

$$\int_{-\infty}^{\infty}|f(t)|^2dt<\infty\tag{1}$$

(i.e., that the signal $f(t)$ has finite energy) is very restrictive when we try to model signals, even though obviously any actually occurring signal must have finite energy. Modeling signals as random processes means that we ignore condition $(1)$. Models are always unrealistic to a certain degree, but many signals can be described very well by random processes even though the signals have finite energy and their models do not. This aspect of the model is often irrelevant.

One example which may serve to understand this fact a bit better is the frequently used model of a (wide-sense) stationary process. Certain statistical properties of such a process do not change over time, and, consequently, realizations of such a process will generally not decay as $t\rightarrow\pm\infty$, and $(1)$ will generally not be satisfied, even though we are only interested in the properties of that process during a certain finite time window. However, power and the power spectrum can be defined for such processes, and most practically useful processes have finite power (or can easily be made to have finite power).


OFDM symbol synchronization In an OFDM receiver, how could we find the start of the symbol for FFT operation.

I know the CP-based method to find the start of symbol. However, what if the sampling is off by some duration. We will not get the same digital signal at the transmitter. How serious is that issue? Is it simply addressed in the channel estimation in the frequency domain?

Thanks in advance.

CP based Symbol Boundary detection algorithm which is fundamentally based on correlation will still work. You just need to do a timing recovery before FFT operation.

A small offset in sampling times will cause a small phase rotation in pilots and these small rotations are used to compute the timing offset.

Once you have a good timing offset estimate, you do timing recovery on the received samples of OFDM symbol. Timing recovery is nothing but doing some low cost polynomial interpolation to get the actual values of the samples at correct sampling times.

Or in OFDM once you have these timing offsets as phase rotation, derotate rest of the data symbols after FFT before slicing to get the complex QAM symbols.


Help in understanding if the maximum likelhood estimation is working properly I am learning estimation theory and need  help in understanding for educational purpose how the concept of ML works with the help of a step by step implementation. I am trying to find out the ML estimates of the coefficients of a moving average model. I am unable to understand what I need to do in order to get the log-likelihood curve, so that I can select the coefficients based on the curve. There are implementations in the internet but I wanted to make my own in order to get a clear picture. The way I have done is for the model

 y(i) =  0.6*x(i-1) ;
 z(i) = y(i) + w(i)

where the true coefficient is h = [ 0.6], the input is x and w is zero mean additive white gaussian noise with variance 1. I have used 4 different values of the estimates out of which one value is the true. If the implementation is correct, then I should get estimated coefficient as 0.6 because I am using the true known coefficient in the estimation procedure. However, this is not so. Where am I going wrong?

The maximum likelihood estimator will give you precisely the true value for the parameter for a vanishingly small number of instances --- unless the noise is zero.

Your noise instance w is not zero-mean, which will probably give more of a bias away from the true value.

If I take your numbers are just change the last entry in w from -0.65 to 0.1 then I get the "right" answer. Note that this changes the mean of w form -0.16 to -0.042 (i.e. closer to zero).

Your values:

The new values: