Does MatLab's `ftt` function gives a Fourier Transformed signal in terms of frequency $f$ or angular frequency $\omega$?

There are two things to clarify here before the answer becomes self-evident.

The trigonometric functions $\cos, \sin, \dots$
The relationship between phase $\phi$ and frequency $f$.

The trigonometric functions $\cos(\phi)$ and $\sin(\phi)$ trace the periphery of a circle. Where we are on the circle is denoted by the argument of the function \phi. But, in general, once we are around the circle once, the values keep repeating periodically.

Now, suppose that we say $\cos(4)$ ... Where is that on the circle? Is this on the "4" or is it further down its "multiples" sequence?

If we said, $\cos(2 \pi 4)$, then we know this represents 4 times around a circle of $2 \pi$. If we say $\cos(2 \pi 0.5)$ then we know that's halfway around a full circle.

Please note that so far we are talking about instantaneous values, akin to saying "What is the cos of 60 degrees? What is the cos of (360 + 60) degrees? What is the cos of (2*360 + 60) degrees? What is the cos of (k*360 + 60) degrees?" and so on.

Let's bring time in to this. Let's now say, not "Where am I, around the circle?" but "How fast am I going around the circle?".

$$ y = \cos(2 \pi t) $$

Where $t$ is time and time is a variable we cannot stop, it starts (usually) at $0.0$ at some arbitrary set off point (i.e., $t=0$ when we power on a system) and increases for ever beyond our control. It also has units of seconds.

In $1.0$ second, we have been around the $2 \pi$ circle once, in $2.0$ seconds we have been around it twice...and so on.

This represents periodic events with a period of $1.0$ sec, the time it takes to go around the circle.

What if we have things that happen more frequently than once every second?

This is easy to deduct, we just have to go around the circle faster. How much faster? I don't care, some factor $f$ faster than 1 rev per second.

$$ y = cos(2 \pi t f + \theta )$$

So, now, the phase at which we find ourselves moving around the circle is $\phi = 2 \pi t f$ plus some $\theta$ term to account for any starting point around the circle.

If we differentiate phase with respect to time we get frequency. That is, the rate at which we move around the circle.

The rate has a $2 \pi$ factor, yes, but that offsets everything and means number of times around the circle. The key factor here is $f$.

So, $2 \pi 440$ radians per second is another way of saying 440 Hz.

Let's now turn to the Fast Fourier Transform (FFT).

We have a system that samples at an $Fs = 44.1$ Hz. We use it to acquire a buffer of 441 samples which we then pass through the FFT. The FFT returns a result that contains 441 complex values, symmetrical about the middle point. The question is, what is the frequency of "bin" 44?

It is $44 \times \frac{44100}{441} = 4400$ Hz. And in this formula we used Hz throughout.

Hope this helps.


Why doesn't a sudden loud noise sound high pitched?

Pitch is a human psychoperceptual phenomena, which does not work the same way as a Fourier transform.

IIRC: Even if an FFT shows strong high frequency components in a short impulse or burst, the human ear-brain-system needs to detect roughly somewhere on the order of around 6 repetitions/periods/cycles of a periodic or pseudo-periodic acoustic waveform before communicating "I hear a pitch" to consciousness.  (depending somewhat on, but nowhere near linearly with, frequency). A single or duple period of even a pure sinewave won't do it.


LPC reduces to AR modelling only if the stochastic time process is stationary (does not change distribution parameters over time) and ergodic (average over time is equivalent to mean of ensemble average).

This connection between auto-regressive coefficients and the autocovariance of the process is described by the Yule-Walker-Equations (mentioned in the book you linked).

The filter coefficients in your book are defined as conjugates to make definitions more easier (as mentioned below equation 2.5, e.g. equation 2.10).


Is HSV Color Space Sufficient for Rudimentary Color ID and Edge Detection

The approach seems reasonable.
Indeed doing edge detection in weighted RGB channel is the classic approach (Though you could also employ more advance methods, See Edge Detection on a Color Image).

I think you could achieve great results if you also look specifically for oval shapes then you reduce the chances for false positives.

Color identification in the H and S channels will give most information (Well, up to differences in luminosity of the color).


How to sketch the following discrete time signal $\sum_{k=0}^{n} \delta[n-k]$?

Your signal is:
$$
x[n]  = \sum_{k=0}^n \delta[n - k]
$$
and you're saying you want to sketch it.

Based on your comment, I'm going to assume that
$$x[-1] = \left . x[n] \right  |_{n=-1} =  \left . (\delta[n] + \delta[n+1]) \right  |_{n=-1} = 0 + 1 = 1$$
and that
$$x[+1] = \left . x[n] \right  |_{n=+1} =  \left . \left(\delta[n] + \delta[n-1]\right) \right  |_{n=+1} = 0 + 1 = 1$$

I believe that means
$$
x[n] = 1, ~~ \forall n \in \mathbb{Z}
$$
because
$$
x[n]   = \cdots + \delta[n - (n-2)] + \delta[n - (n-1)] + \delta[n - n] + \cdots
$$
if summing up or
$$
x[n]   = \cdots + \delta[n - (n+2)] + \delta[n - (n+1)] + \delta[n - n] + \cdots
$$
if summing down.

It's a little odd having the upper limit of the summation depend on the current time and, as I said in the comment, having the summation go down to the upper limit if the upper limit is below the lower one.


You can, in a sense, if you can also sample your signal at an incommensurate sampling rate for the same period.

Any periodic signal can be decomposed into its constituent tones with a DFT and a proper algorithm.  (I have an algorithm I developed about 10 years ago that I have in the pipeline for my blog that I recently found out that McLeod described about 20 years ago.)

So, reduce your question down to a single pure tone.  The frequency for a pure tone, or an alias tone, can be found exactly (see my blog if you doubt this) for any frequency.

So, reduce your question down to a pure tone with a whole number of cycles in your sample frame.

The limitation with a DFT due to the discrete sampling is the aliasing effect.  If a tone is higher than the Nyquist frequency it will appear in the DFT as a frequency less that the Nyquist limit.  So, if a tone appears in bin $k$, where $k < N/2$, the will be a conjugate value in bin $N-k$, which is over the Nyquist limit.  In actuality, you can't tell with the frequency is $k$ cycles per frame, or $N-k$.  The frequency could also be any multiple of N over those.

The possible frequency can be: $k$, $N-k$, $k+N$, $2N-k$ ...

By sampling at two different truly incommensurate rates, the sets of possible frequencies from each of the DFTs will only have one value in common, and that will be the actual frequency.  If the rates are nearly incommensurate, the next possible common value will be well out of range.

Here is a concrete example:

60 Hz Sampling  N = 60

49 Hz Sampling  N = 49

             49                   60

 5 Hz   5, 44, 54, 93 ...    5, 55, 65, 115 ...

47 Hz   2, 47, 51, 96 ...   13, 47, 63, 107 ...

97 Hz   1, 48, 50, 97 ...   23, 37, 63, 97  ...

Now, suppose you thought you could be clever and use every other point of your 60Hz sampling, and every third point as your two different sampling rates.

One second sample frame

30 Hz Sampling  N = 30

20 Hz Sampling  N = 20

Possible Frequencies:

             20                                30

 5 Hz   5, 15, 25, 35, 45, 55, 65 ...   5, 25, 35, 55, 65 ....

47 Hz   7, 13, 27, 43, 47, 63, 67 ...  13, 17, 43, 47, 73 ....

You can't tell whether the 5 Hz is 5, 25, 35, 55, etc.

You can't tell whether the 47 Hz is 13, 43, 47, etc.

So that doesn't work so well.

Hope this helps.

Ced

Followup:

This technique assumes you have "ideal" instantaneous sampling, so it is more theoretical than practical.


Filter design with a phase delay specification?

Least square errors works well with FIR filters. General IIR filters are more difficult and typically require an iterative search algorithm. One specific type of IIR, named warped FIR filters, can also match arbitrary amplitude and phase response with a least square errors approach.


Well, your channel model might give you an explicit PDF, but:

You need to realize that a channel is not a scalar, ie. your PDF doesn't describe a random variable $\in \mathbb C$ or $\in \mathbb R$, but usually something like an impulse response, that is, something like a function of delay (so, your random variable is $\in \mathbb R \times \mathbb R$ or $\in \mathbb C \times \mathbb R$).

I'm not an acoustics guy myself, but if I've learned one thing it's that at least for indoor acoustics, the "flat channel assumption" (that lets us wireless communication folks model channels as simple complex number) doesn't hold, and that your channel is very frequency-selective. Often, it's time-variant, and sometimes even non-linear.

So, I'm not quite sure what to recommend here:

If you're looking for soft decision in wireless communications, I personally think that Proakis is quite OK. That's one of the standard digital communications textbooks.
That's up to your channel model, and you'll have to specify that.
An AWGN channel is certainly not an adequate approximation of an audio channel.


What does boundary discontinuity in DFT imply?

Because the left (top) edge of an image is unlikely to be a reflection of its right (bottom) edge, there are discontinuities all along the edges of an image when it is viewed as N-periodic. These discontinuities are represented in the frequency domain by high-frequency coefficients. From an image compression point of view, using the 2-D DFT as a transform would require to store these high frequency coefficients (or non-zero terms in the decomposition of the image as a sum of elementary components), which do not describe anything meaningful in the image itself, but are a nasty side-effect of an inappropriate choice of transform.


Factor $|a|^{-1/2}$ in definition of mother wavelets

My answer is for real scale $a$ and the fact that wavelet transform is usually defined in $L_2$ with norm

$$||\Psi(\tau)|| = \int_\mathbb{R} \Psi(\tau)\Psi^*(\tau)\mathrm{d}\tau $$

So

$$||\Psi_{a,t}(\tau)|| = \int_\mathbb{R} \frac{1}{|a|}\Psi(\frac{\tau-t}{a})\Psi^*(\frac{\tau-t}{a})\mathrm{d}\tau$$

Set $\tau' = \frac{\tau-t}{a} \implies d\tau' = d\tau / a \implies d\tau = a \times d\tau'$

$$||\Psi_{a,t}(\tau)|| = \int_\mathbb{R} \frac{1}{|a|}\Psi(\tau')\Psi^*(\tau') \times a \times \mathrm{d}\tau' = \int_\mathbb{R} \Psi(\tau')\Psi^*(\tau') \mathrm{d}\tau' = ||\Psi(\tau)||$$

Update for the question "How does it not change the norm considering that it could introduce a negative sign with the norm being sign-less?"

We develop stuff for the case that scaling factor $a < 0$.

with $A > 0$
$$||\Psi_{a,t}(\tau)|| = \lim_{A \to \infty}\int_{(-A-t)/a}^{(A-t)/a} \frac{1}{|a|}\Psi(\tau')\Psi^*(\tau') \times (-|a|) \times \mathrm{d}\tau' \\
= \lim_{A \to \infty}\int_{(A-t)/a}^{(-A-t)/a} \Psi(\tau')\Psi^*(\tau') \times \mathrm{d}\tau' \\
= \lim_{A \to \infty}\int_{(-A+t)/|a|}^{(A+t)/|a|} \Psi(\tau')\Psi^*(\tau') \times \mathrm{d}\tau' = \int_{-\infty}^{+\infty} \Psi(\tau')\Psi^*(\tau') \times \mathrm{d}\tau'$$


Filter design for frequencies in range 0-60Mhz

As you stated there are some parameters associated with filters such as type of filter (butterworth, elliptic, chebyshev I,chebyshev II, bessel,...), order of filter, etc. (Note that there are lots of considerations in filter design bases on the application)

Let's first talk about the order of filters:
Filters with higher orders let you have a filter with better quality. For example in lowpass filter if you want to have a sharper transition band you need to design a filter with high order filers. Order of a passive analog filter is determined by the number of elements you use in filter (resistor, inductor, capacitors). This was all I knew about the order of the filter.

Now let's talk about the most applicable kinds of filter types. The most important difference between filter types is the band they have ripples. Note that the priority of filters is due to the application and when we say one filter is better than the other we must state that in what aspects (order, ripples, poles, zeros, stability, group delay, phase response, minimum phase, ...) this is better than that one. Picture below compares 4 important types of filters designed with the same order. As you can see, Elliptic filter has ripples both in pass-band and stop-band while it has the sharpest transition, and butterwoth filter has no ripples but has least sharpness among filters of the same order.


Complex Data type not accepted in Simulink

The data type of the signal must be fixed in Simulink. Since you use a square root of a real signal, the output will be real if the data is positive, imaginary if negative. Simulink cannot allow that.
The workaround is to use force your real signal to be a complex signal by setting its imaginary part to 0. Then feed the signal to the square root block. The output will always be complex.


What is the advantage of applying a Hamming window for MFCC calculations?

Ideally you would only see that frequency as a narrow peak but because of the finite-length window you also get that other artificial crud: Windowing equals time domain multiplication by the window function (here rectangular or Hamming). Time domain multiplication equals frequency domain convolution: All frequencies will be replaced by the Fourier transform of the window function. This is what those look like:

The convolution spreads the true frequency to frequency bins around it. By using a window function other than the rectangular one there is less of the sidelobes compared to the main lobe, so less artificial long-distance spread making the result cleaner and better suited for frequency-selective analysis. You seem to have by accident chosen a frequency where you don't see the sidelobes with the rectangular window because the dips between the sidelobes land between your frequency bins. So you just see slopes.


Suggest me some good references for sonification

1)

Proceedings of annual International Conference on Auditory Display. The annual ICAD conference is the prime forum for academia and industry to exchange ideas and discuss developments in the field of auditory display.

https://smartech.gatech.edu/handle/1853/49750

2)

http://sonification.de/handbook/

This book is a comprehensive introductory presentation of the key research areas in the interdisciplinary fields of sonification and auditory display. Chapters are written by leading experts, providing a wide-range coverage of the central issues, and can be read from start to finish, or dipped into as required (like a smorgasbord menu).


Taking the FFT of a sinusoidal signal and going back

Part of the problem is that you are using sine waves.  The imaginary components of an FFT result contains all the information about pure sine waves (e.g. that are equivalent to a sin() function that starts with phase of zero at the start of the FFT aperture, and that are exactly integer periodic in aperture).  If you use only the real part of an FFT result, you can only see cosine waves (not sine waves, periodic in aperture), thus your silent result.

For shorter FFTs, the real component of the result might contain windowing artifacts if the frequency of your sine wave is between FFT result bins for the given FFT length (e.g. not integer periodic in aperture).  So, in synthesis, you are hearing an artifact, quantized in frequency.


How to improve the FIR filter performance degradation when the amplitude of the input signal decreases?

A linear filter just attenuate some frequencies more and others less. If you alter the SNR, the output of the FIR filter should change accordingly.
Many applications use more elaborate noise reduction methods in stead of (or in addition to) linear filtering. If you know that your signal is going to be a sine of known frequency, you may estimate its phase and amplitude embedded in noise using a matched filter. If you know that your noise is shot-noise, you could try to exploit that.
You noise seems to be uniform (non-zero mean) noise. Did you intend to call the normal distributed randn()?
-k


Data fusion using 2d discrete wavelet transform (DWT)

If I remember my wavelet transform correctly, this is what happens. Let's imagine that we are working with only the LOW image at the moment. I perform the first-scale decomposition which gives me the LOW_A0 (approximation) and the three detail coefficient images which are LOW_H0, LOW_V0 and LOW_D0. Then, to perform the second scale decomposition, you take the approximation image LOW_A0 and obtain LOW_A0_A1 (approximation) and the corresponding three detail coefficient images LOW_A0_H1, LOW_A0_V1 and LOW_A0_D1. Finally, in the last scale decomposition, you will split LOW_A0_A1 into its corresponding approximation coefficients LOW_A0_A1_A2 and the final detail coefficient images LOW_A0_A1_H2, LOW_A0_A1_V2 and LOW_A0_A1_D2.

So at the end you are left with the following:
1) The lowest approximation coefficients LOW_A0_A1_A2
2) The detail coefficients at each level which are 9 in number (LOW_H0, LOW_V0, LOW_D0, LOW_A0_H1, LOW_A0_V1, LOW_A0_D1, LOW_A0_A1_H2, LOW_A0_A1_V2 and LOW_A0_A1_D2).

Now, also consider that in parallel you have also performed the wavelet decomposition for the HIGH image and you have the same end stage products which are:
1) The lowest approximation coefficients HIGH_A0_A1_A2
2) The detail coefficients at each level which are 9 in number (HIGH_H0, HIGH_V0, HIGH_D0, HIGH_A0_H1, HIGH_A0_V1, HIGH_A0_D1, HIGH_A0_A1_H2, HIGH_A0_A1_V2 and HIGH_A0_A1_D2).

You want to now reconstruct a FUSED image which contains information from both the LOW and the HIGH images. At this stage, you would apply the fusion algorithm. Average the lowest approximation coefficients k*LOW_A0_A1_A2 + (1-k)*HIGH_A0_A1_A2, where k= 0.5 for equal weighting (you may want to weight one more than the other). Then, sum the detail coefficients at each of the corresponding LOW and HIGH detail coefficients. For example,
LOW_H0 + HIGH_H0.

Now take these approximation coefficient and detail coefficients and do the wavelet synthesis to get the final result. Does this make sense? I haven't read the paper/thesis you pointed out to though.


How did Scientists Capture First Images of Molecules Before and After Reaction?

The seminal paper you're looking for is:

Binnig, G., & Quate, C. F. (1986). Atomic force microscope. Physical Review Letters, 56(9),    930-933.

Here's the short version, along with some speculation since I don't have their manuscript.

Atomic force microscopy (AFM) is a scanning probe microscopy technique, meaning a sharp probe is scanned in a raster pattern over/in contact with the sample. A probe my look like this.

I will refer to the point as the tip, and that is where the carbon monoxide molecule mentioned in the article resides. The substrate the tip is attached to will be referred to as the cantilever.

A laser is bounced off the top of the cantilever, onto a dicrotic mirror, and finally into a photodetector, as shown below.

The researcher has control of the laser and mirror alignment. Alignment procedure is to place the laser centered, and as close to the end of the cantilever as possible, for the maximum deflection signal. The mirror is then adjusted to put the reflected laser on the most sensitive part of the photodetector.

The article mentions using non-contact AFM. That is about as meaningful as saying "we can the signal through a digital filter." There are multiple non-contact modes, but my guess (without having access to the manuscript) is they used a form of tapping mode.

The cantilever has a self resonant frequency, typically around 300 kHz. The cantilever is driven to resonance, and scanned over the sample. The deflection signal generated is an amplitude modulated version of the topography. It is also possible to obtain information about material properties this way.

That is the short version. If I get access to the manuscript or someone asks questions, I would be happy to add more information about how this system works.


Differentiation of sine in Fourier domain

I assume you mean the derivative with respect to $t$. In that case, the derivative of $\sin(\omega_0t)$ is not $\cos(\omega_0t)$ but $\omega_0\cos(\omega_0t)$. And luckily, this is also obtained via the Fourier transform relation you mentioned in your question:

$$\begin{align}\mathcal{F}\left\{\frac{d}{dt}\sin(\omega_0t)\right\}&=j\omega\cdot \mathcal{F}\left\{\sin(\omega_0t)\right\}\\&=\pi\omega[\delta(\omega-\omega_0)-\delta(\omega+\omega_0)]\\&=\pi[\omega_0\delta(\omega-\omega_0)-(-\omega_0)\delta(\omega+\omega_0)]\\&=\pi\omega_0[\delta(\omega-\omega_0)+\delta(\omega+\omega_0)]\\&=\omega_0\mathcal{F}\{\cos(\omega_0t)\}\end{align}$$

where I've used the fact that $f(\omega)\delta(\omega-\omega_0)=f(\omega_0)\delta(\omega-\omega_0)$ for any function $f(\omega)$ that is continuous at $\omega=\omega_0$. Consequently you have $\omega\delta(\omega+\omega_0)=-\omega_0\delta(\omega+\omega_0)$.


Synchrosqueezing Wavelet Transform explanation?

Synchrosqueezing is a powerful reassignment method. To grasp its mechanisms, we dissect the (continuous) Wavelet Transform, and how its pitfalls can be remedied. Physical and statistical interpretations are provided.
If unfamiliar with CWT, I recommend this tutorial. SSWT is implemented in MATLAB as wsst, and in Python, ssqueezepy. (-- All answer code)

Begin with CWT of a pure tone:

A straight line in the time-frequency (rather, time-scale) plane, for our fixed-frequency sinusoid over all time - fair. ... except is it a straight line? No, it's a band of lines, seemingly centered about some maximum, likely the "true scale". Zooming,

makes this more pronounced. Let's plot rows within this zoomed band, one by one:

and all superimposed, each for samples 0 to 127 (horizontal zoom):

Notice anything interesting? They all have the same frequency. It isn't particular to this sinusoid, but is how CWT works in correlating wavelets with signals.
It appears much of information "repeats"; there is redundancy. Can we take advantage of this? Well, if we just assume that all these adjacent bands actually stem from one and the same band, then we can merge them into one - and this, in a nutshell, is what synchrosqueezing does. Naturally it's more nuanced, but the underlying idea is that we sum components of the same instantaneous frequency to obtain a sharper, focused time-frequency representation.
Here's that same CWT, synchrosqueezed:

Now that is a straight line.

How's it work, exactly?
We have an idea, but how exactly is this mathematically formulated? Motivated by speaker identification and Empirical Mode Decomposition, SSWT builds upon the modulation model:
$$
f(t) = \sum_{k=1}^{K} A_k(t) \cos(\phi_k (t)), \tag{1}
$$
where $A_k(t)$ is the instantaneous amplitude and
$$
\omega_k(t) = \frac{d}{dt}(\phi_k(t)) \tag{2}
$$
the instantaneous frequency of component $k$, where we seek to find $K$ such "components" that sum to the original signal. More on this below, "MM vs FT".
At this stage, we only have the CWT, $W_f(a, b)$ (a=scale, b=timeshift); how do we extract $\omega$ from it? Revisit the zoomed pure tone plots; again, the $b$-dependence preserves the original harmonic oscillations at the correct frequency, regardless of $a$. This suggests we compute, for any $(a, b)$, the instantaneous frequency via
$$
\omega(a, b) = -j[W_f(a, b)]^{-1} \frac{\partial}{\partial b}W_f(a, b), \tag{3}
$$
where we've taken the log-derivative, $f' / f$. To see why, we can show that CWT of $f(t)=A_0 \cos (\omega_0 t)$ is:
$$
W_f(a, b) = \frac{A_0}{4 \pi} \sqrt{a} \overline{\hat{\psi}(a \omega_0)} e^{j b \omega_0} \tag{4}
$$
and thus partial-diffing w.r.t. $b$, we extract $\omega_0$, and the rest in (3) gets divided out. ("But what if $f$ is less nice?" - see caveats).
Finally, equipped with $\omega (a, b)$, we transfer the information from the $(a, b)$-plane to a $(\omega, b)$ plane:
$$
\boxed{ S_f (\omega_l, b) = \sum_{a_k\text{ such that } |\omega(a_k, b) - w_l| \leq \Delta \omega / 2} W_f (a_k, b) a_k^{-3/2}} \tag{5}
$$
with $w_l$ spaced apart by $\Delta w$, and $a^{-3/2}$ for normalization (see "Notes").
And that's about it. Essentially, take our CWT, and reassign it, intelligently.

So where are the "components"? -- Extracted from high-valued (ridge) curves in the SSWT plane; in the pure tone case, it's one line, and $K=1$. More examples; we select a part of the plane and invert over it as many times as needed.

Modulation Model vs Fourier Transform:
What's $(1)$ all about, and why not just use FT? Consider a pendulum oscillating with fixed period and constant damping, and its FT:
$$
s(t) = e^{-t} \cos (25t) u(t)\ \Leftrightarrow\ S(\omega) = \frac{1 + j\omega}{(1 + j\omega)^2 + 625}
$$

What does the Fourier Transform tell us? Infinitely many frequencies, but at least peaking at the pendulum's actual frequency. Is this a sensible physical description? Hardly (only in certain indirect senses); the problem is,  FT uses fixed-amplitude complex sinusoid frequencies as its building blocks (basis functions, or "bases"), whereas here we have a variable amplitude that cannot be easily represented by constant frequencies, so FT is forced to "compensate" with all these additional "frequencies".
This isn't limited to amplitude modulation; the less sinusoidal or non-periodic the function, the less meaningful its FT spectrum (though not always). Simple example: 1Hz triangle wave, multiple FT frequencies. Frequency-modulation suffers likewise; more intuition here.
These are the pitfalls the Modulation Model aims to address - by decoupling amplitude and frequency over time from the global signal, rather than assuming the same (and constant!) amplitude and frequency for all time.
Meanwhile, SSWT - perfection:

Is synchrosqueezing magic?
We seem to gain a lot by ssqueezing - an apparently perfect frequency resolution, violating Heisenberg's uncertainty, and partial noise cancellation ("Notes"). How can this be?
A prior. We assume $f(t)$ is well-captured by the $A_k(t) \cos(\phi_k (t))$ components, e.g. based on our knowledge of the underlying physical process. In fact we assume much more than that, shown bit later, but the idea is, this works well on a subset of all possible signals:

Indeed, there are many ways synchrosqueezing can go awry, and the more the input obeys SSWT's assumptions (which aren't too restrictive, and many signals naturally comply), the better the results.

What are SSWT's assumptions? (when will it fail?)
This is a topic of its own (which I may post on later), but briefly, the formulation's as follows. Firstly note that we must somehow restrict what $A(t)$ and $\psi(t)$ can be, else, for example, $A(t)$ can simply cancel out the cosine and become any other function. More precisely, the components are to be such that:

More info in ref 2.

How would it be implemented? There's now Python code, clean & commented. Regardless, worth noting:

For very small CWT coefficients, phase is unstable (just like for DFT), which we work around by zeroing all such coefficients below a given threshold.
For any frequency row/bin $w_l$ in SSWT plane, we reassign from $W_f(a, b)$ based on what's closest to $w_l$ according to $\omega (a, b)$, and for log-scaled CWT we use log-distance.

Summary:
SSWT is a time-frequency analysis tool. CWT extracts the time-frequency information, and synchrosqueezing intelligently reassigns it - providing a sparser, sharper, noise-robust, and partly denoised representation. The success of synchrosqueezing is based in and explanied by its prior; the more the input obeys assumptions, the better the results.

Notes & caveats:

What if $f$ isn't nice in $\omega(a,  b)$ example? Valid question; in practice, the more the function satisfies aforementioned assumptions, the less of a problem this is, as authors demonstrate through various lemmas.
In the SSWT of damped pendulum, I cheated a little by extending signal's time to $(-2, 6)$; this is only to prevent boundary effects, which is a CWT phenomenon that can be remedied; here's directly 0 to 6.
Partial noise cancellation? Indeed; see pg 536 of ref 1.
What's the $a^{-3/2}$ in $(5)$? Synchrosqueezing effectively inverts $W_f$ onto the reassigned plane, using one-integral iCWT.
"Fourier bad?" My earlier comparison is prone to criticism. To be  clear, FT is the most solid and general-purpose basis that we have for a signals framework. But it's not an all-purpose-best; depending on context, other constructions are more meaningful and more useful.

Where to learn more?
The refernced papers are a good source, so are MATLAB's wsst and cwt docs and ssqueezepy's source code. I may also write further Q&A's, which you can be notified of by subbing this thread.

References:

A Nonlinear Squeezing of the CWT Based on Auditory Nerve Models - I. Daubechies, S. Maes. Excellent origin paper with succinct intuitions.
Synchrosqueezed Wavelet Transforms: a tool for Empirical Mode Decomposition - I. Daubechies, J. Lu, H.T. Wu. Good followup paper with examples.
The Synchrosqueezing algorithm for time-varying spectral analysis: robustness properties and new paleoclimate applications - G. Thakur, E. Brevdo, et al. Further exploration of robustness properties and implementation details (including threshold-setting).


How to design a a filter and how to determine the poles and zeros on the basis of following informations?

I guess it's for a homework.

1 - You want to attenuate the 50-Hz component. Zeroes placed at 50 Hz will cancel the 50-Hz component, thus you need to place your zeros at 50 Hz.

2 - You want your 300-Hz signal to be unaffected by the zeros, i.e. you want it to keep the same amplitude and the same phase. If you only added the complex pair of zeroes at 50 Hz, it would increase the amplitude of your 300-Hz signal and affect the phase at 300 Hz also. Therefore, you need to "cancel" the effect of your complex zeroes pair. If you add a complex pole-pair at 49 Hz, then the combined effect of the pole pair at 49 Hz and the zero pair at 50 Hz will almost cancel each other at 300 Hz.

Intuitively, poles and zeros are opposites. At 300 Hz, the distance to the zero pair and the pole pair will be approximately the same, thus they will almost cancel each other out.

3 - I told you at which frequency to place your pair of poles and pair of zeroes, but now you need to select a damping factor for both pairs. This issue is a bit more complicated, but you can safely put the pair of zeroes with a damping factor of 0, i.e. on the complex circle. However, you need to put your pair of poles inside the unit circle, but close to the unit circle since poles outside the unit circle are unstable. Zeroes outside the unit circle are stable, but are not desirable.

4 - If you use this method, you will get a gain of approximately 1 at 300 Hz, but not exactly 1... Since I think this is a homework, this is good enough. In real-life, you would do best to work with the transfer function directly, whether in the s-domain or z-domain. https://en.wikipedia.org/wiki/Band-stop_filter


Frequency modulation of $m(t)$ in MATLAB

Integration was incorrect. I was integrating it with $dx= 1$. Doing m(j)*1/Fsoperation in the for loop, I have made $dx= 1/Fs$ now.

As can be seen, frequency changes for different amplitudes now.

Moreover, resultis a constant, it should be changed toresult2 so that it gives result with respect to time $t$.


Bit depth is typically applied to a fixed point representation. 16-bit means the number can go from -2^15 to 2^15-1 so from -32768 to +32767. For 8 bits that would be -128 to +127.

If you normalize to a maximum of 1, your least significant bit (LSB) value is simply the inverse, 1/32678 = 3.0518e-05 for 16 bit or 1/128 = 0.0078125 for 8 bit

Floating point works completely differently. A 32-bit floating point number has a 24 bit mantissa and an 8 bit exponent. So it's basically a 24-bit integer number combined with an 8 bit scale factor that gets "optimized" and recalculated on every floating point operation. See https://en.wikipedia.org/wiki/Floating-point_arithmetic


How to keep the same power (sonority) between pure sines?

From the comments on the question I gather that you want the individual sine waves to appear to be at the same volume. Unfortunately, there is no single EQ setting you could apply to achieve this. As for the reasons, let's look at the equal-loudness contours:

Source

What these contours tell us is that for a specific level of phones there exists an EQ setting that would make your sine waves sound at equal loudness. Unfortunately, the human perception of loudness not only depends on frequency but also on SPL—what that means is that your EQ settings that you meticulously set up for a certain level of phones won't work any more if you change the overall gain of your playback system. To compensate for the SPL-dependency of the human hearing you would need an adaptive EQ that takes into consideration the SPL of the signal at 1000 Hz (or any other reference frequency) by measureing it with a sound level meter and setting the EQ parameters accordingly.


If the noise is a constant background/remains in the same place, you could take various photos in the absence of any light to create an image for background subtraction. My image processing expertise is limited to grayscale but i would imagine the same could be done for RGB. Produce several background images, from that make an average/low pass filter. Then in your desired programming language do some simple matrix manipulation.


Controlling a signals amplitude in the frequency domain

You can use Parseval's theorem for DFT.
$$
\sum_{n=0}^{N-1} |x[n]|^2  = \frac{1}{N} \sum_{к=0}^{N-1} |X[n]|^2
$$
Where $$x[n]$$ - n-th signal sample,
$$X[n]$$ - n-th value of the DFT of the signal


You are on the right track.  You can use decimation and interpolation (or resampling, which is a combination of both decimation and interpolation) to change the sample rate.  For instance, if you start with 48 kHz samples you can decimate by 6 to get down to 8 kHz, do some processing, then interpolate by 2 to get to 16 kHz, do some processing, then interpolate by 3 to get back to 48 kHz.

"Lining up the processes" is not in and of itself difficult.  You can decimate and interpolate whenever you want.  Doing it in a way that makes sense and is useful, on the other hand, requires knowing more about what you start out with and what you want to accomplish.


1)Are these notations referring to delay/shift in input and impulse response?

Yes

  2)Also it mentions that **Periodic or circular convolution is also called as fast convolution

That's a bit of a misrepresentation. For a non-trivial length of signal and/or impulse response the fastest way to implement convolution is to multiply in the frequency domain. That does indeed implement circular convolution. That's typically NOT what you actually want. If you want to leverage frequency domain multiplication for linear convolution, you need to add the right amount of zero padding and/or overlap handling to your algorithm. Read up on "Overlap Add" or "Overlap Save" algorithm.


Is it correct to call a Savitsky-Golay filter of degree 0 a simple moving average?

Yes you can consider a zeroth (or first) order SG filter as a moving average filter. Below MATLAB / Octave code computes the impulse response of a SG filter of order $N$ and length $2M+1$ :
% Savitzky-Golay Filter
%
clc; clear all; close all;

N = 0;                      % a0,a1,...,aN : Nth order polynomial
M = 3;                      % x[-M],...,x[M] : 2M + 1 data

A = zeros(2*M+1,N+1);
for n = -M:M
    A(n+M+1,:) = n.^[0:N];
end

H = (A'*A)^(-1)* A';        % LSE fit matrix

h = H(1,:);                 % S-G filter impulse response (non-causal symmetric FIR)

figure,subplot(2,1,1)
stem([-M:M],h);
title(['Impulse response h[n] of Savitzky-Golay filter of order N = ' num2str(N), ' and window size 2M+1 =  ' , num2str(2*M+1)]);

subplot(2,1,2)
plot(linspace(-1,1,1024), abs(fftshift(fft(h,1024))));
title('Frequency response magnitude of h[n]');

figure,plot(linspace(-1,1,1024), 20*log10(abs(fftshift(fft(h,1024)))));
title('Frequency response magnitude of h[n]');

The impulse response for $N=0$ is :

as can be seen, it's a moving average (constant) impulse response.


Is it possible to Denormalize a matrix of complex numbers?

Obtain the Bitmap-image of a magnitude plot of an FFT.

  De-normalize the Bitmap-image to obtain the original FFT of the image.

This can't work. Preserving only the magnitude throws away half of the information contained in the complex values, and thus, you can't get the original FFT back, ever.


It is all about defining the support of the DFT.
The DFT is periodic with period of $ 2 \pi $.
Hence on any support of the length $ 2 \pi $ you see all the data.

On default, when using MATLAB's fft() the result is on the range $ \left[ 0, 2 \pi \right] $, namely the Low Frequencies are at the beginning of the output vector.
Using fftshift() the DFT is shifted to match $ \left[ -\pi, \pi \right] $ support which means the DC Frequency is centered.

The same logic works on 2D DFT.


Use of complex numbers to calculate Doppler effect?

Bottom Line First

What you suggest will only work if your signal is at baseband and not at a carrier, which is usually the best way to do any simulation work regardless. If your signal is not at baseband first multiply by $e^{-j2\pi f_c t}$ where $f_c$ is the carrier frequency and low pass filter out the higher image (which will now be centered at $-2f_c$) to create the baseband equivalent waveform. With this you can simply multiply by $e^{\pm j2\pi f_d t}$ to introduce a Doppler shift. If you need to have your waveform at a carrier frequency in simulation, then you can shift the Doppler-shifted signal back to the carrier after the Doppler shift operation at baseband described above by multiplying the resulting shifted waveform by $e^{+j2\pi f_c t}$ and then taking the real part if a real signal is desired. Usually however ALL simulation work can be done at baseband with a complex envelope that is equivalent to what is happening at any carrier frequency and much simpler to model.

(Note if the concept of baseband vs passband mentioned above is confusing, see this post that goes into further details on that: modulate complex time domian signal to a carrier frequency)

DETAILS

Simply multiplying any signal by a complex rotation will not work for simulating a Doppler shift in all cases; that applies to when the signal is represented by its complex equivalent at baseband, or if you are only concerned with the complex single-sided equivalent of your signal.

This is because the Doppler will shift a signal either higher in frequency or lower in frequency depending on the relative direction (as you have properly shown), and a real signal is comprised of both positive and negative frequency components, each described by $e^{j2 \pi f t}.$  and multiplying by a single complex rotation will shift EVERYTHING in the same direction.

Multiplying a single complex tone at carrier frequency $f_c$ ( $e^{j2\pi f_c t}$)  by a Doppler frequency $f_d$ ($e^{j2\pi f_d t}$) is simply summing the exponents which results in the expected frequency shift:

$$e^{j2\pi f_c t}e^{\pm j2\pi f_d t} = e^{j2\pi (f_c \pm f_d) t}$$

(And the same would occur is we were dealing with a modulated signal as a band of frequencies around $+f_c$.)

However, look what happens if the signal of interest is a real cosine wave which is comprised of TWO complex conjugate symmetric complex tones:

$\cos(2 \pi f_c t) = \frac{1}{2}e^{j2\pi f_c t}+ \frac{1}{2}e^{-j2\pi f_c t}$

$\cos(2 \pi f_c t)e^{j2\pi f_d t} = \frac{1}{2}e^{j2\pi (f_c+f_d) t}+ \frac{1}{2}e^{-j2\pi (f_c-f_d) t}$

Thus we see that for a cosine wave (or any real signal) which has a positive frequency spectrum and a complex conjugate negative frequency spectrum (same magnitude, opposite phase), all components will shift in the same direction-- NOT what occurs for a Doppler shift.

I would highly recommend doing all such simulations at baseband (remove the carrier frequency and replace with 0 Hz and deal with complex signals expressed as $e^{j\omega t}$, not sines and cosines) in which case you can do such an operation which is significantly simpler. However if it absolutely needed to simulate with a carrier, to shift a real signal at a carrier frequency by a Doppler frequency, you need to first take the Hilbert transform of your signal, which when summed in quadrature with the original signal, effectively removes all the negative frequencies (in practical implementation, this is a quadrature splitter and the resultant complex signal is implemented as two separate real channels), and then with this new complex signal, you can multiply by $e^{\pm j2 \pi f_d t}$ as suggested to shift with the sign chosen dependent on the direction desired. After multiplying you can take the real of the result to re-establish your original real signal with the Doppler shift induced. Note this is exactly how a single-sideband modulator is structured as demonstrated by the slides pasted below:

Note, the operation above is showing how to shift a single-sideband signal at $\omega_{\Delta}$ away from the carrier down to baseband, but the same approach can be used to show how to shift a signal at $\cos(\omega_c)$ to $\cos(\omega_c+\omega_{\Delta})$ where $\omega_{\Delta}$ represents the Doppler shift. Also to note in the last diagram, the Hilbert transform will convert the cosine to a sine, and must be added to the sine with a $j$ term to become the desired exponential rotation:

$e^{j\omega_c t}= \cos(\omega_c t)+ j \sin (\omega_c t)$

In the implementation the required j term shown in front of the sine is implied by carrying that path separately (as is the case with having I and Q paths in a digital modem implementation for example); we need two real numbers to represent a complex number, as is done by having the two signal paths. (either as real and imaginary components such as $I + jQ$ or as magnitude and phase components suchs as $A\angle \phi$, which by the way is equivalent to $Ae^{j\phi}$)

The above approach based on a single-sideband modulator gives excellent insight into how to shift a signal at a carrier frequency in one direction (as happens with Doppler shift), and it would be exactly how a hardware implementation could be done at RF to simulate a Doppler shift. For simulation work my favored approach would be to do all simulation at baseband using the "Bottom Line" approach first described.


Does an audio equalizer consist of a perfect-reconstruction filterbank?

A filter bank system is generally composed of:

an analysis filter bank, that split an input signal into components with filters and may reduce their rate
a synthesis filter bank (FB), that takes the components, potentially increase their rate, and feed them into filter elements

so that the novel outputs can be combined into a signal. It is perfect reconstruction (PR) when, for every input, the output is a potentially delayed or scaled version of the input. Filters can be nonlinear, the rate reduction can be different for each component. Min-Max or morphological wavelets are examples of the above. There are even hybrid filter banks, where the analysis (A) or synthesis (S) FB are either analog or digital. Their design is quite involved, so often A and S are linear, with upsampling and downsampling operators, and an output obtained by summation, as explained by @Marcus Müller. However, the filters may have different bandwidths and the subsampling can be different for each filter.

Many signal processing systems, such as (most) equalizers, are made of an analysis FB, a processing unit (often nonlinear) and a synthesis FB (which could be a mere summation). The non-perfection arises, for PR FB, because of the middle filtering component.

Assuming an equalizer has linear A and S filters, the non PR aspect may arise from two aspects. First, what are the filters? It seems to me that the traditional filters (first or second order) do not allow easy PR. However, with sufficient frequency overlap, they together can form a near PR filter banks, with an error below the hearing threshold.
The second aspect is the processing. If you suspect saturation, or nonlinear phase alteration, then PR is not achieved in  general.

Phase distortion is sometimes considered harmless for audio in some   textbooks. However, picky sound engineers are cautious about audible phase distortion. The article Plugging into EQ discusses those issue, around FIR or IIR filters. The first ones can have linear phase, but have less sharp frequency transitions, and require a handful of computations. The second ones do not have linear phase in general. It also provides a comparative evaluation of several EQ.

So, as Markus, I'll say no in general, with a modest yes (in cases) if below the hearing sensitivity.


Filtered Signal vs File Compression Paradox

+1 on very interesting and insightful experiment.

Some thoughts:

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


Independent Subspace Analysis and Topographic ICA

Read this paper about ISA: http://www.cs.helsinki.fi/u/ahyvarin/papers/ESANN06.pdf

and this one about topographic ICA: http://www.cs.helsinki.fi/u/ahyvarin/papers/NC01_TICA.pdf

Hyvärinen is the author of the FastICA algorithm.

PS.

By the way, you might want to look at ICASSO (http://research.ics.aalto.fi/ica/icasso/), if you are not already familiar with it.


Accurately measuring relative distance between a set of fiducials (Augmented reality application)

You can use structure from motion type of algorithm to estimate the camera pose from the environment, not from the markers and then fuse this camera pose with the marker poses in order to detect the locations of the markers accurately. Knowing the extrinsic pose of your camera (by SFM), you can triangulate all the 3D positions.

For pose estimation, 5-point methods are usually more accurate than 8-point algorithm.

Presumably, you should do further bundle adjustment so that the overall accuracy increases.


Perhaps the standard is Ogata's "Modern Control Engineering".

I also recommend Franklin/Powell/Emami-Naeini's "Feedback Control of Dynamic Systems".

For a more basic introduction you might consider Nise's "Control Systems Engineering", but I'm not a fan of his writing style.


Frequency response of $\mathrm{sinc}[n]$

First of all, there is a mistake in that image: there should be an $n$ in the denominator instead of a $t$.

Check that your signal can be expressed as

$$h[n]=2\frac{\sin\left(\frac{\pi}{2}n\right)}{\pi n}=\frac{\sin\left(\frac{\pi}{2}n\right)}{\frac{\pi}{2}n} \tag{1}$$.

Using Euler's formula we can express that sine as the sum of two exponentials:

$$h[n]=\frac{1}{\pi/2}\frac12\frac{e^{jn\pi/2}-e^{-jn\pi/2}}{jn}=\frac{1}{\pi}\frac{e^{jn\pi/2}-e^{-jn\pi/2}}{jn}$$

Let me multiply and divide that expression by $2$, which doesn't change anything but will be useful later:

$$h[n]=\frac{1}{2\pi}2\frac{e^{jn\pi/2}-e^{-jn\pi/2}}{jn} \tag{2}$$

Notice that there exists an integral that returns that exact function of $n$ and, on top of it, that integral represents the inverse Fourier transform of some function:

$$\frac{1}{2\pi}\int\limits_{-\pi/2}^{\pi/2}2e^{j\omega n} \ \mathrm{d}\omega=\mathcal{F}^{-1}[H(e^{j\omega})]=h[n] \tag{3}$$

where $H(e^{j\omega})$ is periodic in $2\pi$ and, for each period:

$$H(e^{j\omega}) =
\left\{
 \begin{array}{ll}
  2  & \mbox{if } |\omega| \leq \pi/2 \\
  0 &   \ \mathrm{otherwise}
 \end{array}
\right.$$

So the Fourier transform of $h[n]$ is $H(e^{j\omega})$, because the inverse Fourier transform of the latter is the former. And that's the proof.

However, most of the times one does not do all these calculations and just looks up the transform of $\mathrm(1)$. (In the link's table, look up for your signal and replace $W=2$. Then use the linearity of the Fourier transform to get to the result you are looking for.)


Why low-pass FIR Equiripple filter is designed with 3 dB of ripple in the pass-band and at least 60 dB of attenuation in the stop-band?

Stop band attenuation is a good thing right? After all it's the primary reason why you want a frequency selective filter: stop the unwanted frequencies. So a 60dB attenuation in the stop band of a filter is considered to be the norm in many typical applications. But not always so, as the comments indicate; you must always refer to your application requirements for determining how much attenuation is enough (or necessary) for your case.

That being said, there exist tradeoffs in filter design, just like in many other design procedures. And the most fundamental of those tradeoff parameters, for digital FIR filter design, are namely the passband ripple, stopband attenuation, transition band width and the resulting filter's length.

Among these parameters transition band width and stop band attenuation affects the selectivity quality of your filter. But a short transition band with a large stop band attenuation means a (very) high order (long impulse response) FIR filter, which creates opposition due to computational burden.

To overcome this computational burden, due to narrow transition band, once can consider relaxing the filter's definition in, for example, the pass-band region, so that instead of a monotonic or an un-equal ripple filter, you distribute the ripples equally by allowing as large as possible but fixed amplitude ripples in every where through the passband, so that the FIR filter length will be kept short enough.

3dB passband ripple can be considered a reasonable choice for some applications in which pass band distortions due to those ripples may be tolerable.


Partial derivatives in discrete domains are computed using differences (finite difference in more proper terms). There are many methods to do that such as central finite difference or forward, backward and so on.

The convolution operation you mentioned computes forward distances in horizontal direction and results in x-derivatives (namely the $x$ component of the gradient).

To also compute in $y$ direction you might do:

fx = conv2(im1,0.25* [-1 -1; 1 1])

0.25 is there just to scale the values to appropriate range (weights each pixel the same).

For more information on finite differencing you might check:
http://en.wikipedia.org/wiki/Finite_difference#Forward.2C_backward.2C_and_central_differences

Convolution is a local-masking operation where at each pixel you compute the dot product of the kernel and matrix part that corresponds to the kernel. Then, kernel is shifted. I won't go to the details of convolution as it is already explaind in http://en.wikipedia.org/wiki/Convolution. However, let's say we have a matrix

$
   I=
  \left[ {\begin{array}{cc}
   1 & 2 & 4 \\
   6 & 14 & 17\\
   16 & 41 & 18\\
  \end{array} } \right]
$

If we convolve this with the kernel mentioned and discard for now the boundary conditions, we end up with

$
   G=
  \left[ {\begin{array}{cc}
   -2.25 &  -1.25 \\
   -8.25 &  5.0 \\
  \end{array} } \right]
$

(this is part in the middle of the result).

So, you should notice that basically, per each element, we compute a sum of 2 neighbors vertically on the left, sum of 2 neighbors vertically on the right and subtract these intermediate results. As subtraction is done horizontally, it computes the horizontal derivatives. This is what convolution intrinsically does.


Frequency warping when integrators are replaced with backward-euler and forward-euler integration

To the extent this helps, here are some interesting magnitude and phase plots of the backward Euler (here using the method of impulse invariance mapping which for 1/s has the same result as backward Euler but as Matt L points out in the comments is not generally the case) vs Tustin (method of Bilinear Transform mapping) integrators. Notice the Euler on the left has dominant phase distortion while the Tustin on the right has more dominant amplitude distortion (no phase distortion from an analog integrator at all!).

The analog frequency response for a pure integrator as $1/s$ is -20 dB/decade in magnitude ($1/f$) and -90° for phase.


Have you ever seen logarithmic burn-in on a Kalman filter before?

The authors cite this work in their response to the review:

  Zhao H, Lu L. Adaptive recursive algorithm with
  logarithmic transformation for nonlinear system identification in
  α-stable noise [J]. Digital Signal Processing, 2015:S1051200415002535.

They take their inspiration from equation (14):

  $J_p(n)=\sum_{i=1}^n \lambda^{n-1} \cdot \log^p (1+|e(i)|)$

It is notable that $e(i)$ is (as one might expect) the error in the signal, and not linearly increasing with the number of steps since the filter was initialized.


Which Machine or deep learning algorithm is appropriate of this issue?

Your question can be interpreted in many ways. First, convolution with a number is not a common operation. Convolution is defined by:

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

Assuming a constant signal $g(t)=a\in \mathbb{R}$ this is reduced to:
$$(f*g)=\int_{-\infty}^{\infty}af(\tau)d\tau=a\int_{-\infty}^{\infty}f(\tau)d\tau$$

I am assuming, therefore, that you meant multiplication and not convolution. Correct me otherwise and I will edit the answer accordingly.

Secondly, assume for your example the given parameters:

$$Y=[1,2,3,4],J=\left[
 \begin{matrix}
  1 & 2 & 3 & 4\\
  1 & 2 & 3 & 4\\
  0 & 0 & 0 & 0\\
  0 & 0 & 0 & 0
 \end{matrix}
\right], X=\left[\begin{matrix}4 \\ 4\\0\\0\end{matrix}\right]$$

It is easy to see that a couple of use-cases could give us this result. For example, it could be that the chosen $y$ was $y_4=1$ and the chosen column was the most right column $j_3$. On the other hand, the combinations $\{y_1,j_1\},\{y_3,j_0\}$ are also plausible. I have purposedly gave a very bad, linearly dependent, matrix $J$ to illustrate that the choice of $J$ in the design process or formulation of the problem is of major interest and is fundamental for answering you question. The following example illustrates the complete opposite:

$$Y=[1,2,3,4],J=\left[
 \begin{matrix}
  1 & 0 & 0 & 0\\
  0 & 1 & 0 & 0\\
  0 & 0 & 1 & 0\\
  0 & 0 & 0 & 1
 \end{matrix}
\right], X=\left[\begin{matrix}x_0 \\ x_1\\x_2\\x_3\end{matrix}\right]$$

Obviously, $x_0=y*j_0$,$x_1=y*j_1$, and so on. Here, it is easy to deduce $y$ and only a small number of $X$ are even possible.

Generally, if we call $J$ a dictionary, you are looking at something called Orthogonal Matching Pursuit (OMP) problem, which is not a machine or deep learning problem. Otherwise, you can try and use a deep learning model to try and learn the dictionary, though this would be unrecommended as such a solution is not compatible with such a problem in my opinion.

Hope this answers your question. If not, or I misinterpreted your question, feel free to comment and I will try to supply more information if I can.


Measure object planar displacement on slight moving or vibrating plane

...turns out that the plane might move/vibrate slightly (~ +-5mm), would I still getting good results?

I do not see why not. And this is making reasonable assumptions about the equipment.

You have a monocular camera, tracking motion from a viewpoint that is vertical to the direction of motion.

Now, if the plane moves purely in the Z direction (that is, on an axis that is vertical to the lens), then the only thing that will change is the size of the object.

The problem this might create is in the Signal To Noise ratio of the tracking. In other words, if you are tracking a marker against a background and the marker's size "pulsates", then the tracking estimate of the frames where the target is relatively larger will be better (within reason).

By how much will your marker / object size change, depends on the focal length of your lens.

If the plane does not move purely in the Z direction but sort of oscillates back and forth as if it was a pendulum, then it introduces an additional perspective error to the measurement that is proportional to $\sin(\theta)$ where $\theta$ is the angle of oscillation.

In other words, a monocular camera imaging a vertical line segment cannot know if the imaged size is true or the projection of a much larger vertical line segment that stands at some angle $\theta$ towards the lens.

In general however, the "magnification error" (first case) grows smaller with the length of the lens (lower in long lenses) and vice versa. There of course, you have to balance the field of view too so that it includes the variation you are trying to measure.

If the equipment you are using has specific characteristics regarding focus and image formation, you are going to have to take those into account of course.

Hope this helps.


"Fourier Transform can localize signals in frequency domain, but not in time domain." -- What does it mean in layman's terms?

In the Fourier transform, the basis functions are complex exponentials. These functions are perfectly localized in the frequency domain, i.e., they exist at one frequency, but they have no time localization because of their infinite duration. The localization of a function depends on its spread in time and frequency. A complex exponential has zero spread in frequency, but it has infinite spread in the time domain. Consequently, any signal analyzed by the Fourier transform is only localized in frequency (e.g., frequencies of sinusoidal components can be identified), but not at all in time (e.g., temporal changes of signal properties cannot be localized).


Comparing band-pass filters and convolutional kernels?

CNNs sort of irrelevant here because they are applicable when you have some invariant transform - shift or rotation. The idea of CNN is that weights are shared for different data points, not that weights are different. You are thinking about simple DNN.

CNNs are also applied in speech recognition, but for a bit different task. For example, you can consider multiple frames jointly and tie filter shape for different frequencies. This way system learns Gabor-like filters. See Robust CNN-based Speech Recognition With Gabor Filter Kernels.

You can try to train your own kernels from raw data, they will eventually resemble bandpass filters. See Acoustic Modeling with Deep Neural Networks Using Raw Time Signal for LVCSR. This gives a marginal improvement.

You can also start from default bandpass filters as initial values and see if network learns an improved version.

You can also restrict network architecture to consider only the region around central frequency and just learn the weights in that region. Such network will unlikely improve anything.


Instantaneous amplitude is usually done in conjunction with analytic signals and the Hilbert transform. See this answer on this site.

  For signal $s(t)$, given its Hilbert Transform $\hat{s}(t)$ it is defined as a composition:

  $$s_A(t)=s(t)+j\hat{s}(t) $$

  The Analytic Signal that we obtain is complex valued, therefore we can express it in exponential notation:

  $$s_A(t)=A(t)e^{j\psi(t)}$$

  where:

  $A(t)$ is the instantaneous amplitude (envelope)

  $\psi(t)$ is the instantaneous phase.


Change all data types to "float"
Split the filter into two second order sections
If you can't afford the floating point, you need to learn fixed point math, which will require significant effort


Gain at given frequency from $z$-plane zero-pole plot. Two methods gives different results

I think you lost the square root somewhere in your graphical method: Considering the nominator (product of lengths of vectors from every zero given point), I end up with

$$N=2*2$$

i.e. you have two vectors of length 2 (two zeros at z=-1).

Then, for the denominator you have

$$
D = \sqrt{(A\sin\theta)^2+(A\cos\theta-1)^2}^2=A^2-2A\cos\theta+1
$$

And finally $G=N/D$ yields the same result as your analytical solution.


What is the effect of the natural logarithm in the frequency domain?

First note that when you use a logarithmic function you shall avoid negative arguments if it's output should be real valued. Then consider the following relation:
$$ y[n] = \ln( 1 + x[n] )$$
where the input $x[n]$ is conditioned such that $ 0 < 1 + x[n] \leq2$ is maintained. Then a Maclaurin series expansion of this expression will be

$$ y[n] = \ln(1 + x[n]) = \sum_{k=1}^{\infty} \frac{(-1)^{k-1}}{k} x^k[n] .$$

For very small $x[n]$ the following approximation can also be taken:
$$ y[n] \approx x[n] - x^2[n]/2 $$

Applying DTFT on this last line yields:

$$ Y(e^{j\omega}) \approx X(e^{j\omega}) - \frac{1/2}{2\pi} \left( X(e^{j\omega}) \star X(e^{j\omega}) \right) $$

The spectrum, $Y(e^{j\omega})$, of the logarithm of $1+x[n]$ will include multiple convolutions of $X(e^{j\omega})$ with itself. The relation becomes exact (with added terms which were ignored) as $k$ approaches infinity, showing that the spectrum will be broadened.

Indeed for any nonzero bandwidth input $x[n]$, the exact output $y[n]$ would have infinite bandwidth (or equivalently unavoidable aliasing) and no finite upsampling is enough to properly represent it. However in practice as $|x[n]|<1$ is maintained, its high powers, $x^k$ will tend to go zero quickly and can be ignored compared to some suitable threshold...


Finding the Amplitude of a Sinusoid in Noise

If you know you have a single sinusoid in Gaussian noise, then a least squares type parametric estimator might be one method of estimating the amplitude. See  DSP Related, specifically:


Frequent update of the Kalman Filter covariance matrices

The term "Kalman Filter"  tends to be used in genetically and there are a number of variations suitable for different applications.

The documentation in OpenCV for the Kalman Filter describes a strictly Linear type which they say can be modified to be an Extended Kalman Filter (EKF). Your Jacobian implies that your problem is nonlinear, and the OpenCV version is exclusively a linearization.  Many problems are solvable using this approach, but it works better, the smaller your time step (in other terms updated more often). You can also inflate your state covariance.

Another approach, where the state transition $\bf{x_k} \; \Rightarrow \;\bf{x_{k+1}}$ is nonlinear, and using an integration technique, typically Euler Integration is used on the State Evolution equation. The nonlinear state covariance update uses the Jacobian linearization. It doesn't look like you can use the OpenCV KF, for this kind of EKF.

One issue with the KF in general is choosing the number of states you use in your filter. A very common problem in applications is when there is a large change in state, like turns, the Filter exhibits a lot of lag.  If an additional state was present to account for the change (acceleration in Kinematic models) the lag would be much lower, but the the higher dimension can make the Filter more sensitive to small disturbances.  In some cases, this isn't a problem and in others there is the Interacting Multiple Model (IMM) or as a friend of mine called it "Dueling Kalman Filters". The OpenCV doesn't look like it provides the innovations, which probably means that IMM is not a candidate with an all linearized EKF.

The essential problem when you have a jump and don't account for the possibility of a jump, is that the Filter places more trust on the predicted state and the change looks like noise.  It wants (as an anthropomorphic analogy) to keep the state level and steady.  In many cases, the error accumulates slow enough for the Filter to adjust to the jump.  In other cases, track is lost. State Covariance inflation is one way to have the Filter weigh the new measurements over the predictions based on the old state.

One way to avoid linearization problems and approximate integration is to used the Unscented KF, which is often a good solution. The state dimension problem can creep in with this filter a well.

The sensitivity to jump lag manifests itself differently to different forms of the KF and update rate dependency is present in all forms.  Kalman Filters need a lot of tuning and testing in any case.

Step size or update rate implicitly effects time integration error. Small steps, less error. Nonlinear integration, less error.


Do Causal Discrete-time systems have proper transfer functions?

Note that a stable and causal continuous-time transfer function does not need to be strictly proper but only proper, i.e. the degree of the numerator does not exceed the degree of the denominator, but numerator and denominator degree can be equal. E.g.

$$H(s)=\frac{as^2+bs+c}{s^2+ds+e}$$

can represent a causal and stable system, as long as its poles are in the left half of the complex $s$-plane.

For discrete-time systems the same is true. A transfer function for which the degree of the numerator is greater than the degree of the denominator has at least one pole at infinity. Since causal and stable discrete-time systems must have all their poles inside the unit circle, such a system can't be causal and stable. However, as is the case with continuous-time systems, equal degrees of numerator and denominator are possible.

Assume a causal and stable discrete-time system is described by the following linear difference equation with constant coefficients:

$$y[n]+a_1y[n-1]+\ldots+a_Ny[n-N]=b_0x[n]+b_1x[n-1]+\ldots+b_Mx[n-M]\tag{1}$$

The corresponding transfer function is

$$H(z)=\frac{\sum_{m=0}^{M}b_mz^{-m}}{1+\sum_{n=1}^Na_nz^{-n}}\tag{2}$$

Note that there are of course no restrictions on the choices of $M$ and $N$. Nevertheless, the transfer function $H(z)$ is always proper. Assuming $M>N$ gives

$$H(z)=\frac{b_0z^M+b_1z^{M-1}+\ldots+b_M}{z^M+a_1z^{M-1}+\ldots+a_Nz^{M-N}}\tag{3}$$

Obviously, the transfer function given by (3) is proper. $M>N$ just means that the last few coefficients of the denominator polynomial are zero. A very similar argument holds for $M<N$.


Polyphase Filter decomposition. It is not working

The following is a working code that uses 32-component polyphase decomposition of the associated 32-channel anslysis and synthesis filterbanks. As I have already commented, the speed gain is not dramatic in this cae due to short signal and filter lengths. However further architectural improvements as well as coding optimizations can provide better results.

% S0 - Load the prototype lowpass filter impulse response h0[n]:
% --------------------------------------------------------------
load h2.mat;         % h[n] is the prototype lowpass filter of length 512
L = length(h);

% S1 - Create the 32 x 512 analysis filter-bank hha[k,n] by cosine modulation from protoype :
% -----------------------------------------------------------------------------------------
numbands = 32;                 % number of banks (channels)
n=0:L-1;

hha=zeros(numbands,L);         % bank of filters hha[k,n] = 32 x 512 array.
for k=0:1:numbands-1
   hha(k+1,:) = h.*cos( ( (2*k+1)*pi*(n-16) ) / (2*numbands) );
end

% S2 - Create the 32-polyphase components hhap[k,m,n] , for each one of 32 analysis filters hha[k,n]:
% ---------------------------------------------------------------------------------------------------
numpoly = numbands;             % polyphase component number = decimation ratio = number of channels
hhap = zeros(numbands,numpoly, L/numpoly);  % hhap = 32 x 32 x 512/32 , 3D ANALYSIS filter bank array

M = numpoly;                    % polyphase system decimation ratio
for k=1:numbands
    for m = 1:numpoly
        hhap(k,m,:) = hha(k,m:M:end);       % create the m-th polyphase component of k-th channel filter
    end
end

% S3 - Design the 32 x 512  synthesis (complementary) filter bank :
% -----------------------------------------------------------------
numbands = 32;                  % number of banks
n=0:L-1;
hhs = zeros(numbands,L);        % bank of filters
for k=0:1:numbands-1
   hhs(k+1,:) = h.*cos( ( (2*k+1)*pi*(n+16) ) / (2*numbands) );
end

% S4 - Obtain the 32-polyphase components hhsp[k,m,n] , for each one of 32 synthesis filters hhs[k,n]:
% ----------------------------------------------------------------------------------------------------
numpoly = numbands;             % polyphase component number = interpolation ratio = number of channels
hhsp = zeros(numbands,numpoly, L/numpoly);  % hhap = 32 x 32 x 512/32 , 3D ANALYSIS filter bank array
M = numpoly;                    % polyphase system decimation ratio
for k=1:numbands
    for m = 1:numpoly
        hhsp(k,m,:) = hhs(k,m:M:end);       % create the m-th polyphase component of k-th channel filter
    end
end

% S5 - Generate the test input signal
% -----------------------------------
N = 2*1024;
wav_in = cos(0.01791*pi*[0:N-1]);        % pure sine tone

% S6 - Apply test signal to the filterbank: ANALYSIS STAGE :
% -----------------------------------------------------------
yyd = zeros( numbands, floor(N/numbands));   % decimated outputs..
M = numbands;
for k=1:1:numbands
    temp = conv([wav_in(1:M:end),0] , hhap(k,1,:));
    for m=2:M
        temp = temp + conv([0,wav_in(M-m+2:M:end)],hhap(k,m,:));
    end
    yyd(k,:) = temp(L/(2*M)+1 : L/(2*M)+N/numbands);
end

% S7 - Apply SYNTHESIS filterbank on the decimated signal :
% ---------------------------------------------------------
ys = zeros(1, N);

for k=1:numbands
    temp = zeros(1, N+L-1);
    for m = 1:numpoly
        temp(m:numbands:end-31) = conv( yyd(k,:) , hhsp(k,m,:) );
    end

    ys = ys + temp(L/2+1:L/2+N);
end
ys = numbands*ys;

% SX - DISPLAY RESULTS:
% ---------------------
L = length(h);
figure,subplot(2,1,1)
stem([0:L-1],h);title('The Prototype Lowpass Filter');
subplot(2,1,2)
plot(linspace(-1,1,4*L),20*log10(abs(fftshift(fft(h,4*L)))));
grid on;

figure
plot(linspace(-1,1,4*L),20*log10(abs(fftshift(fft(hha(1,:),4*L)))));
hold on
for k=2:numbands
    plot(linspace(-1,1,4*L),20*log10(abs(fftshift(fft(hha(k,:),4*L)))));
end
title('32 CHANNEL FILTERBANK');

figure,subplot(2,1,1)
plot(wav_in);title('input signal')
subplot(2,1,2)
plot(linspace(-1,1,4*N),20*log10(abs(fftshift(fft(wav_in,4*N)))));

figure,subplot(2,1,1)
plot(ys);title('Synthesized Back');
subplot(2,1,2)
plot(linspace(-1,1,4*N),20*log10(abs(fftshift(fft(ys,4*N)))));


After the AWGN block you need to divide with the channel coefficient value before you apply the BPSK demodulator.

If you click on the Rayleigh block it gives the option to enable a "Complex Gain Path Port". You can use this output to accomplish the division mentioned above.

Instead of division you can also multiply the output of the AWGN block with the complex-conjugate of the channel gain.


What is DC level in white gaussian noise?

DC level in signal processing refers to the average or the mean value of a signal. So a zero-mean signal will have an average value of zero over its domain of definition.

When this concept is extended to the stochastic processes, we shall distinguish between time average and ensemble averages. A stochastic process can be viewed as an ensemble of time functions and a sequence of random variables. When a stochastic process is assumed to be ergodic, then its time averages over a sample function equal its ensemble averages over the random variables. In such a case the zero-mean random process creates time functions which have zero averages, hence zero DC levels.


What you have to show is that
$$W_3^2+W_3+1=0\tag{1}$$
with $W_3=e^{-j2\pi/3}$. It's a straightforward exercise to prove $(1)$.
In general we have
$$\sum_{k=0}^{N-1}W_N^k=0,\qquad N>1\tag{2}$$
with $W_N=e^{-j2\pi/N}$. Eq. $(2)$ can be shown by using the formula for the geometric sum.


Get phase information from bursty time series with amplitude variability and sharp edges

I recently responded to a question elsewhere that shared some similarity with this one, so parts of that are re-used here but you can refer to it for the general idea anyway.

Since most of the information you are trying to capture exists in the phase spectra, you don't really "care" about the amplitude information as long as it doesn't produce any distortions.

Therefore, I would like to clarify that my response at this point is focusing on removing the amplitude modulation...carefully.

The main idea is rather simple: Detect the local maxima and reshape each pulse proportionally so that the amplitude modulation vanishes but the timing (and shape) of the pulses is preserved. The following code snippet contains a lot of inline comments on how this is done. It is in Python but it's not incredibly hard to port it to another language (e.g. MATLAB (?)).

from pandas import DataFrame
from scipy.interpolate import interp1d
import pylab
import numpy

def getEnvelopeModels(aTimeSeries):
    '''Detects peaks and returns a suitable interpolation model
    Parameters:
        aTimeSeries: An Nx2 DataFrame that contains the timings and the actual samples in respective columns (col 0, col 1)
    '''
    #Detect peaks and mark their location
    #Get the derivative of pulse train to get rid of "flat tops"
    df = numpy.diff(aTimeSeries[1])
    #Append first data point
    res = [(aTimeSeries[0].iloc[0],aTimeSeries[1].iloc[0])]
    #Simple threshold detector on the derivative to recover peaks and their location
    for k in xrange(0,len(df)):
        if df[k]>10:
            res.append((Q[0].iloc[k+1],Q[1].iloc[k+1]))
    #Append last data point
    res.append((aTimeSeries[0].iloc[-1],aTimeSeries[1].iloc[-1]))
    #Convert to array
    res = numpy.array(res)
    #Fit a model, here it is a 'zero' model which basically maintains
    #the y value until the next (x,y) interval. This is how the pulse width is preserved.
    u_p = interp1d(res[:,0],res[:,1], kind = 'zero',bounds_error = False, fill_value=0.0)
    #Return the model and the peak data
    return (u_p,res)

#This "if" is required here so that Python doesn't execute
#the whole script even one attempts to "borrow" the getEnvelopeModels
#from another script with a simple 'from myFile import getEnvelopeModels'
if __name__ == "__main__":
    #Load the dataset
    #This is the one that was uploaded earlier to pastebin
    Q = DataFrame.from_csv("file.csv", index_col=False, header=None)

    #Clamp min to zero
    #Q[1] = Q[1]-min(Q[1]);

    #Grab the peaks' model
    F = getEnvelopeModels(Q)

    #Evaluate model to get vector
    u = numpy.array(map(F[0],Q[0]))

    #Workout the scaling coefficient by finding the
    # 'weakest' pulse and scaling everything to that one
    # The weakest pulse was chosen to avoid distorting the limited
    # range of the weak pulses by "stretching" them to higher amplitudes.
    coef = numpy.min(F[1][1:,1]) / numpy.max(F[1][1:,1])

    #Modify the evaluated model to workout a scaling vector
    u = coef/(u/numpy.max(u));

    #Plot everything
    #Original signal
    pylab.plot(Q[0],Q[1]);
    pylab.xlabel('x');
    pylab.xlabel('y');
    pylab.title('Original');
    pylab.grid(True)
    pylab.figure();
    #Scaled signal
    pylab.plot(Q[0],Q[1]*u);
    pylab.xlabel('x');
    pylab.xlabel('y');
    pylab.title('Scaled');
    pylab.grid(True)
    pylab.figure();
    #Scaling coefficient
    pylab.plot(Q[0],u,'r');
    pylab.xlabel('x');
    pylab.xlabel('y');
    pylab.title('Scaling Vector');
    pylab.grid(True)
    pylab.show()

    #Apply scaling to the original DataFrame
    Q[1]=Q[1]*u;
    #Save the DataFrame
    Q.to_csv('processed.csv',sep=' ',index=False,header=None)

The plots are as follows:

Original Signal

Scaling Vector

Processed Signal

By the way, processed.csv can be downloaded from here.

There is a tiny little problem here with the signal having a floating minimum. I tried to clamp it to zero by subtracting the min from the vector (it is commented in the code at the moment) but it requires more work than this.

So, although the local maxima envelope effect was diminished there is still the local minima which might still carry some effect.

I did not want to do anything more drastic to clamp the minima to zero before checking that this would be alright (?). It would basically take the same "scaling" idea but now applied to the "lower part" of the signal.

Hope this helps.


confidence interval estimate for coherence

You might look at

G. Carter, C. Knapp and A. Nuttall, "Statistics of the estimate of the magnitute-coherence function," in IEEE Transactions on Audio and Electroacoustics, vol. 21, no. 4, pp. 388-389, Aug 1973.

doi: 10.1109/TAU.1973.1162487

Abstract: Expressions for the statistics of the estimate of magnitude coherence between two random stationary Gaussian processes are presented. These statistics include the probability density function, the cumulative distribution function, the bias, and the variance of the estimator. The expressions presented are in convenient and accurate forms for digital computer evaluation. Graphical examples of the bias and variance are included. Simple approximations are also given for the maximum bias, variance, and mean-square error.

keywords: {Density functional theory;Distribution functions;Gaussian processes;Mean square error methods;Probability density function;Random processes;Random variables;Statistical distributions;Statistics;Upper bound},

URL


Highest rate achievable by current DSP and A/D chips?

This is a pretty benign requirement. The Open Source software defined radio community does significantly more than 100 MS/s complex samples in software running on stock PC hardware. The samples come from stock SDR hardware.

However, you wouldn't be doing a 100s of taps long channel estimator / LMS filter at that rate.

It's impossible to estimate the computational complexity of your problem without at least knowing what kind of interferer you want to cancel. However, 10 MS/s really doesn't sound dramatic at all.


Invertibility of Room Impulse Response: Reproducing Research Paper

I wrote the following code for the same purpose and thought of sharing it. For simplicity I assumed impulse response to be simple averaging type(Low Pass Filter).

$h[n]= [0.5,0.5]$

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

hlp         =   np.array((0.5,0.5))
#hhp        =   np.array((0.5, -0.5))
f, T        =   signal.freqz(hlp)

fig1 = plt.figure()
plt.title('Frequency response')
ax1 = fig1.add_subplot(111)
plt.plot(f, 20 * np.log10(abs(T)), 'b')
plt.ylabel('Amplitude [dB]', color='b')
plt.xlabel('Frequency [rad/sample]')
ax2 = ax1.twinx()
angles = np.unwrap(np.angle(T))
plt.plot(f, angles, 'g')
plt.ylabel('Angle (radians)', color='g')
plt.grid()
plt.axis('tight')

fig2 = plt.figure()
plt.title('Nyquist Plot')
plt.plot(np.real(T), np.imag(T))
plt.ylabel('Im(H(w))')
plt.xlabel('Re(H(w))')
plt.show()

This code resulted the following graphs.

 I plotted this frequency response just for a sanity check(to make sure the numerator and denominator coefficient is taken correctly).

Now I understand it the paper better :-))).


Your DC blocker transfer function is

$H(z) = (1-a)\frac{1-z^{-1}}{1-(1-a)z^{-1}}$

and an alternative (equivalent) difference equation is

$y_n = (1-a)(x_n - x_{n-1}+y_{n-1})$

Although this filter provides no gain it is possible to find inputs that exceed +1. For instance this input [1 -1 -1 -1 1].

Because you seem to work with a floating point realization you should not need to scale down before filtering.

Edit:
The output of your filter is bounded this way: $|y_n|<2(1-a)$. So an input gain of 0.5 guarantees that the output of the filter can't exceed +1 if the input is bounded between -1 and +1. Can you come up with an input that produces an output value that exceeds +1.5 or an input that produces an output values that is less than -1.5?

Edit2:
Filters can provide amplification as well as attenuation, otherwise they wouldn't be that interesting. Your filter does not provide amplification in the sense that a sinusoid passing through your filter will not have an increased amplitude. You can see that by plotting the amplitude/frequency response of your filter. However, for more complex inputs the non-linear phase response of your filter can cause 'overshoots' in the output. An allpass filter for instance, even though it has a completely flat amplitude response can also provide 'overshoots' for fullscale inputs.

How to deal with the spikes depends very much on what the other modules in your signal path are doing. Maybe some of the modules create headroom maybe some of them consume headroom. I don't know. Considering your DC blocker in isolation then you can apply a headroom gain of say 0.7 or 0.8 and then saturate your output. Although this setting is likely to perform some saturations my guess is that they will be completely inaudible. You will have to confirm this by experiment.


Inverse mapping (affine transformation) for singular transform matrix

Usually, when an inverse is required but the matrix you are taking the inverse of is not guaranteed to be invertible, then the (Moore-Penrose) pseudoinverse is used:

In R, if $T$ is

T
        [,1] [,2] [,3]
[1,]    1    1    0
[2,]    1    1    0
[3,]    0    0    1

Then the pseudoinverse calculated using ginv is

ginv(T)
     [,1] [,2] [,3]
[1,] 0.25 0.25    0
[2,] 0.25 0.25    0
[3,] 0.00 0.00    1

Clearly, this won't give the identity matrix when multiplying $T$, but it's as close as you can get with $T$ singular.

T*ginv(T)
     [,1] [,2] [,3]
[1,] 0.25 0.25    0
[2,] 0.25 0.25    0
[3,] 0.00 0.00    1


Advantages of Splines in Compressed sensing

That's the old question:

You say it's better, and wonder why, and everyone can only ask you to mathematically define better.

As soon as you define a measure for processing quality, you will be able, typically with not too much linear algebra, derive based on the properties of your signal, your measure and your basis, why one is better than the other.


How to difference a loud sound (playing notes), from a soft noise in spectrum. (Pitch detection with bass)

The intuitive answer is to threshold with respect to the magnitude.

You are performing a short-time Fourier transform (STFT), that is, segmenting the original into short time windows and performing Fourier transform on each of them. The output of STFT is a frequency-time (t-f) bin table of complex values.

There is an inherent tradeoff in STFT. The shorter your windows are the higher is your time resolution and the lower is your frequency resolution. If you want a high-frequency resolution you will have to use longer windows at the cost of poor time resolution. This paper gives a clue on window length optimization. otherwise, you can simply use the common music window length from papers or online examples. For bass, the frequencies are low = long wavelength, hence, you will want longer windows.

Basically, you can get the magnitude and\or the angle of your t-f table. Assuming your bass has much more energy than the noise gathered by the system, you can look at the magnitude of bins in the required frequency range and drop time segments below a certain threshold. Play a bit with the result and present it as an image. It will give you a science of where to cut. This is a MATLAB example of presenting the magnitude as an image (also called spectrogram)


What is the difference between continuous, discrete, analog and digital signal?

A signal is indeed a function. Given a signal $f(x)$, according to whether continuous or discrete for both the variable $x$ and the function $f(x)$, there are four types of combinations:

(1) $\mathbf{continuous}$ $x$ and $\mathbf{continuous}$ $f(x)$

This is the most common $\mathbf{analog}$ signal.

(2) $\mathbf{continuous}$ $x$ and $\mathbf{discrete}$ $f(x)$

For this one, we can imagine the ideal base-band waveform used in digital communication. As this one:

(3) $\mathbf{discrete}$ $x$ and $\mathbf{continuous}$ $f(x)$

This is indeed the signal in most the "digital signal processing" textbooks. An example, as others have pointed out, is the output of the CCD sensor.

(4) $\mathbf{discrete}$ $x$ and $\mathbf{discrete}$ $f(x)$

This is the $\mathbf{digital}$ signal. Digital signals are used in practical implementation aspects, and they actually exist in a conceptual manner.

If we concern the discrete feature of the function, the problem will be more complex, therefore, in most "digital signal processing" textbooks, the signals are $\mathbf{not}$ digital indeed. An interesting fact is that, for the the classical textbook by A. V. Oppenheim, the name was "digital signal processing" in the 1st edition, but the name was changed to "discrete-time signal processing" for the later editions.


How to estimate taps required for Parks-McClellan filters

If you're using MATLAB, the function firpmord exists to help you with that. Like some other MATLAB functions, it doesn't link to any libraries or mex files, it's simply MATLAB code that runs. The only reason I mention it is that when you open this function (open firpmord), it has a subfunction remlpord that was written by (ta-daaaaa!) J. H. McClellan himself. It's using a matrix of hardcoded numbers and references Rabiner & Gold, Theory and Appications of DSP, pp. 156-7. The method therefore must be somewhat empirical, though I won't argue one way or the other. In any case, you can study this function (it's very short) and write your own based on it. I failed to find any specific papers that address the problem though, but perhaps there are books.


FIR filter design with time and frequency domain constraints

The most straightforward way to impose constraints in the frequency domain as well as in the time domain for the design of (linear phase) FIR filters is to use linear programming.

Constraints on the step response or on the impulse response are naturally linear, and constraints on the frequency response can also be formulated as linear constraints in the case of linear phase FIR filters.

Take a look at this document for some examples. More details on the linear programming formulation of the FIR filter design problem can be found in this paper.


Selecting the Low Pass FIR Filter, Decimation, DSP

When you are decimating a signal x(n), to x(an) (a>1),in the frequency domain the spectrum expands. This might lead to aliasing of the adjacent spectral bands.So to avoid it we use a "anti aliasing filter".

Now choosing an FIR low pass filter depends on many factors like, order of the filter, attenuation in stop band, range of pass band frequencies, etc..

Till my knowledge, sampling is not at all a point of concern while designing the filter. you only decimate the signal after filtering.

To quote an example, in the base band part of digital communication, at the transmitter side, you upsample the signal and you use a anti imaging filter, to avoid multiple copies of the same spectrum being repeated. where as in the receiver you first filter out the signal and then decimating it for further processing.

Cheers,
phani tej


Finding the Correlation Matrix of Wiener Filter using Matlab

Here is a simple example of Wiener filter in MATLAB. First, the input data, denoted as $\mathbf{x}(n)$. Note that we are considering only at time instant $n$ and an extension for $N$ time samples can be made to incorporate more time domain samples by concatenating these vectors together to get a matrix $\mathbf{X}=[ \mathbf{x}(n) \ \ \mathbf{x}(n+1) \ \ \cdots \ \mathbf{x}(n+N-1)]$ (that is a vector for each time instant). But for the single time instant case, lets say that $\mathbf{x}(n)$ is IID Gaussian with some input mean and varaince, then this input vector can be generated in MATLAB as:
x = sqrt(inputVariance)*randn(10, 1) + inputMean;

The 10 comes from your problem description. Now that you have some example input data, you can go ahead and compute the correlation matrix as you defined $\mathbf{R} = \mathbb{E}\big[ \mathbf{x}(n) \mathbf{x}(n)^H \big]$:
R = x*x'/size(x, 2); % size(x, 2) = # of time domain samples

Since our example input data is easy to understand (Gaussian), we can know what to expect to find in the correlation matrix. We generated the input data to be uncorrelated, so the off diagonal elements should be small and the diagonal elements should be close to the variance. One useful way to look at these matrices in MATLAB is as an image using:
figure
imagesc(R)

After doing this, you might be confused though. Because you expect to see something like this:

But instead you might see something closer to the next figure. This is because you are only taking a single time domain sample! Not enough samples for the statistical properties to emerge from the correlation matrix. In fact, it is an estimate of the true correlation matrix and you can't expect an estimator to be that good with only one sample (the first image was generated by following the steps but using x = randn(10, 100) to generate 100 time domain samples).


Detecting Arabic verses spanning multiple lines

This can be solved fairly straightforwardly with simple template matching. I don't know exactly how you have it set up, so I'll just describe the algorithm generally and use illustrations.

Observe that the verse numbers have a distinctive border that can easily be used to detect the start and end of a verse. So create a binarized template for that pattern and store it. Something like this:

Since the number of lines in a screen are known in advance (you're formatting the page) and each verse has a constant height, you can easily infer (algorithmically) where the Y coordinates for the centerlines of the verses should be on the screen. This demonstrates the idea:

When the user touches a verse, get the X-Y coordinates and snap the Y coordinate to the nearest verse center.
Then starting with the X coordinate, perform a simple template matching (cross-correlation) across that row. The first match (peak in the cross-correlation) in the forward direction (to the left), will be the end point for the verse. If there are no matches in the reverse direction (to the right), then move up one verse (which you can do, because you know the Y coordinate of the centerline) and repeat. The first match from the left end will be the start point of the verse. Similarly, if there is no forward match on the line, move down one line and repeat.

Here's a short illustration of the idea. The yellow box is where the user touches the verse. You then do the cross-correlation with your template and the blue circles will be the match.

I also use template matching in this answer, if you're interested in seeing it in action.
Once you've determined the start point for the verse, then use an Arabic text recognizer to infer the verse number inside that border and play the corresponding audio file.

Simpler solution:

A simpler solution, if you don't want to go through this is to store the X-Y coordinates of the verse starting points (keep it simple and use the center points) and once you get the coordinates of the user input, you can again snap it to the centerline and then walk backwards to see where the verse starts. This might have the advantage of being faster.

I didn't put this forward as the first solution because you seemed to reject a similar idea in the comments. In the end, it depends on your constraints — would you rather do computational work (template matching — which, by the way, also requires you to store the template) or using memory (storing coordinates).

If I were you, I'd probably go with this one, but the image processing solution can be fun to try.


First, you have to extract the Room Impulse Response (RIR) of the space from your sine sweep recordings. You can do this by taking the ratio of the FFTs of the  measured response and sine sweep input and then taking its IFFT, rir = ifft(fft(measured, Nfft)./fft(sin_sweep,Nfft)). I generally use Nfft = 2^nextpow2(length(measured))
You can truncate the RIR to $L$ samples, depending on how long the reverb tail is.
Break your music source into frames of $M$ samples.
For each frame, you want to convolve a signal of length $M$ with the RIR of length $L$. The output of the resultant convolution will be $M+L-1$ samples. So you have to pick an FFT size accordingly. Then multiply the FFT of the current frame of the signal, with the FFT of the RIR.
The final step is overlap add.

Here is some MATLAB style code.
% sig - Input signal
% rir - Room Impulse Response extracted from sine sweep
% M - length of each frame
% L - length of the RIR

nframes = round(length(sig)/M));
Nfft = 2^nextpow2(L+M-1);
conv_sig = zeros(nframes*M,1);
for frame_no = 1:nframes
    sig_frame = sig((frame_no-1)*M+1:frame_no*M);
    conv_frame = ifft(fft(sig_frame,Nfft).*fft(rir,Nfft));
    conv_frame = conv_frame(1:M+L-1);
    conv_sig((frame_no-1)*M+1:frame_no*M+L-1) = conv_sig((frame_no-1)*M+1:frame_no*M+L-1)+conv_frame
end


how to compute the signal passing from the low pass filter?

Low pass filters get the name because they pass signals at low frequencies and attenuate signals at higher frequencies. The definition of low and high depend on the the cut-off frequency of the filter. In your case, the cut-off frequency is $\frac{\pi}{5}$. Your filter is the ideal low pass filter where signals above the cut-off frequency are completely zeroed out. The input signal is a perfect cosine signal so all of its energy is concentrated at one specific frequency. It happens to be at a frequency above the cut-off frequency, so what should happen?


1) No, if you have a particularly low level in one of the FFT points it will indicate that you have low energy in that particular frequency, but you can have several other frequencies at the same time.

2) In theory yes, but in practice the numerical precision will make the resulting signal slightly different than the original one.

Note: before applying the FFT you'll need to decode the mp3 file in order to get the audio data.


What do we achieve by using PAM symbols symmetrical with respect to 0?

The probability of error for PAM depends on the distance between symbols. This means that the following two PAM constellations, both with distance $2$, will have the same probability of error:
\begin{align}
\mathcal{A}_1 &= \lbrace 2, 4 \rbrace \\
\mathcal{A}_2 &= \lbrace -1, 1 \rbrace
\end{align}

When using $\mathcal{A}_1$, the average energy spent per transmitted symbol is $(4+16)/2=10$ joules per symbol. When using $\mathcal{A}_2$, however, the average energy is $(1+1)/2=1$ joule per symbol: ten times less.

Even if transmitting a signal with constant DC (as pointed out by @AlexTP) were not a problem, you still want to use constellations that are symmetrical around zero because they minimize the energy expenditure for a given BER.

A good and simple excercise is to prove that this is true in general, and not just for the example I provided.


Why I get another spectrum when upconvert an o-qpsk baseband signal with a bladerf compared to a hackrf device?

One source of this is spectral regrowth (overdriving a device into saturation such that the pulse shaping that limits the spectrum is degraded), but compared to typical spectral regrowth spectrums I have seen, this looks more like your signal was greatly attenuated and then amplified together with noise, and then that combined signal passed through a bandpass filter. Check if there is a 5 MHz filter earlier in the Tx chain to further support this theory, and then also try to reduce the power level at the input to the device.


Pitch detection to find music notes (from voice)

1) Get some practical experience on what you are studying.  By this I mean is take your recordings of your voice and look at them close up on you computer.  You can use a program like Audacity.

2) Yes, you will have to learn about the Discrete Fourier Transform (DFT, the FFT is a computationally efficient version).  I write blog about it explaining it in my own way.  You can find the link on my profile page.

3) Pick a different topic that is worthy of a Master's Thesis.  There are many techniques for pitch determination, which is a pyscho-acoustic phenomenon vs frequency estimation which is basically a physics problem.  Since you are a singer, you might be interested in how sound waves are attenuated as they travel through air.  You simply have to set up some microphones at different distances from the source and look at their waveforms to see what I am talking about.

4) Depending on your math chops, you may consider coming up with data visualization techniques for signal data.  There is lots of room for creativity here.

Hope this helps.

Ced


Keyword spotting - constant phoneme length

There is no relation between fixed frame length and varying phoneme length. Since frames are just 0.01 seconds and average phoneme length is about 0.1 seconds, there would be a frame located around a center of the phone. And in that point of time DNN detector that spans 20 frames around central frame will detect the phone of interest even if it does not have any time-varying architecture. The fully connected DNN they use could do that pretty reliably.

Of course they could use something more advance and adapted to varying length sequences, for example, multilayer bidirectional LSTM, that would result in better detection accuracy but at the same time it would be much more computation-consuming and thus battery-consuming, it would also affect reaction time, thats why they prefer simpler architecture.

Their approach would not work for some significant phoneme length variations. For example you can test that "heeeeeeey siiiiiiiiri" does not work that reliably. At the same time one could argue that "heeeeeeey siiiiiiiiri" is not the same as "hey siri" for a human either and should not activate Siri. For a balance between energy consumption and detection accuracy such simple solutions are a good choice.


Noise removal in medical segmented image

One solution I have found is as follows:

Thresholded on grayscale value.
Remove objects on the basis of size.
Some more morphological operations.