Using the same state transition information as this answer to another question, but using:

y(t) = round(H*x_truth(:,t) + rand(1,1,"normal")*sqrt(R));

as the signal model's output equation, we can apply the same Kalman filter.

This is not really accurate, because the round function is a nonlinearity sort of like quantization. However, quantization can also be modeled as an additive noise, so we'll proceed.

The results are shown in the plot below.

Here,

the black line is the true position,
the red + signs are the quantized, noisy position measurements, and
the green line is the Kalman filter's estimate of the position.

Is this the sort of "smoothing" you're interested in?

And the error between the true position and the Kalman filter's estimate is:

Including full scilab script for reference.

// Signal Model
DeltaT = 0.1;
F = [1 DeltaT; 0 1];
G = [DeltaT^2/2; DeltaT];
H = [1 0];

x0 = [0;0];
sigma_a = 0.1;

Q = sigma_a^2;
R = 0.1;

N = 1000;

a = rand(1,N,"normal")*sigma_a;

x_truth(:,1) = x0;
for t=1:N,
    x_truth(:,t+1) = F*x_truth(:,t) + G*a(t);
    y(t) = round(H*x_truth(:,t) + rand(1,1,"normal")*sqrt(R));  // <== Changed line!!
end

// Kalman Filter
p0 = 100*eye(2,2);

xx(:,1) = x0;
pp = p0;
pp_norm(1) = norm(pp);
for t=1:N,
    [x1,p1,x,p] = kalm(y(t),xx(:,t),pp,F,G,H,Q,R);
    xx(:,t+1) = x1;
    pp = p1;
    pp_norm(t+1) = norm(pp);
end

// Plots
figure(1);
clf;
plot(x_truth(1,:),'ko');
plot(y,'r+');
plot(xx(1,:),'g.');

figure(2);
clf;
plot(x_truth(1,:) - xx(1,:));


How to change a continuous signal to form a discrete signal?

When your continuous-time signal $$x(t)= 2\cos(200\pi t)+8\cos(250\pi t)$$ is sampled at a rate of $F_s = 400$ Hz, (or sampled with a period of $T_s= \frac{1}{F_s} = \frac{1}{400}$ seconds per sample), you can use the following to find the mathematical representation of the discrete-time signal $x[n]$ :
$$ x[n] = x_c(t)|_{t=nT_s} $$
Applying this to $x(t)$ leads to:
$$x[n] = x(nT_s)= 2\cos(200\pi n T_s)+8\cos(250\pi n T_s)$$
$$x[n] = x(n/400)= 2\cos(200 \pi n /400)+8\cos(250\pi n /400)$$
$$x[n] = 2 \cos(\frac {2\pi}{4} n) + 8 \cos(\frac {5\pi}{8} n)$$
Now observing the pattern of sinusoidal signals in discrete time as $\cos(\omega_0 n)$ where $\omega_0$ is the discrete time radian frequency between $-\pi$ and $\pi$ (or $0$ to $2\pi$ if you prefer so) you can deduce your discrete time frequencies as $\pi/2$ and $5\pi/8$
Note that $\omega_0 = \pi$ will be the highest frequency discrete time sinusoid.


Implementing an analog approximation of an FIR filter given its impulse response

FIR filters typically have a sinc^x type impulse response spectrum or something that looks like it, in some but not all cases. The sinc^x shape is mainly seen when the impulse response h(n), if plotted on cartesian coordinates, look like an even function time shifted where the coefficients are placed symmetrical about the origin. This is a simplistic answer based on a particular type of FIR filter such as a Type I or II. You have not given enough information to give a more specific answer.

So your analog filter would have to look like a sinc^x or some type of comb filter with lobes falling off faster than a typical comb filter. But again this is a simplified answer.

In short it is theoretically impossible to do what you want exactly, or even approximately, but I can't tell you that you cannot somehow come up with a crude approximation, say using delay lines


What is special about the frequency $\omega_0=\pi$ that suddenly causes rate of oscillation decrease?

As a consequence of the sampling theorem, the frequency spectrum of a discrete random variable is periodic, but it turns out that signals that are inherently discrete in time and not derived from an A/D converted continuous time signal are also periodic in frequency.  As a consequence, as the frequency is increased, it goes from the positive portion of one periodic interval, to the next, i.e. it transitions to the negative image of the next interval.

In the old days, DSP was mostly a sub discipline of Electrical Engineering, and DSP was focused on doing things digitally that were previously done analog.  In order to talk about DSP, without having to reference, a sampling frequency, a convention called unit sampling was commonly used, and a period of the spectrum is from -1/2 to 1/2.

As time progressed it was noted that DSP wasn't directly tied to processing analog derived signals and unit sampling was somewhat artificial, so it was recognized that $2\pi$ was a natural periodic interval, an $2\pi$ fits between $-\pi$ and $\pi$, which is why $\omega=\pi $ is a transition point in the spectrum.


How to demonstrate the orthogonality of OFDM subcarriers?

Orthogonality is defined as "the inner product of two vectors equals zero".
Now, in OFDM, the transmit vector for a single subcarrier is exactly one row vector $\mathbf D_k$ of the DFT Matrix $\mathbf D$, multiplied by the complex value of a symbol $c_k$, i.e. $c_K \mathbf D_K$.
Two different subcarriers $k, l, k\ne l$ hence have the inner product $\langle c_k\mathbf D_k,c_l\mathbf D_l\rangle$; inner products are linear things, hence that's
\begin{align}
\langle c_k\mathbf D_k,c_l\mathbf D_l\rangle &= c_kc_l \langle \mathbf D_k,\mathbf D_l\rangle\\
&= c_kc_l \begin{cases}0&k\ne l\\\|\mathbf D_k\| & k = l\end{cases} &\text{q.e.d.},
\end{align}
because the DFT matrix is unitary.


Frequencies of complex exponentials in Discrete Fourier Transform

It's not necessary for a discrete-time signal to have the same period as the continuous-time signal that it represents. As long as the continuous-time signal is band-limited and as long as the sampling frequency is greater than twice the bandwidth of the continuous-time signal, the discrete-time signal can perfectly represent the continuous-time signal.
In your example, the continuous-time signal is
$$x(t)=\cos(2\pi f_0t)\tag{1}$$
and the discrete-time signal is
$$x_d[n]=x(nT)=\cos(2\pi f_0Tn)\tag{2}$$
$x_d[n]$ perfectly represents $x(t)$ if $f_0T<\frac12$ is satisfied. In general, $x_d[n]$ is not even a periodic sequence. It is only periodic if $f_0T$ is rational.
With the choice $f_0T=3/16$ we make sure that $x_d[n]$ perfectly represents $x(t)$, and we also make sure that $x_d[n]$ is periodic with period $16$. Note that there can't be $3$ periods of the discrete-time signal in the DFT frame, simply because $16/3$ is not an integer.
Since $x_d[n]$ has period $16$, taking a length $16$ DFT results in a single frequency component. With $f_0T=3/16$, that component occurs at frequency index $k=\pm 3$, which exactly represents a continuous-time sinusoid at frequency $3/16\cdot f_s$ (with $f_s=1/T)$.


Efficient frequency offset estimation of known signals

Time aligned
If the signals are time aligned, you can conjugate-multiply the received signal with the reference signal and look for the dominant frequency.
Max-Likelihood estimation of the frequency shift and complex gain of a known signal in Gaussian Noise
Let $\vec r$ be the noisy observation of a frequency translated version of a sent $\vec x$ with unknown gain.
Here, let us consider the channel $\mathbf A \in \Gamma$
$$
\vec r = \mathbf A \vec x + \vec w  \in \mathbf C^N \\
[\mathbf A]_{nn} = \rho \exp(\omega n) , \rho \in \mathbf C^N , \omega \in \mathbf R \\
[\mathbf A]_{mn} = 0, \forall m \ne n  \\
p(\vec w|\tau) \propto \exp\left( \frac {-\| \vec w \|^2} {\tau}  \right)
\\
\mathbf A^* \triangleq
\operatorname*{argmax}_{\mathbf A}
\max_{\tau}
\underbrace{p(\vec r| \mathbf A,\vec x ,\tau)}_{  \propto \exp \left( \frac {-\| \mathbf A \vec x - \vec r \|^2} {\tau}  \right) }
\\
\mathbf A^* =
\operatorname*{argmax}_{\mathbf A}
\underbrace{
-\| \mathbf A \vec x - \vec r \|^2
}_{
\triangleq f(\mathbf A)
}
\\
f(\mathbf A)
=
- (\mathbf A \vec x - \vec r )^H
(\mathbf A \vec x - \vec r )
\\
= - \vec x^H
\underbrace {\mathbf A^H\mathbf A }_{|\rho|^2 \mathbf I
}
\vec x
+ 2\mathop{re}
(\mathbf A \vec x)^H \vec r
-
\underbrace
{\vec r^H \vec r}_{
\|\vec r\|^2}
$$
Notice the first and third terms have no dependence on frequency $\omega$ nor the phase of $\rho$.  So we simply need to find the frequency that best aligns $\vec x,\vec r$.  That will imply the phase.
$$
\omega^* = \operatorname*{argmax} \sum_n |\exp(-j\omega n) \bar {x_n} r_n|^2
$$
One can use an FFT or some such frequency estimator to recover the dominant frequency
Not Time-Aligned
If the signals are not time-aligned, then you've got a time and frequency search problem on your hand.  With low SNR, this is essentially the problem solved by the cross-ambiguity function (CAF).  If you've got several dB of margin you can take a shortcut.
If you have a reference signal with some AM variation; you can time-align first by AM detecting both signals and then correlating your input by the reference AM.  If your  reference signal doesn't have much amplitude variation, but its spectrum is "spikey" you might be able to correlate the magnitude spectra of the received and reference.  This is the same trick applied in the frequency domain, rather than time domain.
With the exception of the CAF, these are all fairly cheap operations.


Question About Linear and Circular Convolution - 1D and 2D

Let me answer you:

For a signal of size $ m $ and a filter of size $ n $ the output of Linear Convolution is $ n + m - 1 $. In case of 2D signal of size $ \left( m, n \right) $ and filter of size $ \left( p, q \right) $ the output size is $ \left( m + p - 1, n + q - 1 \right) $.
You can read about Circular Convolution in Wikipedia. Basically when a convolution is applied on finite discrete signals one should take care of the boundaries. In most cases the default is assuming the signal i padded with zeros which results in Linear Convolution. If you use padding which build a periodic / circular signal and then apply convolution you will get Circular Convolution. It turns out that frequency domain multiplication of discrete signals is equivelnt of Circular Convolution in spatial domain.
You need to pad it with zeros and line the axis origin to match the image.
Have a look at my answer for Kernel Convolution in Frequency Domain - Cyclic Padding. I also shared a MALAB code which shows how to build the kernel correctly.


How to extend hough lines over the image space

Fit a line through the end points using standard line model $ax+by+c=0$. Then draw them on the image from -10000, 10000. The points remaining on the image will give you what you are asking for.


How to find inverse of convolution integral?

The integral of a signal $x(t)$

$$f(t)=\int_{-\infty}^tx(\tau)d\tau$$

can be written as the convolution of $x(t)$ with the unit step function $u(t)$:

$$f(t)=(x\star u)(t)$$

Consequently, the integral of a convolution

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

can be written as

$$g(t)=((x\star y)\star u)(t)$$

And due to associativity of convolution we have

$$g(t)=((x\star u)\star y)(t)=((y\star u)\star x)(t)$$

So options $(c)$ and $(d)$ in your question are both correct.


Note that the squared magnitude of the frequency response is given by
$$\big|H(j\omega)\big|^2=\frac{1}{1+\epsilon^2T^2_N\left(\frac{\omega}{\omega_c}\right)}\tag{1}$$
In the $s$-domain we have
$$H(s)H(-s)=\frac{1}{1+\epsilon^2T^2_N\left(\frac{s}{j\omega_c}\right)}\tag{2}$$
Computing the zeros of $(2)$ does not only result in the zeros of the filter's transfer function $H(s)$, but also in the zeros of $H(-s)$. The zeros of $(2)$ lie on an ellipse, and for each zero $s_0$ in the left half-plane, there's also a zero $-s_0$ in the right half-plane. Since we want a stable filter, we assign all zeros in the left half-plane to $H(s)$.


You can run the signal through a low pass filter to get rid of the high frequency content. After that, if you know the expected length of the flatness of the signal, you can perform matched filtering with a template of that length. If the cross-correlation result exceeds some threshold, consider the signal detected.


Is there an example of an infinite signal?

After the proposed heat death of the universe, no active life form will be available to measure the rest of any given signal. That's not a real world to anyone.

One could define a signal not as function of time but as function of location, and wrap around the location. An example is the present ambient temperature as function of distance traveled to the east from Pontianak, Indonesia, on the Equator. We don't know if the universe is infinitely large so measurements on a straight line would be a bit iffy.


The frequency response of a single FFT bin filter looks like a Sinc function, which has a massive amount of overshoot or ripple at frequencies between FFT bins.  So your filter is only useful if you can strictly guarantee that the input to the FFT only contains unmodulated frequencies that are strictly and exactly periodic in the FFT aperture length (e.g. not a window on longer data of unknown periodicity or a non-periodic stream).

A notch filter is equivalent to adding the inverse of the equivalent bandpass filter to the signal.  If the bandpass has ripples, they will get added to the result of the notch, thus potentially adding a lot of noise.

So you end up with a not-so-simple (as it requires pairs of FFT computations) that likely produces terrible results.


What is the definition of a "complementary filter"?

Usually, a complementary filter (like a complementary function) complements another filter.  The two filters that are complementary to each other add to one.  Or, at least, add to an all-pass filter (which is what Linkwitz-Riley crossovers do.

so either

$$ H(f) + G(f) = 1 $$

or

$$ H(f) + G(f) = A(f) $$

for $H(f)$ and $G(f)$ being complements of each other and $|A(f)|=1$ is an APF.


Similarity in signal plots but dissimilarity in sensed content

First, the original $x$ and reconstructed  $\hat{x}$ signals have a peak amplitude around $0.35$. Their difference peaks above $0.40$. That happens between two time-shifted like-alike signals.

Second, the difference seems to relate to the amplitude. Processing maybe be non-linear, and it could be interesting to look at relative differences like $2(x -\hat{x})/(x +\hat{x})$ as well.

To enhance understanding, what would  3D plots $(x[n],x[m],x[n]-x[m])$ and $(x[n],x[m],2\frac{x[n]-x[m]}{x[n]+x[m]})$ like like?


What is the reason for filtering out the negative frequencies of a signal?

I'll answer your questions out of order:

3. Negative frequencies do exist. There is nothing controversial about them.

For a cosine signal: since $\cos(2\pi(-f)t)=\cos(2\pi ft)$, a cosine of negative frequency $-f$ is equal to a cosine of frequency $f$.
Since $\sin(2\pi(-f)t)=-\sin(2\pi ft)$, a sine of frequency $-f$ is $\pi$ radians out of phase with respect to a sine of frequency $f$.
More importantly, a complex exponential $e^{2\pi ft}$ can be represented by a point in the complex plane, which rotates counter-clockwise if $f$ is positive and clockwise if $f$ is negative.

5. It turns out that, for all physical signals (whose imaginary part is zero), the magnitude spectrum is even (the negative frequencies are a mirror image of the positive frequencies). There is no need to display or calculate them. That's the reason a spectrum analyzer will only display positive frequencies. Complex signals do not have an even magnitude spectrum, and you need to calculate it for negative as well as positive frequencies.

2. Filtering out the negative frequencies is just that: remove the negative frequencies of a signal. The result is a complex signal, because the resulting spectrum is not even. The Hilbert transform you mention at the start of your question is an easy way to implement such a filter. Without it, you'd need a complex filter, which are not trivial to implement and use.

1. There can be many reasons to filter out the negative frequencies of a signal. In digital communications, it is used as one step in obtaining the complex envelope of a modulated signal. In general, a modulated signal looks like $s(t)=\Re[a(t)e^{2\pi f_ct}]$, where $a(t)$ is a quadrature signal, and is complex, and where $f_c$ is very large (up to several gigahertz). However, the information is contained in $a(t)$, so we would like to recover it from $s(t)$. This is accomplished by:

Filtering out the negative frequencies of $s(t)$
Downconverting the result to baseband


Difference between single channel and single microphone separation

It really depends on what you mean by a channel. Normally a channel in a source separation scenario would the number of sensors detecting your unknown/partially known signal sources. In the audio domain sensors are referred to as microphones.

If this viewpoint is taken single channel separation and single microphone separation are equivalent. However depending on the domain the processing behind these algorithms is often different.

For the one channel (sensor) case in any signal domain, the assumptions of independent component analysis (ICA) do not hold (there are more signal sources than sensors and the number of sensors is not greater than one). Computation auditory scene analysis (CASA) can be used to address this case. Though note that the other answer is incorrect in saying that this is for only one signal source. CASA can be expanded to separate more than one source.


In general LTI System is invertible if it has neither zeros nor poles in the Fourier Domain (Its spectrum).
The way to prove it is to calculate the Fourier Transform of its Impulse Response.

The intuition is simple, if it has no zeros in the frequency domain one could calculate its inverse (Element wise inverse) in the frequency domain.

Few remarks for the practice world:

If for a certain system we know the input signals are limited to certain band (Or zones) in frequency domain it is enough the LTI system has no zeros in this zone only.
In practice we have noise hence even values which are not mathematically zero but close to 0 (Or even low relative to the gain of the LTI system, Say ~-60 [dB]) are making the system not invertible.

Look for many question in this community about Deconvolution and Inverse Problems.


Why is bandwidth always limited in a real (physical) channel?

Information bandwidth is dependent on signal to noise ratios.  At absolute zero, quantum level signal quantization and quantum noise will limit the lower bound on the noise floor.  At higher temperatures, thermal noise creates a higher noise floor in any information receiving equipment.  There may also be an upper limit on power density in a signal before its mass equivalent creates a black hole, which may be problematical in terms of information conservation.  In terms of transmitters, at some point you will melt the wires.

But in terms of pure math, there is no (aleph zero) finite limit on the amount of information that a single real number (one sample) can contain.


Kalman Filter | Why Shift the Belief Distribution

The idea is that the robot is trying to figure out where it is.  The only information it can use is:

Its door sensor which says either you're next to a door or you're not next to a door.
Its map of where the doors are.

That's how we go from the "maximum confusion" (uniform) plot to the one with three humps. This uses the first you're next to a door measurement.

Then, the robot knows it moves to the right.  Therefore it must update its probability with a similar movement to the right.  However, its motion sensors are a little uncertain, so this has the effect of smoothing the three bumps (well, the whole distribution) a little bit. So the three bumps move to the right and they spread out a little.

The peaks of the three spread out humps are still the robot's best guess as to where it is on its map.

Now we get the second you're next to a door measurement.  We multiply our prior (the slightly spread out three-humped distribution) with a new three humped distribution (from our known map).  The prior second-door location and the map second-door location will align, so that will generate the largest peak.


How to get the orientation value of a rotated square?

Disclaimer:  I am not a MATLAB user, so I had to search for this and have no way to test it.

You should be able to find two opposite corner points by using the "Extrema" property from "regionprops".  The tangent of an angle is a slope, so calculate the slope and take the inverse tangent.  Do note that the line may be vertical, i.e. undefined slope, and you have to deal with that case.  If you choose the left-top and right-bottom, you should get a more horizontal line and not have to worry about it.

You also have to keep radians and degrees straight.

As math:

$$ \Delta y = y_{LT} - y_{RB} $$
$$ \Delta x = x_{LT} - x_{RB} $$

$$ m = \frac{\Delta y}{\Delta x} $$

$$ \theta = tan_{deg}^{-1} (m) $$

$$ d = \theta - 45 $$

By finding the leftmost and rightmost $x$'s, there shouldn't be a problem with the slope calculation.  The line between the two extreme points should be a diagonal of the square.  A diagonal should be at 45 degrees, so $d$ is how much you need to rotate the figure to get it there.


How we can represent a bandpass signal equation in matlab

Usually (almost always), simulation is done using the discrete time baseband signal. Lets say you want to simulate 32 BPSK symbols (either +1 or -1). The MATLAB code for this is:

 symbols = sign( rand(32, 1) - 0.5 );

Now that you have the symbols, the next step is to upsample. You can either use the MATLAB function or do it yourself:

 symbols_up = upsample(symbols, samplesPerSymbol);

Then create and use a pulse shaping filter (just use rcosdesign):

 pulseFilter = rcosdesign(rolloff, span, samplesPerSymbol, 'sqrt');
  txSig = filter(pulseFilter, 1, symbols_up);
  figure, plot(txSig)


How to match zero-pole diagrams to their frequency responses (Discrete Time)

These are all frequency selective filters. Note the zeros on the unit circle, which correspond to zeros of the frequency response. The filter's passband is determined by range of the angles of the poles, and the stopband is determined by the range of the angles of the zeros.

As an example, let's consider PZ-diagram $\#1$: since the stopband is around DC (i.e., angle zero), and there are no more zeros, just poles, at angles greater than approximately $\pi/3$ (and, by symmetry, smaller than approximately $-\pi/3$), this must be a high pass filter. The cut-off frequency must be close to $\pi/3$, so all frequency responses can be eliminated except for $\#4$ and $\#5$. However, frequency response $\#5$ has a zero at DC, which is not present in the PZ-diagram. Consequently, we can match PZ-diagram $\#1$ to frequency response $\#4$. I'm sure you can match the other diagrams yourself.


How is a signal multiplied with the window functions?

Matlab window functions such as hamming(N) returns the $N$ point Hamming window sequence $w[n]$ in a $N \times 1$ matrix (a column vector).

Therefore if your sequence $x[n]$ is represented in a row vector of size $1 \times N$ then for obtaining the sample-by-sample product $x[n]w[n]$, you should either transpose $w[n]$ or the other $x[n]$, so that their dimensions permit the correct sample-by-sample multiplication in between and produce the output sequence $v[n]=w[n]x[n]$ represented in a matrix of the size $N \times 1$ or $1 \times N$ depending on which one you have transposed.

v0 = x.*w      % x: Nx1, w: Nx1, v0: Nx1
v1 = x .* w'   % x: 1xN, w: Nx1, transpose w => v1: 1xN
v3 = x' .* w   % x: 1xN, w: Nx1, transpose x => v3: Nx1

Note that if the transposed signal is complex valued, then you should use the conjugate transpose instead :

v1 = conj(x').* w  % x[n] was complex valued...

Note also that the sample by sample product is perfomed by the $.*$ operation under Matlab;

  v = x .* w    % is a sample by sample product implementing v[n] = x[n]w[n]

If you omit the period and use the $*$ alone, you get matrix product or the dot product in case of vectors involved.

For this case when you have two vectors, $x$ of size $N \times 1$ and $w$ of size $1 \times N$, then their matrix product would produce a matrix of size $N \times N$

v = x * w    % x: N x 1, w: 1 x N => v: N x N


What is the transfer function $H_{xy}$ of this flow graph?

If you mark with $v[n]$ the signal between summing nodes for input $x[n]$ and output $y[n]$, then:

$v[n] = x[n] -8v[n-2] -4v[n-3]$

Re-arranging sides:

$v[n]+8v[n-2]+4v[n-3]=x[n] $

Using the Z-transform:

$V(z)[1+8z^{-2}+4z^{-3}]=X(z)$

Now let's calculate the transfer function between $v[n]$ and $x[n]$:

  $H_{x\to v}(z)=\dfrac{V(z)}{X(z)}=\dfrac{1}{1+8z^{-2}+4z^{-3}}$

Now moving to the relation between output $y[n]$ and $v[n]$:

$y[n]=v[n]+2v[n-1]$

Again calculating the Z-transform:

$Y(z)=V(z)[1+1+2z^{-1}]$

Now converting the transfer function $H_{vy}(z)$:

  $H_{v\to y}(z)=1+2z^{-1}$

Final step is to merge these two together. We know that in time domain it is a convolution: $h_{x\to y}[n]=h_{x\to v}[n]\star h_{v\to y}[n]$, thus in Z-transform domain it will be the multiplication, yielding the final result:

  $H_{x\to y}(z)=H_{x\to v}(z)H_{v\to y}(z)=\dfrac{1+2z^{-1}}{1+8z^{-2}+4z^{-3}}$


Spectrum in real life with Spectrum in theory

why does this power measurement (which is an instant function of time)

It's NOT a instant function of time.
For example instantaneous power of a sine wave is
$$ p(t) = x(t)^2 = A \cdot cos^2(t) = \frac{A^2}{2}(1+cos(2t))$$
You have to apply some amount of time averaging to get rid of the $cos(2t)$ component and that's what all spectral analyzer do.

: if I stop the input signal, the spectrum becomes 0) ...

Not instantaneously. Since there is time domain averaging involved the power measurement can only change at a finite rate.
In general there is an "uncertainty principle" at work here: the higher the spectral resolution, the slower the system can react or track. For your first implementation, that is determined by the bandwidth of your demodulation filter. The smaller the bandwidth, the slower the output will change.
If you use a Fourier Transform based analyzer, this will be determined by the length of the (finite) time domain window you apply the transform over.

why does this power measurement ... coincide with the theoretical power spectral density |X(f)|2

It doesn't. All physical measurements are approximate in nature. They all have a margin of error, tolerance, repeatability variance etc. You can't measure anything exactly, but only "good enough" for your specific application".
You also need to define what exactly you want to measure within the context of your specific application and requirements. You don't need the exact power spectrum of an infinitely long signal, since infinitely long signals do not exist in nature.


Why is not Fourier Transform Good for Non-linear Processes

Because complex exponentials $e^{\jmath \omega t}$, which are results of Fourier transform, are the eigenfunctions for linear, time invariant (LTI) systems. See eigenfunction of LTI. Also see this answer on SP.SE.

Thus, Fourier transform is useful for analyzing linear (not suitable for non-linear one) time invariant (can be intepreted as stationnary) system.


The use of FDM for a single message signal

You can simply send the message multiple ($N$) times and if all N messages were received at the same SNR you would coherently average the messages for a processing gain in SNR equal to $10\log_{10}(N)$ in dB.  This is effectively trading bandwidth for SNR as you are using more resources to send the same message. To coherently add you remove the complex carrier phase for each message prior to adding in the average. If the messages were not received at the same SNR (such as if it was a fading channel) you would optimally weight each message by the SNR of each message prior to averaging.
This latter point is similar to what occurs in a matched filter in that each sample within a symbol duration is optimally weighted by the SNR for that sample prior to averaging over the symbol duration.

See @MBaz's good comments under the OP's question clarifying that there really is no actual SNR gain if you consider the total signal power of all messages sent, since the total signal power would need to increased to realize the gain listed above.


Why is my time domain interpolation via zero-padding in frequency domain wrong?

The answer above is correct. Just to clarify a bit further, using x = np.linspace(0,10,5) will produce 5 numbers from 0 to 10 inclusively
np.linspace(0,10,5)
array([ 0. ,  2.5,  5. ,  7.5, 10. ])

You don't want the last number because in your example the last number is the first number of the next period. A correct implementation would be:
periods = 4.*np.pi
x    = np.arange(0., periods, periods/n)
xpad = np.arange(0., periods, periods/npad)


Lightweight audio listener for detecting blip

One lighter-weight approach might be to use a narrow-band IIR filter centered at your blip frequency, rectify or square, feed this result to a low-pass filter, and look for a value passing some threshold.  A Goertzel filter is another possibility.

Note that the most current iOS Accelerate framework includes a reportedly power-efficient biquad IIR filtering function.  Also note that on recent iOS devices, even a more heavyweight approach (FFTs, autocorrelation, and such) on a single channel of real-time audio input might use only a tiny percentage of total device power.


Fourier coefficients of discrete difference of a square wave

It's quite straightforward to express the Fourier coefficients of a shifted version of a function in terms of the Fourier coefficients of the original function.
If the Fourier series of a $T$-periodic function $f(t)$ is given by
$$f(t)=\sum_{k=-\infty}^{\infty}c_ke^{j2\pi kt/T}\tag{1}$$
then it follows that the Fourier series of a shifted version must be
$$\begin{align}f(t-t_0)&=\sum_{k=-\infty}^{\infty}c_ke^{j2\pi k(t-t_0)/T}\\&=\sum_{k=-\infty}^{\infty}c_ke^{-j2\pi kt_0/T}e^{j2\pi kt/T}\\&=\sum_{k=-\infty}^{\infty}\tilde{c}_ke^{j2\pi kt/T}\tag{2}\end{align}$$
Hence, the Fourier coefficients of $f(t-t_0)$ are given by
$$\tilde{c}_k=c_ke^{-j2\pi kt_0/T}\tag{3}$$
And, by linearity, the Fourier series coefficients of $g(t)=f(t)-f(t-t_0)$ are
$$d_k=c_k-\tilde{c}_k=c_k\big(1-e^{-j2\pi kt_0/T}\big)\tag{4}$$

For periodic discrete-time sequences we can use the Discrete Fourier Tansform (DFT), which is just a discrete Fourier series representation:
$$x[n]=\frac{1}{N}\sum_{k=0}^{N-1}X[k]e^{j2\pi nk/N}\tag{5}$$
with
$$X[k]=\sum_{k=0}^{N-1}x[n]e^{-j2\pi nk/N}\tag{6}$$
In a completely analogous manner to the continuous-time case, we can derive an expression for the DFT coefficients of the sequence $y[n]=x[n]-x[n-n_0]$, $n_0\in\mathbb{Z}$:
$$Y[k]=X[k]\big(1-e^{-j2\pi n_0k/N}\big)\tag{7}$$


How do multirate systems introduce 'efficiency' in computation?

I will try to explain from my understanding.

For interpolation, the whole process can be viewed as two parts: up-samping with padding zeros, and the low-pass filtering. The zero padding will increase the length and thus the computational costs. For the low-pass filtering, there are lots of zeros, thus these multiplications could be omitted.

For decimation, the first part is anti-aliasing low-pass filtering, and the 2nd is the down-sampling. Since only part of the low-pass filtering output is needed, it is apparent that the filtering for other points can be omitted.


Steady state variance of a stochastic differential equation - relation between the frequency and time domains

BOTTOM LINE UP FRONT: I think the exponential decay growth in $\left<|x(t)|^2\right>$ can be shown in the frequency domain only if the "boundary terms" are nonzero when we compute the Fourier transform of $dx(t)$ from the original SDE.

I provide only a start in the work below.

Since these processes seem like they could possibly be complex-valued, I will consider $\left<x(t)\overline{x(t')}\right>$ and $\left<y(t)\overline{y(t')}\right>$, where the $\texttt{overline}$ indicates complex conjugation.

I will write out some steps, because I need to see some things that you did not write explicitly.

\begin{equation}
dx(t) = ax(t)dt + by(t)dt.
\end{equation}
Fourier transform:
\begin{equation}
\begin{split}
\int e^{-i\omega t}dx(t) &=~ \int e^{-i\omega t}ax(t)dt + \int e^{-i\omega t} by(t)dt\\
\int\left[d\left(e^{-i\omega t}x(t)\right) - x(t)d(e^{-i\omega t})\right] &=~ a\widehat{x}(\omega) + b\widehat{y}(\omega)\\
\underbrace{\int d\left(e^{-i\omega t}x(t)\right)}_{\textrm{Boundary terms: $-c$, a constant}} - \int x(t)(-i\omega)e^{-i\omega t}dt &=~ a\widehat{x}(\omega) + b\widehat{y}(\omega)\\
i\omega\widehat{x}(\omega) &=~ a\widehat{x}(\omega) + b\widehat{y}(\omega) + c
\end{split}
\end{equation}
Result:
\begin{equation}
\widehat{x}(\omega) = \frac{b\widehat{y}(\omega) + c}{i\omega - a}
\end{equation}

\begin{equation}
\begin{split}
x(t) &=~ \frac{1}{2\pi}\int \widehat{x}(\omega) e^{it\omega}d\omega\\
&=~ \frac{1}{2\pi}\int \frac{b\widehat{y}(\omega) + c}{i\omega - a} e^{it\omega}d\omega\\
& \\
x(t') &=~ \frac{1}{2\pi}\int \widehat{x}(\omega') e^{it'\omega'}d\omega'\\
&=~ \frac{1}{2\pi}\int \frac{b\widehat{y}(\omega') + c}{i\omega' - a} e^{it'\omega'}d\omega'
\end{split}
\end{equation}
\begin{equation}
\end{equation}
\begin{equation}
\begin{split}
\left<x(t)\overline{x(t')}\right> &=~ \left<\frac{1}{2\pi}\int \frac{b\widehat{y}(\omega) + c}{i\omega - a} e^{it\omega}d\omega\frac{1}{2\pi}\int \frac{\overline{b\widehat{y}(\omega') + c}}{-i\omega' - a} e^{-it'\omega'}d\omega'\right>\\
&=~ \frac{1}{4\pi^2}\int\int\frac{b^2\left<\widehat{y}(\omega)\overline{\widehat{y}(\omega')}\right> + bc\left<\overline{\widehat{y}(\omega')}\right> + b\overline{c}\left<\widehat{y}(\omega)\right> + |c|^2}{(i\omega - a)(-i\omega' - a)}e^{it\omega}e^{-it'\omega'}d\omega d\omega'
\end{split}
\end{equation}

\begin{equation}
\begin{split}
\left<\widehat{y}(\omega)\overline{\widehat{y}(\omega')}\right>
&=~
\left<\int y(\tau)e^{-i\omega\tau}d\tau\overline{\int y(\tau')e^{-i\omega'\tau'}d\tau'}\right>\\
&=~
\int\int\left<y(\tau)\overline{y(\tau')}\right>e^{-i\omega\tau}e^{i\omega'\tau'}d\tau d\tau'\\
&=~ \int\int\delta(\tau-\tau')e^{-i\omega\tau}e^{i\omega'\tau'}d\tau d\tau'\\
&=~ \int e^{-i(\omega - \omega')\tau}d\tau,
\end{split}
\end{equation}
after integrating in $\tau'$. This integral does not converge, but in the sense of distributions,
\begin{equation}
\int e^{-i(\omega - \omega')\tau}d\tau = 2\pi\delta(\omega - \omega').
\end{equation}

What to do about the cross-terms? We must consider the integrals in $\omega$ and $\omega'$ separately for those. Consider the one for $\omega$:
\begin{equation}
\int\frac{\left<\widehat{y}(\omega)\right>}{i\omega - a}e^{i t\omega}d\omega
\end{equation}
If we assume analyticity of $\left<\widehat{y}(\omega)\right>$ (as a function of complex-valued $\omega$) and assume that it decays rapidly enough for large $|\omega|$ in the complex-$\omega$ plane, then we can appeal to contour integration and Cauchy's integral formula. The only pole in the complex-$\omega$ plane is found at  $\omega = -ia$, where I assume $a$ is real. The integral is
\begin{equation}
\int\frac{\left<\widehat{y}(\omega)\right>}{i\omega - a}e^{i t\omega}d\omega
~=~
\lim_{N\to\infty}\oint_{\gamma_N}\frac{\left<\widehat{y}(\omega)\right>}{i\omega - a}e^{i t\omega}d\omega
~=~ 2\pi i\left<\widehat{y}(-ia)\right>e^{at},
\end{equation}
where $\gamma_N$ is a path that includes the real interval $[-N,N]$ and a semi-circular arc connecting $N$ and $-N$. You can find examples of this path in almost any undergraduate complex analysis book. We assume that the integral on this arc decays to zero as $N\to\infty$.

I think this is how we pick up the exponentially decaying growing behavior in $\left<|x(t)|^2\right>$.


Dominating frequency in a "clear" 50 Hz signal

The size of $y$ is 201. So frequency bin of fft is not 1 Hz. Change code for calculating $t$:

Fs=200.
F=50.
t = [i*1./Fs for i in range(200)]

Now you have only 1 peak on 50 Hz.

P.S.

Plot magnitude spectrum
Try to modify F (49, 49.5, 49.67). Enjoy!

Welcome to DSP and new questions!


MATLAB's documentation doesn't recommend the use of firrcos, instead use rcosdesign as recommend. Type doc firrcos  in your command window, here are the first lines:

firrcos Raised Cosine FIR Filter design.
WARNING: firrcos is not recommended. Use RCOSDESIGN instead.

My guess is that firrcos is actually discontinued in new releases. To better help you, try using rcosdesign as recommended. You can actually get the actual MATLAB code for this function by typing the following lines in your MATLAB's command window: edit rcosdesign (you can also check edit firrcos). This opens a new editor window with the actual code for the function. This can help you.
Note that in your EDIT1, Fc is actually the cutoff (not carrier) frequency if you double-check the documentation, and that with three input arguments as you have the default $F_s = 2$. The roll-off factor is a measure of the excess bandwidth expressed as a fraction. The span of the filter in symbols is for you to define.
To add to what has been mentioned by @MBaz, what you need at the transmitter is the root-raised cosine filter. This is so that at the receiver side you also put a  root-raised cosine and the two together result in a raised-cosine filter which is a Nyquist filter with zero ISI while the root-raised cosine doesn't exhibit zeros (in amplitude) at symbol times.
EDIT:
To see the equivalence between the two functions you'll have to give the correct input parameters to both functions. As discussed, first things is having a root-raised cosine. You can do this by setting the 6th and 4th input of firrcos and rcosdesign respectively to sqrt. The equivalents of the remaining parameters are as follows:
beta = R;
sps = Fs/(2*Fc);
span = N / sps;

Then you have the unnormalized two filters as:
h1 = firrcos(N, Fc, R, Fs, 'rolloff', 'sqrt');
h2 = rcosdesign(R, N/(Fs/Fc/2), Fs/Fc/2, 'sqrt');

More on this together with the needed normalization can be found in MATLAB's Communications Toolbox™ Release Notes for R2020b under the chapter on R2013b (I guess it got discontinued after then) titled Functionality Being Changed or Removed.
For your Simulink, click on Simulink Library Browser then on the list on your left choose Communications Toolbox then pick Comm Filters. There you have all your filters, you can set the chosen raised-cosine filter as discussed up.


Pre-echo on very short and percussive sounds - for example crisp solo drums, or better, electronic drums. The TR-808 rimshot and clave are very short and impulsive sounds, about 5ms long (for the rimshot), and at less than 160kbps their attacks get "swallowed". At less than 96kbps some notes appear to be duplicated or out of place rhythmically because of pre-echo. Worse, each note in a repeated pattern will sound different depending on how it falls within the analysis window. Here is an example.

Try solo vocals too. I think casual listeners are better at detecting anomalies in voices (a kind of auditory "uncanny valley" effect). A classical example used in the audio coding community is the intro from Suzanne Vega's "Tom's Diner".


Why having a ridge like Ambiguity Function implies lack of Doppler resolution?

If you examine a contour plot of the Ambiguity function - the lines connect equal values of the Ambiguity function. This means that if you travel along one of these lines (i.e. keeping the same value of the function) then your system will not be able to distinguish between the signal that have this combination of Doppler and time delay. Note that the values for time/Doppler resolution depend on what value of the Ambiguity function you choose i.e. which contour line you decide on. I'm saying just pick one for the purposes of illustration.

If you look the Ambiguity function for a CW pulse of 2 secs duration:

You see that if you choose the lowest contour value, that we would be unable to distinguish between time delay of -1.5 secs and 1.5 sec, and our Doppler resolution ranges from -0.5 Hz to 0.5 Hz. If you were use a pulse of longer duration the Ambiguity function would be wider (left to right) and shorter (top to bottom). This means your time resolution would be worse, but your Doppler resolution would be better. The opposite is true if you used a shorter pulse.

If we examine the Ambiguity function for a LFM chirp pulse:

Now we see that the ambiguity exists for a combination of Dopplers and time delays (chose one of the contours that is essentially an ellipse).Previously, for the CW pulse, across a range of time delays the ambiguity was mostly constant for the same Doppler value. The advantage of the LFM pulse over a CW pulse, is that I can use the same pulse length but achieve better time resolution. You could also use a longer pulse and achieve better time resolution. The longer pulse means you can get more energy on the target and improve the signal to noise ratio.


Implementing a phase-neutral band-pass filter using SciPy

How about scipy.signal.filtfilt and scipy.signal.butter filter? See the example at the bottom of the documentation.

Here's another example with and without phase neutrality: Applying filter in scipy.signal: Use lfilter or filtfilt?


Apply "non-causal" filter buffer-wise, a.k.a "soft real-time filtering"

Overlap-add or overlap-save/scrap with zero-padded data are the common methods of using block based convolution on streaming data.  Pad the convolution (FFT/IFFT fast, or linear) by at least the length of the impulse response above your desired noise floor, minus 1.  The basic idea is that these methods save the remainder of the impulse response that doesn't fit in the current block's filtered results, and applies that remainder to blocks that are processed later.  For linear phase FIR filters, this still incurs as least half the length of the impulse response in latency.


Implementation of Overlap Save Filtering with Decimation

Attempt 2
I am not quite sure what is the interpretation of your filtering, in special, the IFFT part, you are taking a slice of the spectrum shifted by 3*Nfft/4, I think it should be (circularly) shifted by 7*Nfft/8.
  % your method
  xa = ifft(X((3*Nfft/4)+1:4*Nfft/4),Nfft/4) / 4;
  % what I would do
  xb = ifft(fftshift(fftshift(X)(3*Nfft/8+1:5*Nfft/8))) / 4;

Notice that you factor is (Nfft/2)^2 because you multiplied the spectra by Nfft.
To make sure I am not missing any constant I plotted the histogram of the correlations between the fft(trimmed(signal)) and decimated(ifft(signal)) for your method and the modified formula (correlation of 1 is exact up to an affine transformation).

The complete code I used to check is this
cora = []
corb = [];
for ii = 1:10000
  x = zeros(Nfft, 1);
  x(1:4:end) = randn(Nfft/4, 2)*[1;1j];
  x = cumsum(x);
  X = fft(x);
  % your method
  xa = ifft(X((3*Nfft/4)+1:4*Nfft/4),Nfft/4) / 4;
  % what I would do
  xb = ifft(fftshift(fftshift(X)(3*Nfft/8+1:5*Nfft/8))) / 4;

  % make sure this is not for a constant factor
  cora(ii) = abs(corr(xa, x(1:4:end))); % 0.66
  corb(ii) = abs(corr(xb, x(1:4:end))); % 0.98
end
figure
subplot(211),hist(cora, 25),xlim([0, 1]);
title('correlation histogram for your method'), xlabel('correlation')
subplot(212),hist(corb),xlim([0,1])
title('correlation histogram for the proposed method'), xlabel('correlation')

Working your example.
A few modifications that makes your frequency and time domain filtering consistent
clear all;close all;clc
fs=1e6;
sim_time=20e-3;
Nfft=256;

t=0:1/fs:sim_time-1/fs;
freq = fs/Nfft*(0:Nfft-1);

%create composite signal
fsig1= 50e3;
fsig2= 100e3;
fsig3= 300e3;
fsig4= 872e3;

sig1=2*cos(2*pi*(fsig1-5e3)*t)+4*cos(2*pi*fsig1*t)+8*cos(2*pi*(fsig1+5e3)*t) + j*(2*sin(2*pi*(fsig1-5e3)*t)+4*sin(2*pi*fsig1*t)+8*sin(2*pi*(fsig1+5e3)*t));
sig2=8*cos(2*pi*(fsig2-5e3)*t)+4*cos(2*pi*fsig2*t)+2*cos(2*pi*(fsig2+5e3)*t) + j*(8*sin(2*pi*(fsig2-5e3)*t)+4*sin(2*pi*fsig2*t)+2*sin(2*pi*(fsig2+5e3)*t));
sig3=2*cos(2*pi*(fsig3-5e3)*t)+8*cos(2*pi*fsig3*t)+4*cos(2*pi*(fsig3+5e3)*t) + j*(2*sin(2*pi*(fsig3-5e3)*t)+8*sin(2*pi*fsig3*t)+4*sin(2*pi*(fsig3+5e3)*t));
sig4=10*cos(2*pi*fsig4*t)+j*10*sin(2*pi*fsig4*t);

sig=sig1+sig2+sig3+sig4;
sig=randn(size(sig)); # white noise for testing
P=41;%filter length
h1 = exp(-((1:P)/10).^2); % (Octave) without designfilt

%create complex band pass by multiplying with complex exponential
fc=800e3;
cmp_exp=exp(j*2*pi*fc*t);
h1_c=h1.*cmp_exp(1:length(h1));
%take FFT of filter impulse response
H1=fft(h1,Nfft);% low pass prototype
H1c=fft(h1_c,Nfft); %complex bandpass filter
%implement overlap save
wLen = Nfft - (P - 1);
BLOCKS = 10;
y_osa = zeros(1, wLen * BLOCKS);
for i=1:BLOCKS % This code snippet compute the OSA filter with decimation by 4
    sig_win = sig(wLen*(i-1)+(1:Nfft));
    X_sig=fft(sig_win,Nfft);
    Y_sig=H1c.*X_sig;%multiply with complex bandpass filter
    y_sigd = ifft(Y_sig);
    y_osa(wLen*(i-1)+1:wLen*i)=y_sigd(P:Nfft);
end

y1 = conv(h1_c, sig(1:(BLOCKS+1)*wLen));
y1 = filter(h1_c, 1, sig);
y1 = y1(P:wLen*BLOCKS+P-1);

figure
hold on
plot(real(y1))
plot(real(y_osa))
```


Filter to speed up convergence of geometric series

The problem is that in the case with positive $c$ you have a DC desired signal ($lim$) and a low frequency noise ($c^n$). In the case when $c$ is negative you have a high frequency noise signal, which is of course much easier to separate from the desired DC value.

If you figure out some way to estimate $c$ (sign and magnitude) you can use a simple FIR filter

b=[1,-c]/(1-c); a=1;

which will exactly eliminate the exponentially decaying noise:

$$y[n]=\frac{1}{1-c}(x[n]-cx[n-1])=\frac{1}{1-c}(lim+c^n-c\cdot lim - c\cdot c^{n-1})=lim$$


What is the ideal type of barcode for paper tickets?

First, about the USB barcode reader: at my previous job (QA at a label printer company), we had one from Scanology. It worked as a keyboard device, the scanned data is sent to the pc as if it were typed. So you don't need any special software.

Code 39 has no error correction, but even a single misread line won't result in wrong data, because it simply won't produce a valid character. So you either get valid data, or you get nothing.
A major advantage of code 39 is that you don't need barcode software, you can simply use a barcode font.
Code 128 is slightly better than code 39 in terms of error detection, but still no fault tolerance. And you can't use a barcode font. In your case, you would use 128A (0-9, A-Z and a few special chars).
If your most critical factor is fault tolerance, then I would suggest QR, with high error correction setting, which gives up to 30% redundancy. This means that up to 30% of your QR code can be unreadable before the scanner gives up. The Wikipedia page on QR shows an example of a torn QR code that is still readable: http://en.wikipedia.org/wiki/QR_code#Error_correction
Another advantage is that it can be scanned in any direction, and it's a quite popular format.
Data Matrix and PDF-417 are similar to QR in terms of fault tolerance, but they are not widely supported for smartphones. I scanned thousands of QR and PDF-417 codes during my time as QA engineer at that label printer company, and I found a lot of high quality QR scanning software for Android and iPhone, but hardly any for PDF-417 and Data Matrix. I assume the situation will be the same (or worse) for Windows Phone.

My advice:

Easy implementation: Code 39. Use a Code 39 barcode font. Works with cheaper barcode scanners.
Very high fault tolerance: QR at level H (High). Use the ZXing library. Handheld USB barcode scanners may cost a tiny bit more.


Difference between Leaking Factor and Forgetting Factor

The play similar role in those algorithms - the ability to forget the past and adapt to current reality.

In the LMS, the classic implementation has $ \alpha = 1 $.
Namely the optimal weights at any point are function of all inputs.

The Leakage factor allows to weigh the past differently in a damped manner which over times means the far past has practically no significance on the current result.
The other side means the ability to integrate data against noise is decreased which is the same balance the user stands when selecting the parameters of the RLS.

Namely, as usual, the balance between high bandwidth filter which adapts fast yet is sensitive to noise or the low bandwidth filter which is slow yet if the data is stationary can handle high energy noise.


According to the link you provided:

$H(error) = H(y-\hat{y}) = H(Var[y]+Var[\hat{y}])$. While $Var[y]$ decreases with increasing of SNR, $Var[\hat{y}]$ stays the same (because it's not dependent of the noise, it depends on the signal alphabet only). So we may guess $H(error)$ decreases with increasing of SNR.


Non invertibility of system $y[n]=x[n]-x[n-1]$ using transform method?

The impulse response of your system is $h[n] = \delta[n] - \delta[n-1]$ and its Z-transform is accordingly

$$H(z)=1-z^{-1}$$

Now, you can evaluate the frequency response of this system:

$$
|H(z)|_{z=\exp(j\omega)}=|1-\exp(j\omega)|
$$

you see, that $H(z)|_{z=\exp(j0)}=0$ (i.e. it vanishes for $\omega=0$). Now, in order for a filter to be invertible, its Z-Transform on the unit circle (i.e. when setting $z=\exp(j\omega)$ must not vanish. If it does not vanish, the inverse is given by $\frac{1}{H(z)}$. If it vanishes, the inverse does not exist.


Finding filaments in high dimensional space

Following my comment, I now have time to post it as an answer.
If you expect that a large enough portion of the vectors belong to that line, you can use RANSAC.
RANSAC works by selecting a small and random subset of the data and fit it to the desired model (say a line). The fitting can be achieved in any way you like, e.g. least-squares. Then, the data is divided to inliers (the points in space that "agree" with the model) and outliers (those who don't "agree"). The process is repeated until you have enough inliers.
This would only work if the probability of selecting the points on the line is high enough.

As far as multiple lines, I don't think RANSAC supports this inherently. You can do it iteratively by removing all inliers from the set once a line has been found.

If you are working in MATLAB\Octave, my favourite RANSAC implementation can be downloaded from here.


Laplace transform of a finite duration signal

The property claimed by Schaum and Oppenheim is also true for the given example. Note that the Laplace transform

$$X(s)=\frac{1-e^{-5(s+2)}}{s+2}\tag{1}$$

has no pole at $s=-2$:

$$\lim_{s\to -2}X(s)=\frac{1-(1-5(s+2))}{s+2}\Big{|}_{s=-2}=5$$

So the ROC is indeed the entire $s$-plane. Even though the ROCs of the two individual signals in your solution have the same right half-plane as ROC, their sum has the entire $s$-plane as ROC because the two signals cancel everywhere except in a finite interval.


Separating signal from noise in similar frequency range

As you have noticed, it is quite unlikely that you can filter your data with a standard frequency filter, thus a linear filter.

Unless the signal and the noise are well separated in time, unless you have a model, you will probably need to resort to other techniques, such as nonlinear filtering, time-frequency or time-scale filtering, or add additional information or assumptions, such as thought signal or noise properties (to name a few: gaussianity, positivity, independence, sparsity), and bend to blind source separation, latent variable analysis, nonlinear matrix separation, etc.


Lets say that $x_1[n]=\delta[n]$ and $x_2[n] = \delta[n-1]$, where $\delta[n]$ is Kronecker delta and that $y_1[n],y_2[n]$ are the output signals for $x_1[n],x_2[n]$.
The downsampling operator discards every other sample. Without the loss of generality, lets assume that we are keeping the even samples and throwing the odd ones. Therefore, after downsampling $x_2[n]$, it is all zeros, which means that $y_2[n] = 0$ for each $n$. Now since you say that $h[n]$ is not all zeros, it is guarantied that $y_1[n]$ will have some non-zeros values.
In other words, $y_1[n]\ne y_2[n]$ for some $n$, and you system is not TI.


Features for audio levels imbalance in films between speeh/music/effects

I would try a sum and difference calculation between the L and R channels to start. Usually when film audio is mixed, the speech is in both channels in phase while the special effects is out of phase.

So if you do:

X = L+R;
Y = L-R;

X should give you the signal with more speech and Y should give you the other effects. Then if you do an RMS calculation on X and Y, you can find how much are they off by.


FIR filter coefficients from frequency responses

The frequency argument for freqz() should be in radians/second, not normalized frequency in [0,1). It has a form that makes it easy for you:

% length of frequency response vector
N = 50;
% cutoff frequency
wc = 0.25;
% design 12th-order lowpass filter
b1 = fir1(12, wc);
% calculate its frequency response at `N` points around unit circle
hd = freqz(b1, 1, N, 'whole');
% inverse transform the frequency response
b2 = ifft(hd);
% `b1` and `b2` should be the same.
plot(1:length(b1), b1, 1:length(b2), b2);

This yields the plot:

You only see one trace because they are right on top of one another, as they should be. b2 is longer in duration because it has the same length as the frequency response vector, which was arbitrarily chosen to have length 50.


My own interpretation:

Sound: a mechanical wave that propagates through the air or water.

Audio: sound in the 20 Hz to 20 kHz range; in other words, sound that is (at least in theory) audible to a large number of humans.

Voice: sound produced by the human vocal tract.

Speech: intelligible voice (i.e. not grunts or screams)

Tone: signal dominated by a single frequency, even if it has harmonics.

Music: Harmonious, pleasant audio produced by man-made musical instruments. Often accompanied by song, which is audio just like music but produced by a human singer's vocal tract.

Musical tone: a tone produced by a musical instrument or singer.

Musical note: one of a particular, finite set of musical tones that, when played in combination, produce a pleasant auditory experience.


Depends on what you are controlling. For DC-motors it is the inertia of the device that acts as a low-pass filter of the PWM modulated signal resulting in a continuous motion. For most LEDs it is the human eyes that do the apparent low-pass filtering.

If the PWM-frequency is not very high you can actually see this by moving your head from left to right while looking at the PWM modulated LED. The light will appear on and off while the LED is moving through your field of vision (this can get very annoying with LED car taillights).


Is the phase difference between recorded sounds from left and right ear the same as the phase difference between their HRTF?

Consider a signal with a particular frequency at phase $\alpha$. Passing this through the two filter paths changes the phase of the frequency to $\alpha + L$ and $\alpha + R,$ where $L$ and $R$ are the values of the phase frequency responses of the two HRTF filters at the particular frequency. Yes, the phase difference $(\alpha + R) - (\alpha + L)$ is equal to $R - L,$ which you call more simply the "phase difference" between the HRTF's.


SNR relation with modulation index and noise power

If $x(t)$ is your message signal, then you can write the AM signal as

$$s(t)=A[1+mx(t)]\cos\omega_0 t$$

where $A>0$ is a real-valued constant, $m$ is the modulation index, and $\omega_0$ is the carrier frequency (in radians). If you assume that the noise added by the channel is white with power spectral density $N_0$, then the SNR after demodulation is

$$SNR=\frac{A^2m^2\overline{x^2(t)}}{4N_0B_x}$$

where $B_x$ is the bandwidth of the message $x(t)$ and $\overline{x^2(t)}$ is its average power.


What Is the Algorithm Behind Photoshop's "Black and White" Adjustment Layer?

Each (color) image is composed of RGB components. when you add (or reduce) a constant value to all pixels only in RED components you will see the effect equivalent to moving the RED tab towards the right, and same way reducing the RED component by a constant will have the reverse effect.

Like wise you can increment/decrement each component by a fixed value as described. If you increase/decrease ALL RGB components by a same value this will be equivalent to change in brightness (basically you are adding/removing WHITE color).

The Cyan, Blue, Magenta - corresponds to such transformation in CMYK color space. (But i guess, Blue in this color space corresponds to mix of Cyan and Yellow. So that's a bit tricky. The transformation is essentially same for all.

The last element Tint: { Hue and Saturation } corresponds to same operations but here, the images is first transformed in HSV model and then HUE and Saturation is added/subtracted independently.

I don't know about the exact relationships of markers of the dial to the corresponding numbers but that you can figure out by trying out practical values.


Finding period of a square wave with varying sampling frequency

Pick a nice center value, say 1.0, scan across looking for crossings from down to up (or the reverse).  Mark the spot of your first crossing, count up some fixed number of crossings, call the count N.  Take your system timer value at your first crossing and your last crossing.  Divide the time interval by N and you'll have a very good estimate of the period.  How uneven the crossings are in the interval is immaterial.  If you can assure the two crossing that you are timing happen in "the same state", you will lose one more level of uncertainty.

You can't do it much simpler than that.  Increased accuracy will come from interpolating the crossing point rather than taking the closest sample's time.

----------------------

Follow up

You'll have to figure the Interpolate routine out yourself.  It's not that hard.

The logic is very straightforward and significantly less computationly heavy than an auto-correlation technique.  More accurate in this case too.

Pseudo-code

    FirstIndex = -1
    N = 0

    For i = 0 to PointCount
      If signal[i] <= 1.0 and signal[i+1] > 1.0 Then
         If FirstIndex < 0 Then
            T1A = System.Time(i)
            T1B = System.Time(i+1)
            FirstIndex = i
         Else
            TNA = System.Time(i)
            TNB = System.Time(i+1)
            LastIndex = i
            N += 1
         End If
      End If
    Next

    T1 = Interpolate( T1A, T1B, Signal[FirstIndex], Signal[FirstIndex+1] )
    TN = Interpolate( TNA, TNB, Signal[LastIndex], Signal[LastIndex+1] )

    Period = ( TN - T1 ) / N

You may not need to interpolate, (TNA - T1A)/N may be accurate enough.


Can someone explain complex mapping as appears in Ogata's textbook

If we denote $s=\sigma+j\omega$, then for the mapping $z=e^{sT}$, the negative real axis in the $z$-plane in the range $z\in[0,-1]$ corresponds to $\omega=\pm\omega_s/2$ (i.e. Nyquist) and $\sigma\in [-\infty,0]$, where $\omega_s$ is the sampling frequency in radians. So if the curve in the $s$-plane moved along the lines $s=\pm j\omega_s/2$ then the curve in the $z$-plane would move on the negative real axis of the $z$-plane. If it doesn't (as suggested in the figure) then this range of the $z$-plane will be excluded. So if with this mapping you want poles on the negative real axis of the $z$-plane, you need corresponding poles with imaginary parts $\pm \omega_s/2$ in the $s$-plane.


A discretized amplitude is directly obtained by the numerical value of an intensity attached to the pixel location. The brightest would the maximum intensity over the whole picture. It seems that in Python with numpy you get gt it as:
brightest = numpy.amax(imageData)

and in Matlab
brightest = max(imageData(:))


Additive but not homogeneous continuous system?

You are misinterpreting the issues entirely.

  A system with complex inputs and complex outputs that maps its input
  $x(t)$ to output $[x(t)]^*$, the complex
  conjugate of $x(t)$, is additive, is not homogeneous for complex scalars and is homogeneous for real scalars.

Consider the system with complex inputs and complex outputs that maps its
input to its complex conjugate.  Thus, it maps
$$x_1(t) \to [x_1(t)]^*, \quad x_2(t)\to [x_2(t)]^*.$$
So, what happens if the input is the sum $x_1(t)+x_2(t)$? well, the
rule is that the input gets mapped to its complex conjugate and so
$$x_1(t) + x_2(t) \to [x_1(t) + x_2(t)]^* = [x_1(t)]^* + [x_2(t)]^*.$$
Since this holds for arbitrary choice of $x_1(t)$ and $x_2(t)$, we
conclude that the system is indeed additive: the response to a
sum of signals is the sum of the responses.

Let $\alpha$ denote a complex
number. Then, this "complex conjugator" system maps
$$\alpha\cdot x(t) \to [\alpha\cdot x(t)]^* = [\alpha]^*\cdot [x(t)]^*
\neq \alpha\cdot [x(t)]^*$$
unless $\alpha$ happens to be a real number.  So, the system is not
homogeneous since the required condition does not hold for arbitrary
(complex) $\alpha$.
However, it is true that if $\alpha$ is restricted to be a real number,
then the scaling property (a.k.a. homogeneity) does hold.

Addendum:  The OP wants

  a system that maps strictly real inputs to strictly real outputs, and that cannot be scaled with a complex (only with real) constants, and that is additive, but not homogeneous. That's my initial question,....

Since the input and output are required to be real, and not complex, the
notion of "cannot be scaled with a complex constant" does not make sense
any more because the scaling of the input makes it a complex-valued signal,
and the system allows only for a real-valued input. So essentially the question
boils down to

"Does there exist a system with real-valued inputs and outputs
that has the additive property but not the homogeneity property?"

to which the answer is No, as pointed out in the
Wikipedia article on linear systems. Most of the ideas there are reproduced
in the dsp.SE answer
mentioned in a comment by MBaz but the dsp.SE answer goes on to make
the additional claim (not in Wikipedia) that

"not much else is needed to extend to $a \in \mathbb C$"

which is false.  Indeed, as MattL notes carefully in a comment,
MattL's answers are applicable to scaling by real numbers only.

An outline of the proof that for real-valued systems, additivity
implies homogeneity is as follows:

Consider a system with real-valued inputs and outputs that
maps any given $x$ to $A(x)$. The function $A(\cdot)$ has the
property that $A(x_1+x_2) = A(x_1)+A(x_2)$ and the system is called
additive since the response to the sum $x_1+x_2$ is the sum of
the responses $A(x_1)$ and $A(x_2)$.  Mathematical induction extends this to
$$A\left(\sum_{i=1}^n x_i \right) = \sum_{i=1}^n A(x_i).\tag{1}$$

Choosing all the $x_i$ in $(1)$ to have value $0$, we see that
$A(0) = nA(0)$ and so $A(0)=0$.
Choosing all the $x_i$ in $(1)$ to be $x$, we see that $A(nx)=nA(x)$
and so the scaling property holds for positive integers $n$.
Choosing $n=2$ and $x_2 = -x_1$ in $(1)$, we have that
$A(x_1+(-x_1)) = A(x_1)+A(-x_1)$. But $x_1+(-x_1)=0$ and so we
see that $A(0) = 0 = A(x_1)+A(-x_1)$, that is,
$$A(-x_1)=A((-1)x_1) = -A(x_1) =(-1)A(x_1).$$ So the scaling
property holds for $-1$, and by mathematical induction for
all negative integers as well.
Scaling holds for rational numbers $\displaystyle \frac 1n$ as well. Choose
$\displaystyle x_i = \frac 1n x$ in $(1)$ to get
$$A(x)=nA\left(\frac 1n x\right) \Rightarrow
A\left(\frac 1n x\right) = \frac 1n A(x).$$
Mathematical induction allows us to extend this to more general
rational numbers $\displaystyle \frac mn$ as well.
Finally, to show that scaling applies to arbitrary (not necessarily
rational) real numbers, Wikipedia says that we need to assume that
the mapping is continuous and use the fact the rational numbers
are a dense subset of the reals.

So, except for this last part which requires some mathematics
beyond high school, we can show straightforwardly that
additivity implies homogeneity.


Can I search Doppler shift in freqency domain

Note my function below showing how I simultaneously searched for code and carrier (Doppler) offset for spread spectrum signals such as GPS in one 2D IFFT computation*. The 2D IFFT is processing intensive, but this approach has a lot of merit if you need very fast ("one-shot") acquisition.    (*Note the total processing involves two 1D FFT's, an element by element matrix multiplication for all Doppler and Code bins, and a 2D IFFT).

The basic premise is using fft's to cross-correlate two sequences (your waveform that has code and carrier offsets, and your reference code):

$\textrm{XCORR} = \textrm{ifft}(\textrm{fft}(a)\textrm{fft}(b)^*)$

Specifically the circular cross correlation between two sequences is equal to the inverse fft of the complex conjugate multiplication of their respective fft's.

If you circularly shift one of the fft's before multiplying, you calculate the result for a different Doppler bin!  It is that easy.  For actual implementation detail I included my Matlab function below. Here is an example plot of the result:

One FFT is computed for the received sequence, and one FFT is computed for the reference code. A matrix of delayed FFT's is created for one of the two that only involves replicating columns with a circular shift of the FFT already computed (which is a Hankel matrix) , and then a matrix of the other is created by simply replicating the FFT computed in each column, matching the number of Doppler columns created. The two matrices are complex conjugate multiplied element by element ($M_{ij}= A_{ij}B_{ij}^*$) Then (this is where all the processing is needed), a 2D inverse FFT is computed of the product matrix resulting in the solution depicted in the surface plot above.

On one axis for the figure above is the normalized frequency offset (Doppler) and the other axis is the relative code delay, and the vertical axis is the correlation where the delay and Doppler were simultaneously determined. So with one capture (whose length depends only on post correlation SNR required of the received sequence) you can complete acquisition processing with the use of two ffts, an element by element matrix multiplication and a 2D IFFT.

function [fdout,freq_axis]= fdcorr(a,b,frange)
%FDCORR Joint frequency and delay correlation
%   Performs an efficient circular cross correlation between
%   two vectors over all possible frequency offsets (Dopplers)
%   between -pi and pi where 2pi is the sampling rate.
%
%[FDOUT, FREQ_AXIS] = FDCORR(A,B)
%
%   A:   vector representing received signal sampled at FS = 2pi
%   B:   vector representing correlation code sampled at FS = 2pi
%
%   FDOUT: 2D Correlation surface vs delay and frequency offset
%   FREQ_AXIS: Frequency values aligned to frequency axis (-pi to pi)
%
%[FDOUT, FREQ_AXIS] = FDCORR(A,B,FRANGE)
%
%   Limits frequency range as defined in FRANGE.    Typically
%   frequency offsets are much smaller than the Nyquist bandwidth, so
%   using a reduced range for FRANGE from -pi to pi will be more
%   efficient.
%
%  FRANGE: Two element vector from lowest to highest frequency to
%          search over within a frequency range from -pi to pi.
%          Example [-pi/10 pi/10]
%
%FDCORR(A,B,<FRANGE>) with no output arguments plots the surface
%       correlation magnitude vs frequency and delay offset in the
%       current figure window.  FRANGE optional.
%
%
% Dan Boschen 11/22/2009

%==========================================================
% VERSION HISTORY
%
% 11/22/2009     Inital Release
%
%
%==========================================================

%==========================================================
% APPLICATION NOTES:
%
% Vector "A" represents a received signal that has been encoded with
% a correlation code (such as GPS and CDMA applications) or a training
% sequence.  It has an unknown frequency offset due to Doppler offsets
% and frequency offsets in the receiver, and an unkown delay shift
% relative to the start of the code sequence represented by vector "B".
% FDCORR will efficiently find the correlation between the received
% sequence and reference code over frequency and delay offsets.
% Vectors A and B are interchangable.
%
% To detect the maximum correlation value and delay and freq location use:
% [temp1, temp2]=max(fdout);
% [maxcorr,freq_index]=max(temp1);
% delay =temp2(freq_index)-1;
% freq= freq_axis(freq_index);
%
% The basic operation uses the Cross-Correlation theorem:
% xcorr= ifft(fft(a).*conj(fft(b)))
%
% Note the detected frequency offset is course and is the
% frequency of the closest Doppler bin where each Doppler bin
% is seperated by 2pi/N
%
%==========================================================

%
%error handling
if nargin==3
    if length(frange)~= 2
        error('FRANGE must be a two element vector')
    elseif abs(frange)> pi
        error('FRANGE values must be between -pi and pi')
    end
end

%force column vectors:
a=a(:);
b=b(:);

%============================================================
%create all Doppler versions of a by shifting the fft:
%============================================================

N=max(length(a),length(b));
fft_a= fft(a,N);

%frequency axis for all doppler bins -pi to pi

%index is the "fft shifted" freq axis. Method below faster than
%using fftshift(0:N-1), but same result.
freq_index=[-floor((N-1)/2):floor(N/2)];
freq_axis= freq_index*2*pi/N;   %convert to -pi to pi range

if nargin==3     %reduced doppler range
    freq_include=find(freq_axis>=frange(1) & freq_axis<=frange(2));
    freq_axis=freq_axis(freq_include);
    freq_index=freq_index(freq_include);
    num_freqs=length(freq_index);
else
    num_freqs=N;
end

%create Hankel matrix for all Dopplers (each column is the fft
%of vector a at a different Doppler offset)
fft_idx=(1:N)';
fft_dopplers = fft_idx(:,ones(num_freqs,1))+freq_index(ones(N,1),:);
fft_dopplers= mod(fft_dopplers-1,N)+1;    %Hankel subscripts
fft_dopplers(:)=  fft_a(fft_dopplers);      %fill data

%alternate method if solving for all Doppler bins (processing intensive):
%fft_all_dopplers= hankel(fft_a,[fft_a(end); fft_a(2:end)]);

%==================================================================
%Circular Correlation using circ_corr= ifft(fft(a).*conj(fft(b))
%==================================================================
fft_code_mtx= fft(b,N)*ones(1,num_freqs);
fdout=(ifft(fft_code_mtx.*conj(fft_dopplers)));

%===============================================================
%plot result
%===============================================================
if nargout ==0
    surf(freq_axis,1:N,(abs(fdout)));
    shading interp;
    xlabel('Normalized Frequency (xpi rad/sample)');
    ylabel('Delay Index');
end


Differences between FFT from acceleration, velocity and displacement data of a dynamic system

I have the same problem and I made some simulation and calculations to understand this phenomenon.

Facts:

velocity is an integral of acceleration
displacement is an integral of velocity
sine wave having frequency f is sin( 2*pi * f * t ) where t is time in seconds
integral of sin( n*t ) is -1/n * cos( n*t )
double integral of sin( n*t ) is -1/n^2 * sin( n*t )
To make equations cleaner, I'll use k = 2*pi

Calculations:

If you have hypothetical acceleration measurement signal having frequencies 1Hz, 10Hz, 100Hz with equal amplitudes, then:

acceleration: sin(k*t) + sin(k*10*t) + sin(k*100*t)
velocity: -1/(k) * cos(k*t) - 1/(k*10) * cos(k*10*t) - 1/(k*100) * cos(k*100*t)
displacement: -1/(k^2) * sin(k*t) - 1/(k^2*10*10) * sin(k*10*t) - 1/(k^2*100*100) * sin(k*100*t)

So, it seems that:

Both acceleration, velocity and displacement signals have the peaks at the same place: 1Hz, 10Hz, 100Hz
the 10Hz and 100Hz peaks of velocity are lower and the peaks of displacement are MUCH lower that the peaks of acceleration
the 100Hz peak of displacement is 10,000 times lower that 100Hz peak of acceleration

Simulations:

I made short R script which created the acceleration signal above. Then I integrated it (cumulative sum) to get velocity and displacement signals and finally took fft of them. The results verified the calculations


What is the variance of DFT of Fourier coefficient of difference of a vector of white noise?

This appears to be a homework/quiz question; we don’t provide full solutions but the following hint should help:
All samples in a white noise process are independent and identically distributed (IID) so it is helpful to understand how the variance grows as you sum N IID random variables: the sum of $N$ IID random variables with variance $\sigma^2$ will have a variance of $N\sigma^2$. Note that the exponential in the DFT equation only changes the phase of each sample without changing the amplitude. Importantly, all the samples were independent from each other and identically distributed; how does this change when we change the phase of those samples?
If we roll 10 fair die, the face value of each die will be independent. If we then turn over each die by same or different rotations, would the new values also be independent?


Hough transform is your friend. Basicaly it deals with straight lines (and almost all fibers in your example are straight). You can find open source implementation in OpenCV (here is doc + result example).

Moreover there is some generalizations that able to detect custom shapes.
P.S. Since you're using C#, OpenCV wrapper could be helpfull for you.


You're right, there usually is a better IIR filter (if you have enough data). The discrete-time Wiener filter is not "by definition" FIR. It is common to constrain the filter to the FIR case because it's often more straightforward to implement, and because such a filter can be made adaptive more easily. Also, in practice you often want to consider only a finite data window, which makes the optimum Wiener filter FIR.

However, you don't need to impose the constraint that the filter length be finite. The most general causal discrete-time Wiener filter for data available over the infinite past is an IIR filter.

A derivation of the FIR and the IIR cases can be found (among others) in Digital Signal Processing by Proakis and Manolakis.

EDIT: If we only have a finite amount of data available, the optimum Wiener filter has a finite impulse response. Assume we have $N$ data points (one current and $N-1$ consecutive past data points), then the best we can do with a linear time-invariant filter is to linearly combine those $N$ data points, by multiplying each data point by a coefficient and adding them up:

$$y[n]=\sum_{k=0}^{N-1}h[k]x[n-k]\tag{1}$$

The coefficients $h[k]$ are optimized such that the mean square error between the filter output $y[n]$ and some desired signal is minimized. The $N$ coefficients $h[k]$ are the impulse response of the corresponding FIR Wiener filter.


Lpcm encoded data kept on an analog source

The theoretical max number of bits you can get is function of the bandwidth and the signal to noise ratio of your media. It's often referred to as "Channel Capacity", see
http://en.wikipedia.org/wiki/Shannon%E2%80%93Hartley_theorem

In practice this depends a lot on the specific properties of your channel (transfer function, noise properties etc.) and there are many technologies you can deploy effectively deal with these properties.


O(1) peak segmentation algorithm in matlab?

I'm not sure why you need to call imreconstruct more than once. The way you normally use morphological reconstruction is described here: You basically subtract some constant (the peak threshold) from your source image, then reconstruct it. What you get is an image that is identical to your original image except at the peaks:

So you simply take the difference and binarize, or compare the two images pixel-wise.


Need intuition on the frequency response of a signal

The frequency response tells you how a system reacts to a sinusoidal excitation of a given frequency, in terms of amplitude and phase.

As any signal, periodic or not, can be described as a sum of sinusoids, this allows you to deduce the response to any signal, by the superposition principle.

Interesting features of the frequency response are the continuous gain, the presence of resonances and the bandwidth. A flat curve reproduces the original excitation without deformations.

You get a feeling of a frequency response when looking at how a slowing spin-drier vibrates.


Implementing real-time signal processing filters

There is a recent (2015) 150+ pages review called Smartphone-Based Real-Time Digital Signal Processing with lab codes. There are chapters like:

2: Android Software Development
3: iOS Software Development
6: Real-Time Filtering
7: Adaptive Filtering

  Real-time or applied digital signal processing courses are offered as
  follow-ups to conventional or theory-oriented digital signal
  processing courses in many engineering programs for the purpose of
  teaching students the technical know-how for putting signal processing
  algorithms or theory into practical use. These courses normally
  involve access to a teaching laboratory that is equipped with hardware
  boards, in particular DSP boards, together with their supporting
  software. A number of textbooks have been written discussing how to
  achieve real-time implementation on these hardware boards. This book
  discusses how smartphones can be used as hardware boards for real-time
  implementation of signal processing algorithms as an alternative to
  the hardware boards that are currently being used in signal processing
  teaching laboratories. The fact that mobile devices, in particular
  smartphones, have now become powerful processing platforms has led to
  the development of this book, thus enabling students to use their own
  smartphones to run signal processing algorithms in real-time
  considering that these days nearly all students possess smartphones.
  Changing the hardware platforms that are currently used in applied or
  real-time signal processing courses to smartphones creates a truly
  mobile laboratory experience or environment for students. In addition,
  it relieves the cost burden associated with using a dedicated signal
  processing board noting that the software development tools for
  smartphones are free of charge and are well-developed. This book is
  written in such a way that it can be used as a textbook for applied or
  real time digital signal processing courses offered at many
  universities. Ten lab experiments that are commonly encountered in
  such courses are covered in the book. This book is written primarily
  for those who are already familiar with signal processing concepts and
  are interested in their real-time and practical aspects. Similar to
  existing real-time courses, knowledge of C programming is assumed.
  This book can also be used as a self-study guide for those who wish to
  become familiar with signal processing app development on either
  Android or iPhone smartphones.


Why is the frequency response different based on the amplitude of the impulse I feed into the filter?

The lower amplitude the input signal has compared to the resolution of the filter coefficients the weaker is its ability to shape the frequency of the input signal. It also depends on the span of the filter coefficients.

You get the same result if you feed the filter with 32-bit coefficients with a signal having an amplitude covering only the lower 16-bits. If you want to work around this effect you need to amplify the input signal.


What happens when the poles of this z-transform function are outside the ROC for a signal?

Each ROC has it’s own impulse response  .  There isn’t a general impulse response or a single ROC in most cases.

Your problem can have 3 ROC with 3 impulse responses and 3 stability conditions

One ROC is the causal ROC, Another is the anti-causal ( running the filter in reverse time) and the last, which isn’t common has a name that requires me to look it up, but this isn’t my homework problem so i’m going to stay put.

I know that in the first edition of Oppenheim and Schaefer, they explicitly say you can’t write down an impulse response from the Z transform without stating the ROC


Power spectral density: where does the difference from decibel and regular representation come from?

If you zoom in on the linear curve, you would see that it's also "gradually decreasing", just much steeper, than the logarithmic curve. The logarithm simply makes it more visible.
And, yes, what you see is indeed spectral leakage. It can be reduced by using a non-trivial window function as opposed to the rectangular window being employed by default in the periodogram function. Try using a Hann window, that should reduce spectral leakage.


Does FFT magnitude of a sound drop to 1/10 of its value if sound level is decreased by 20dB?

-20dB corresponds to $1/10$ of the amplitude or $1/100$ of the energy. That's the same thing, since energy (all else being equal) is proportional to the square of the amplitude.

It also corresponds very roughly to about $1/4$ perceived loudness but that's an entirely different can of worms


I resolved my issue. The problem was with timing. My samples were rendered correctly, but the problem was with waiting for the buffer. I switched to event driven approach instead of manually calculating the buffer length and waiting, and the problem was resolved.

I don't know whether this question might be useful for someone. If not or if the source of the problem is not related to DSP, I will delete it. But for now I'm keeping it here with my answer.


Is Parks-Mcclellan algorithm guaranteed to stop?

You have to distinguish between the algorithm and a specific implementation. The Parks McClellan algorithm is based on the Remez exchange algorithm, which is guaranteed to converge. So in principle any correct implementation of the algorithm should converge. However, convergence of the algorithm is only guaranteed if all intermediate results are computed with infinite precision, which is not the case in a practical implementation. Some design specifications can result in numerically ill-conditioned problems that will prevent correct convergence of the algorithm.

In short, yes, the algorithm is guaranteed to converge. However, in some cases, numerical problems can prevent convergence to the optimal solution. But even in the usual case  of convergence, the number of necessary iterations cannot be determined beforehand.


What's the advantage of adaptive IIR filter against FIR?

These are the key differences between FIR and IIR filters, regarding the feature you wish to control are the following:

$$
\begin{array}{c|lcr}
\text{Feature} & \text{IIR} & \text{FIR} \\
\hline
\text{Implementation} & \text{Poles & Zeros} & \text{Zeros Only} \\
\text{States} & \text{Yes} & \text{No} \\
\text{Phase Delay} & \text{*} & \text{Half Integer} \\
\text{Stability} & \text{*} & \text{Always} \\
\text{Ripple} & \text{Yes} & \text{*} \\
\text{Cut-Off} & \text{Yes} & \text{*} \\
\end{array}
$$

The * indicates the feature can be controlled, by adding orders in most cases.

The standard definitions of FIR and IIR filters are:

FIR:

$$H(z)=b_0z^0+...+b_nz^n$$
$$y(t)=b_0u(t)+...+b_nu(t-n)$$

IIR:

$$H(z)=\frac{b_0+b_1z^1+...+b_nz^n}{1+a_1z^1+...+a_nz^n}$$
$$y(t)=b_0u(t)+...+b_nu(t-n)-a_1y(t-1)-...-a_ny(t-n)$$

$u$ is the input, $y$ is the output, $x$ is the states (below), $t$ is the time, scaled by a sampling time $dt$, $n$ is the number of orders of the filter. Each filter has $n$ size coefficient vectors, plus constant direct-output term $b_0$ (optional), and $a_0$=1. For simplicity assume $\sum{b_i}=1$ and $\sum{a_i}=1$, though this is not required anywhere.

Implementation. By definition, FIR include zeros only, leading to a linear system in the history vector for $u$: $[u(t-1)...u(t-n)]$.

IIR include both poles and zeros, also leading to a linear system in the history vector not only for $u$, but for $y$ too. Because of this, by one side IIR can be unstable; but by other side, they can be designed to have smooth ripple and sharp cutoffs with a minor number of orders.

States. FIR are static systems in the history vectors, meaning the filter is not dynamical, do not have states, is not recursive, no feedback.
IIR are dynamical systems in the history vectors, meaning the filters have states, is recursive, has feedback, hence have "memory" from past inputs & outputs.

Phase Delay. The Phase Delay $\tau_\phi$

$$y(t)=y_0(t-\tau_t) sin(\omega(t-\tau_\phi)+\theta)$$

can be easily controlled in FIR implementations.
If $b_k=b_{n-k}$,$k=0...n$, the phase delay is constant, equal to $n/2$ (the center of the FIR coefficients shape, its impulse response), equal to the group delay, and thus the filter become linear phase, with phase equal to $\omega\tau_phi$.

Because IIR have infinite impulse response, they can be minimum phase instead of linear phase, though the phase achieved can be much less than the phase of a FIR for the same number of orders.

Stability. FIR are always stable, IIR can be designed to be stable, if stability is required.

Ripple. IIR can be designed to be flat-ripple both in pass-band|stop-band|both (butterworth|chebyshev|elliptic), FIR requires a major (tending to "infinite") number of orders for equate this property.

Cut-Off. IIR can be designed to have a sharp cut-off or narrow transition bands, FIR requires a major (tending to "infinite") number of orders for equate this property.

Related Articles:

https://ocw.mit.edu/courses/electrical-engineering-and-computer-science/6-341-discrete-time-signal-processing-fall-2005/lecture-notes/lec08.pdf
https://www.quora.com/Why-are-FIR-filters-preferred-over-IIR-filters
http://iowahills.com/A8FirIirDifferences.html
http://forums.prosoundweb.com/index.php?topic=2045.0
http://www.vyssotski.ch/BasicsOfInstrumentation/SpikeSorting/Design_of_FIR_Filters.pdf


You must have understood the notion of digital linear modulation or discrete time vs continuos time (see Chapter 2).

Another reference.

OFDM can be thought as FDM with sinc pulse whose delay-$T$-shifted versions form an orthonormal basis. In frequency domain, they are seperated by $1/T$ which is denoted $\Delta f$ , i.e. subcarrier spacing. You fix a $T$, then $\Delta f = 1/T$; you fix $\Delta f$, then $T = 1/\Delta f$.

Look in frequency domain, because data is modulated in seperated subcarriers, the final signal is the sum of these subcarriers. Because the Fourier transform preserves the sum operation, in time domain, the final signal is also the sum of inverse-Fourier-transformed subcarrier signals that have the same duration $T$. Thus the final signal, which is sum of several signals having the same duration $T$, has the same duration $T$ that is called OFDM symbol duration or FFT duration or FFT period.


DCT is not used only for going back to time domain, it is also used due to it's "de-correlating properties". By saying that I mean compression of most of the information in smaller number of coefficients. Although in theory we could use the DCT to calculate the spectrum, but I never saw anyone doing that. Usually thing is simple:

want to analyse the spectrum? Use DFT
want to compress and de-correlate data? Use some flavour of DCT

Also one very important thing. In case of DCT, the equivalence property of convolution and multiplication in frequency domain doesn't hold:
$$x[n]\star y[n] \nLeftrightarrow X_{DCT}[k]\cdot Y_{DCT}[k] $$
which is the fundamental property (homomorphic) used when calculating the MFCC's. It allows to perform the Cepstral Mean Normalization and remove the "convolutional distortions", such as vocal tract. For more info please refer to my previous answer.


You can look at these set of notes provided by Analog Devices as part of their "Basic Linear Design" seminar, or just look for notes on second order filters.

The transfer function of a second order high pass filter is
$$H(s) = \frac{A s^2}{s^2+\frac{w_0}{Q}s+w_0^2}$$ where $w_0$ is the cutoff frequency, $Q$ is the $Q$-factor, and $A$ is the gain.


This question is about modelling mechanical response of behavior in terms of signals

I dont know what this axis([-.76 -.16 -15 5])’ means

As far as I am concerned axis([-.76 -.16 -15 5]) set minimal and maximum values for x axis (first and second values) and y axis (third and forth values) in plot.

For example:

x = linspace (-10, 10, 100);
y = x .^ 2;
plot (x, y, 'r', "linewidth", 2);
title ("S");
xlabel ("x");
ylabel ("y");

x = linspace (-10, 10, 100);
y = x .^ 2;
plot (x, y, 'r', "linewidth", 2);
title ("S");
xlabel ("x");
ylabel ("y");
axis ([-5, 5, -5, 20]);


Why do we choose cutoff frequency to be taken for "1/sqrt(2) magnitude reduce"?

I presume that red line I drew represents that components of a signal
  would neither increase or decrease in magnitude.

That would be a wrong presumption. The magnitude of the frequency response is (almost) never truly flat, and you would see this the more you zoom in. That's why you need to pick a convention of what exactly the "band edge" is.

As with most conventions, it's somewhat arbitrary. Half power has the advantage that if adjacent filters have the same cutoff frequency the total power will be flat. So if you add the output of an 1 kHz high pass and a 1 kHz low pass, the energy sum of the outputs will be flat with no dips or peaks in the cross over region.

-3dB is also not the only convention either: a Linkwitz Riley filter typically uses -6dB at the nominal cutoff. -3dB also doesn't make a lot of sense for parametric EQ filters and for high and low shelve filter.

So it all depends on context and one should never automatically assume -3dB.


What is the instantaneous frequency of this signal

Concerning the comment from Juancho, we can assume that the time is normalized to e.g. 1s, meaning that $t=1$ corresponds to $1s$. Then, we dont have a problem with the units (this a very common operation and not specific to your problem).

Now, regarding the "instantaneous frequency" (also check Wikipedia and another question on math.se), it is not rigorously defined. Due to the uncertainty principle you normally cannot assign an accurate specific frequency to an infinitely short interval of time (which would the term "instantaneous" correspond to). So, in general, the question for an "instantaneous frequency" needs to be asked with care and you need to understand that its interpretation might not be as straight-forward as it seems.

However, in case of a simple sine wave, consider the simple expression

$$x(t)=\cos(2\pi \cdot f\cdot t),$$

where we know that its frequency is $f$. The argument to the $\cos$ is $\phi(t)=2\pi f t$. If we calculate $\frac{d\phi(t)}{dt}=2\pi f$, which tells us the phase of the cosine changes with rate $2\pi f$. This term can be considered as the instantaneous angular frequency (i.e. frequency is $f$, as we would expect) of the sine wave.

Let's extend this expression to the general one

$$x(t) = \cos(\phi(t)) $$

where $\phi(t)$ is the phase of the cosine. We can now define the rate of change of the phase as the instantanous angular frequency: $\omega(t)=\frac{d\phi(t)}{dt}$.

With this idea, you can calculate the instantaneous frequencies (there are multiple) of your signal

$$\phi_1(t)=\pi t, \quad\phi_2(t)=\pi t^2,\quad \phi_3(t)=\pi\cos(\pi t).$$

So, the instantaneous angular frequencies are given by

$$\omega_1(t)=\pi,\quad\omega_2(t)=2\pi t,\quad\omega_3(t)=-\pi^2\sin(\pi t).$$

I.e. your signal consists of three different instantaneuos frequencies.


why does butterworth IIR filter for a square pulse show ripple on edges in time domain but not the same for sine wave

This happens because when lowpass filtering, you remove the high-frequency content that gives the square pulses their sharp edges. Now, if you're not yet versed in Fourier analysis of signals this will be a bit tough, but I'll provide some graphs to show what I mean.
Let's take a simple pulse and a sinusoid like the ones you have. We're going to filter them with a simple 4-th order Butterworth filter. Below are the time-domain and frequency-domain plots of the signals. The spectrum of the filter is also included. The spectrum plots are in decibels (dB).

If you don't know, filtering in the time domain is multiplication in the frequency domain. So you can look at the frequency spectrums and multiply them together, seeing where the attenuation takes place. You can see that the pulse has more significant contributions from high frequencies compared to the sinusoid. So after filtering, the sinusoid will be less distorted than the pulse which is exactly what we see below

We can take this a step further and remove more frequency content by using a tighter filter. Below are the same spectrums with the new filter. You see that the pulse is now even more distorted, while the sinusoid remains relatively unchanged.


Different phase response plots using MATLAB's $\tt freqz$ and $\tt fvtool$?

Citing the help output of freqz:

[H,W] = freqz(B,A,N) returns the N-point complex frequency response
vector H and the N-point frequency vector W in radians/sample of
the filter.

So, your output is not amplitude and phase, but (complex) transfer function and frequency. Hence, to plot the phase using freqz you need to use

plt(w, angle(h));


You are trying to design a FIR allpass filter. That's impossible (other than a simple delay) since true allpass filters consist of zero/pole pairs that are inverse of each others and the FIR poles are always at the origin.

However, with enough resources, you can always get "close enough" for your specific application. Of course, you would have to define what "close enough" means.

The actual design is simple enough

Sample your filter target on a very dense FFT grid
Do an inverse FFT
Window/truncate the resulting impulse response until you meet your "close enough" requirements. Typically the impulse response will decay exponentially, so truncation often gives better results that any type windowing (for the same filter length)
In most cases the impulse response will be non-causal, so you will need to add bulk delay.

Below is a piece of Matlab code that does that

%% start with a Butter filter
fs = 44100; % sample rate
bwCut = 1000; % 1000 Hz cut off frequency
bwOrder = 3;  % thrid order
%% design the filter
[b,a] = butter(bwOrder,bwCut/(fs/2));

%% sample on a large FFT grid
nFFT = 16384;
delta = zeros(nFFT,1); delta(1) = 1;
bwIR = filter(b,a,delta); % impulse response
bwTF = fft(bwIR);         % transer fundction

%% design the target
targetTF = exp(-1i.*angle(bwTF)); % that works every where except for Nyquist, where
                                  % the phase is undefined, we'll manually smootht this out
targetTF(nFFT/2+1) = 1;           % Nyquist phase repair
h0 = real(ifft(targetTF));        % target impulse response

%% now we have to cut it to a practical size
% let's see if we can capture all samples that are above a certain
% thershold. That's in essence a rectangular window
thresh = max(abs(h0)).*0.01;      % -40dB re max
tmp = circshift(h0,nFFT/2);       % shift to the middle of FFT window
i1 = find(abs(tmp) > thresh);     % find the area above threshold
h1 = tmp(min(i1):max(i1));        % grab that area

%% analyse the result. This really depends on the application.
% In this example we will look at the difference of the resulting
% transfer function to the "ideal" answer. The difference looks at both the
% real and imaginary parts, so it's a decent metric that combines magnitude
% and phase. It coould be interpreted as the amount of "noise" induced by
% the truncation error.
idealTF = abs(bwTF);

% model the results, rectangular window
h2 = conv(h1,bwIR);                   % convolution
h2 = h2(1:nFFT);                      % cut off the excess length
h2 = circshift(h2,min(i1)-nFFT/2-1);  % time align
diffTFRect = fft(h2)-idealTF;         % calculate difference
% same length hanning window
h1w = h1.*hanning(length(h1));
h2 = conv(h1w,bwIR);                  % convolution
h2 = h2(1:nFFT);                      % cut off the excess length
h2 = circshift(h2,min(i1)-nFFT/2-1);  % time align
diffTFHanning = fft(h2)-idealTF;

%% plot it
clf;
i = (1:nFFT/2);
fx = [diffTFRect diffTFHanning];
plot(i*fs/nFFT, 20*log10(abs(fx(1+i,:))));
xlabel('Frequency in Hz ->');
set(gca,'xlim',fs*[1/nFFT 1/2]);
set(gca,'xscale','log');
set(gca,'ylim',[-120, 10]);
ylabel('Difference in dB');
grid('on');
legend('Rectangle','Hanning');
title('Error');


Effects of sampling rate approaching the Nyquist frequency? (not undersampling)

Sampling at 2*Fn is suitable only for signals that are infinitely long and perfectly bandlimited to below Fn, which conditions don't exist in actual signals being sampled or resampled (given non-perfect low-pass filtering), or in any non-infinite length signals.

A finite length rectangular window will itself introduce higher frequency content to the signal being sampled or recorded.

Therefore, there are major advantages to a sample rate higher than 2*Fn, in terms of not corrupting the signal, especially if any higher frequency (close to Fn) spectral content near the beginning or end of your recorded signal is important.

Added:  Reconstructing a signal represented by data at a sample rate well above twice the highest frequency can greatly reduce the computational (or hardware) requirements of the reconstruction filter or interpolator.


Unless your flow graph contains signal processing blocks that you wrote in python (See the "embedded python block"), no, C++ will not make your flow graph be faster.
The reason is simple: Python is only used to set up the flow graph (read the generated python code!), i.e. to instantiate and connect the blocks together.
However, the GNU Radio scheduler and most blocks are running without interacting with the Python world at all. Hence, the GIL is needed at no point at all. Your Python calls top_block.run() and that's all that Python sees – the processing happens off-thread separately.


Interpolating/Decimating CIC Filter Group Delay

Consider an $D$-tap FIR filter with liner phase, the group delay (measured in samples) is $$g=\frac{D-1}{2}\tag{1}$$
and therefore, if it is measured in seconds it will be $$g=T_s\frac{D-1}{2}\tag{2}$$
where $T_s=1/F_s$.
The CIC filter which is also denoted as recursive running sum filter is indeed a special implementation of a moving-average filter. The moving average filter is
$$y[n]=\frac{1}{D}[1+z^{-1}+z^{-2}+\cdots z^{-D+1}]=\frac{1}{D}\sum_{n=0}^{D-1}z^{-n}=\frac{1}{D}\frac{1-z^{-D}}{1-z^{-1}}$$
which can be implemented by a $D$-tap FIR filter. If you compare carefully, the above response is identical to the frequency response of each stage of a CIC filter.
Therefore, we can use $(1)$ and $(2)$ for calculation of the group delay in each stage.
Comparing the response given for moving average with that of the CIC filter, assuming $R$ is rate change and $M$ is the number of samples per stages, then $$D=RM$$ So for $N$ stages we have $$g_{\text{cic}}=N\left(\frac{RM-1}{2}\right)$$ expressed in samples and
$$g_{\text{cic}}=NT_s\left(\frac{RM-1}{2}\right)$$ in seconds.


Hint: this is a lot easier if you use the Euler formula https://en.wikipedia.org/wiki/Euler%27s_formula

Then the integral simply becomes
$$X(n)=\frac{1}{L}\int_{0}^{2L}e^{(-a+i\omega)\cdot t}, a_{n}=real[X(n)], b_{n} = imag[X(n)]$$

which is pretty straight forward to solve


Is it meaningful to find linear convolution of just two random sequences?

No, if the two sequences are random then convolving them is not useful.