Real signals are "mirrored" in the real and negative halves of the Fourier transform because of the nature of the Fourier transform.  The Fourier transform is defined as the following-

$H(f) = \int h(t)e^{-j2\pi ft}dt$

Basically it correlates the signal with a bunch of complex sinusoids, each with its own frequency.  So what do those complex sinusoids look like?  The picture below illustrates one complex sinusoid.

The "corkscrew" is the rotating complex sinusoid in time, while the two sinusoids that follow it are the extracted real and imaginary components of the complex sinusoid.  The astute reader will note that the real and imaginary components are the exact same, only they are out of phase with each other by 90 degrees ($\frac{\pi}{2}$).  Because they are 90 degrees out of phase they are orthogonal and can "catch" any component of the signal at that frequency.

The relationship between the exponential and the cosine/sine is given by Euler's formula-

$e^{jx} = cos(x) + j*sin(x)$

This allows us to modify the Fourier transform as follows-
$$
H(f) = \int h(t)e^{-j2\pi ft}dt \\
= \int h(t)(cos(2\pi ft) - j*sin(2\pi ft))dt
$$

At the negative frequencies the Fourier transform becomes the following-
$$
H(-f) = \int h(t)(cos(2\pi (-f)t) - j*sin(2\pi (-f)t))dt \\
= \int h(t)(cos(2\pi ft) + j*sin(2\pi ft))dt
$$

Comparing the negative frequency version with the positive frequency version shows that the cosine is the same while the sine is inverted.  They are still 90 degrees out of phase with each other, though, allowing them to catch any signal component at that (negative) frequency.

Because both the positive and negative frequency sinusoids are 90 degrees out of phase and have the same magnitude, they will both respond to real signals in the same way.  Or rather, the magnitude of their response will be the same, but the correlation phase will be different.

EDIT: Specifically, the negative frequency correlation is the conjugate of the positive frequency correlation (due to the inverted imaginary sine component) for real signals.  In mathematical terms, this is, as Dilip pointed out, the following-

$H(-f) = [H(f)]^*$

Another way to think about it:

Imaginary components are just that..Imaginary! They are a tool, which allows the employ of an extra plane to view things on and makes much of digital (and analog) signal processing possible, if not much easier than using differential equations!

But we can't break the logical laws of nature, we can't do anything 'real' with the imaginary content$^\dagger$ and so it must effectively cancel itself out before returning to reality. How does this look in the Fourier Transform of a time based signal(complex frequency domain)? If we add/sum the positive and negative frequency components of the signal the imaginary parts cancel, this is what we mean by saying the positive and negative elements are conjugate to each-other. Notice that when an FT is taken of a time-signal there exists these conjugate signals, with the 'real' part of each sharing the magnitude, half in the positive domain, half in the negative, so in effect adding the conjugates together removes the imaginary content and provides the real content only.

$^\dagger$ Meaning we can't create a voltage that is $5i$ volts. Obviously, we can use imaginary numbers to represent real-world signals that are two-vector-valued,  such as circularly polarized EM waves.


Why do we need DFT when we already have DTFT?

The answer is the same to the question: "Why do we need computers to process data when we have paper and pencil?"
DTFT as well as the continuous-time Fourier Transform is a theoretical tool for infinitely long hypothetical signals.
the DFT is to observe the spectrum of actual data that is finite in size.


Necessary Conditions for stability in z domain?

The one and only condition for BIBO stability of a 1D discrete-time system, in the z-domain, is that its transfer functions's ROC (region of convergence) should include the unit circle : $|z| =1$. Therefore, it's a necessary and sufficient condition for BIBO stability of a 1D SISO system.

There are no other conditions (to my knowledge).

EDIT [based on comments] :
By definition of ROC (region of convergence) there cannot be any poles inside the ROC region. This's related not with stability but with convergence. And also, causality is not related with stability too.

Properties of ROC is a different question to be answered.


FFT seems more accurate than MUSIC algorithm?

I created a signal with two sinusoids and added increasing amounts of random noise to it. The results of running the FFT and MUSIC on this signal are shown in the image below.

You are using pmusic(). MUSIC (MUltiple SIgnal Classification) is a general algorithm that can estimate all kinds of parameters out of a signal. pmusic() is a MATLAB function that utilizes this algorithm to evaluate the locations of peaks' location in the spectrum of the signal. More to the point, it does not mind a better representation of the level in other frequencies (not the location of the peaks).

I thought the whole point of the MUSIC algorithm is that it is supposed to give more accurate results than the FFT.

Not at all! FFT is an analytical method to represent a given signal in the frequency domain. There can be a vast discussion on the accuracy of the method given the parameters it is executed with but it is a mathematical operation with a close solution. MUSIC algorithm, however, is an optimization method. As an optimization method, it is iteratively estimating the optimal solution to a problem. It can actually generate two different solutions given the same problem and two different initialization. Optimization methods are usually used when a close solution cannot be formulated. Each of these two can be better depending on the specifics of the problem. For example, how close are the frequencies? How fast to you want your solution? What is your tolerance for accuracy?


Dynamic Time Warping of inverted binary signals

I had a similar experience, i had to align side channel traces, I've tried to implement my own DTW that do lot of things not only aligning, but the this method was extremely slower than separating the processing and DTW.

So if you are using python the FastDTW will do the job insanely quicker than any handmade DTW, so the laborious way is the quickest way to get the aligned signals, besides that on very large signals DTW need huge amount of memory, so separating the two process is in my experience is the best solution.


How is multiplexing achieved in spread-spectrum modulations like CSS?

disclaimer

Because of the lack of document about LoRa CSS in Internet, the analysis below may be wrong with respect to LoRa system, not to the (linear) CSS principle. Any comment or update are appriciated. Thanks.

LoRa CSS modulation

LoRa CSS uses linear chirp spread spectrum. To transmit $\log_2 M$ bits, LoRa CSS modulation divides bandwidth $[-\frac{B}{2} \dots \frac{B}{2}]$ to $M$ portions.

In the figure, to transmit symbol $-k$ with $-M/2 < k < M/2$, LoRa use the waveform having the instantaneous frequency:
$$
  f(t)=\begin{cases}
               \mu\frac{B}{T}t - k\frac{B}{M} + B \qquad -\frac{T}{2} \leq t \leq -\frac{T}{2} + k\frac{T}{M}\\
               \mu\frac{B}{T}t - k\frac{B}{M} \quad\qquad -\frac{T}{2} + k\frac{T}{M} < t \leq \frac{T}{2}\\
               0 \qquad\qquad \textrm{otherwise}
            \end{cases}
$$

with $\mu = +1$ or $\mu = -1$. The spreading factor is $M = BT$.

By sampling at exactly $B$ Hz, $f(t)$ becomes continuous, we can works directly with $f(t) = \mu\frac{B}{T}t - k\frac{B}{M}$ with $-\frac{T}{2} \leq t \leq \frac{T}{2}$ in linear CSS decoding framework. (In fact it is not simple like that because of the transition phase at the discontinuos point in the figure. However, it does not change the insight of the answer. I will update my answer "later" :). I appreciate the fact that someone interested in math can do it properly)

\begin{eqnarray}
f(t) = \mu\frac{B}{T}t - k\frac{B}{M} = \frac{1}{2\pi}\frac{d\phi(t)}{dt}\\
\phi(t)  = 2\pi\left(\mu \frac{B}{2T}t - k\frac{B}{M}\right)t
\end{eqnarray}

To decode chirp, we use inverse chirp at receiver:
\begin{equation}
c_d(t) = e^{j2\pi\left(-\mu \frac{B}{2T}t^2\right)}
\end{equation}

Decoding by multiplication in $-T/2 \leq t \leq T/2$:
\begin{equation}
z(t) = c(t)c_d(t) = e^{j\phi(t)} \times ^{j2\pi\left(-\mu \frac{B}{2T}t^2\right)} = e^{j2\pi\left(\mu \frac{B}{2T}t - k\frac{B}{M}\right)t} \times e^{j2\pi\left(-\mu \frac{B}{2T}t^2\right)} = e^{j2\pi(-kB/M)t}
\end{equation}

Remember that we have assumed working on the signal sampled at $B$ Hz ? Thus instead of having $z(t)$ what we got is a sequence of sample $z[m] = z(t=m/B)$ with $-M/2 \leq m \leq M/2$. Note that $M=BT$ if you feel lost.

Now DFT of M samples $z[m]$ to have $Z[u]$ with $-M/2 \leq u \leq M/2$, you will have exactly only one peak at the bin $u = -k$ (if no synchronization error, no Doppler, no multipath). Bingo.

Multiplexing

It depends how we define "multiplexing". If we use the term in multiple access level (MAC layer), we can have TDMA and FDMA-like. If we are talking about user communicating simultaneously, the only thing that (I think) we can use is the factor $B/T$, or the spreading factor $M = BT$ if $B$ is fixed.

Let's fix the spreading bandwidth $B$. The chirp duration $T = M/B$ changes if $M$ changes and so do the chirp wavefor that depends on $B/T$ : $c_1(t) = e^{j2\pi\left(\mu \frac{B}{2T_1}t - k_1\frac{B}{M_1}\right)t}$

Thus
\begin{equation}
z(t) = c(t)c_d(t) = e^{j2\pi\left(\mu \frac{B}{2T_1}t - k_1\frac{B}{M_1}\right)t} \times e^{j2\pi\left(-\mu \frac{B}{2T}t^2\right)} = e^{j2\pi(\mu \frac{B}{2}(\frac{1}{T_1} - \frac{1}{T})t - k_1B/M_1)t} = e^{j2\pi(\mu \frac{B^2}{2}(\frac{1}{M_1} - \frac{1}{M})t - k_1B/M_1)t}
\end{equation}

$M_1 \neq M$, $z(t)$ has instateneous frequency varying in time. The DFT result should be something like that, i.e. containing no relevant peak.

LoRa defines 6 spreading factor $2^7$ to $2^{12}$ (7 if we take the special SF $2^6$). Thus, for a fixed bandwidth, LoRa system can have simulatenously 6 (or 7) users.


How do you compute for fc and t in modulating an I and Q component

The answer is: it depends on the pulse shape.

A 16-QAM baseband signal can be written as $$s_{BB}=\sum_{k=0}^{N-1} a_k p(t-kT_p),$$ where $a_k$ is a complex symbol taken from your 16-QAM constellation, $R_p=1/T_p$ is the pulse rate (how many pulses per second you transmit), and $p(t)$ is the pulse shape (ideally, $p(t)$ is a Nyquist pulse).

The bandwidth of $s_{BB}(t)$ is the same as the bandwidth of $p(t)$ (assuming you're transmitting uncorrelated data). Assume this bandwidth is $B$. Then, after mixing with the carrier, the passband signal's maximum frequency is $f_c+B$, and you need a sampling frequency $f_s>2(f_c+B)$ to satisfy Nyquist.

If you don't know what the pulse shape is, then you're probably transmitting rectangular, non-overlapping pulses. As a rule of thumb, I define the bandwidth of a rectangular pulse as five times the pulse rate: $B=5R_p$. Note that many textbook authors, such as Stallings, use $B=R_p$. I prefer to be conservative and define a larger bandwidth.

Then, your sampling frequency is $f_s>2(f_c+5R_p)$.

As an example, let's say that you want to transmit 1000 pulses per second, with a carrier frequency of 1 MHz, and rectangular pulses. Then,

fs = 2*(1e6+5*1000);  % or larger
t = 0:1/fs:Tend       % Tend depends on how many pulses you transmit


That actually would be a great interview questions. Let's derive it from scratch:

$\text{sinc}(x)$: The Fourier Transform of a $\text{sinc}(x)$ is a rectangle. The normalized $\text{sinc}()$ function results in a bandwidth of 1 Hz (from -0.5Hz to 0.5Hz)
$t-10$: The $(t-10)$ term is irrelevant to the bandwidth: time shift
is equivalent to multiplying with $e^{-j\omega t_0}$, where $t_0$ is time shift, in
the frequency domain. That changes the phase but not the amplitude
or bandwidth
$125*t$: that is a time compression which results in frequency
stretching. Since we compress time by 125 the bandwidth increases by
125, i.e. to 125 Hz (from -62.5 Hz to 62.5 Hz)
$(\cdot)^2$: Squaring is multiplication in the time domain which corresponds to convolution in the frequency domain. The convolution of two equal rectangles is a triangle. The bandwidth doubles to 250 Hz (from -125 Hz to 125Hz). In general the bandwidth of the multiplication of two arbitrary signals is the sum of the individual bandwidths. That's a simple consequence of the multiplication laws for sine and cosine: you get only sum and difference frequencies.
$50*(\cdot)$: that's a simple scale factor that doesn't change the shape of phase or amplitude.

So overall we would expect the spectrum to be of triangular shape (on a double linear plot) going from -125 Hz to +125 Hz with a total bandwidth of 250 Hz.

Writing the Matlab code for this is not as simple as it looks. The code itself is trivial, however you have to pick a sample rate and a proper time axis. In order to do this you need to assess the general shape of the signals first. For example the time domain signal technically goes from $-\infty$ to $+\infty$ and is centered around $t$ = 10s. t = 0:.01:5, really won't do the trick. So in essence you have to do the theoretical analysis before you can write the code (which tends to be a recurring theme in DSP :-))

% 1000 HZ sample rate, since we now the BW is smaller than 500 Hz
fs = 1000;
% pick a good time length. Because of the time shift we need significantly more
% than 10s. Technically the sinc() has inifite length so we need to extend
% it far enough to cover the vast part of the energy AND we need to start
% at negative times, not at t = 0;
n = 2^16; % that's about 60s. Should be plenty
t = (-n/2:n/2-1)'/fs;   % time axis in seconds, symmetric around t = 0;
x = 50*sinc(125*(t-10)).^2;
% frequency axis:
f0 = fs/n;  % frequency resoution in Hz
faxis = (-n/2:n/2-1)'*f0;   % frequency axis in Hz seconds, symmetric around f = 0;
fx = fft(x);
plot(faxis,circshift(abs(fx),n/2));
title('Spectrum of 50*sinc(125*(t-10)).^2');
xlabel('Frequency in Hz');


How do KLT and SVD differ from each other?

From the direct implementations of the two, SVD and KLT are different since $S$ in SVD is a diagonal matrix, while $C$ in KLT is not. Besides, KLT is originally used to describe a stationary process. In its expansion, an infinite number of basis functions is required to
form a complete set. Yet for a random vector $X$, it can be rephreased in terms of a linear combination of orthonornal basis vectors, then covariance matrix $C$ is furtherly derived.

For a 2D image, we can view it as $M$ realizations of stochastic process $x$ where $M$ is the row of the image. And the assumption is that all rows have the same row covariance matrix $C$, then the KLT can be somehow converted to the form of SVD. In other words, when considering a multivariate stochastic process, the SVD and the KLT are computationally equivalent. That's why in most cases people view the both as the same. For mathematical details, refer to The KL transform and eigenimages.

In the application of image compression, the global KLT is used for the whole image, while SVD is better applied to each block. In Patrick et al's Hybrid KLT-SVD Image Compression, they combined the two method to implement the adaptive transform on the blocks to achieve a better coding efficiency.


How to calculate sine wave SNR from oscilloscope

For all the reasons above, an oscilloscope is inadequate to this task.  Here's a couple of alternatives:
Method 1: find someone with a really good spectrum analyzer, and borrow time on it.
Method 2: build a really good notch filter (for a one-off you'll want to use analog components, and good ones).  Characterize it (I didn't say this would be easy).  Then run your 1kHz signal through the notch and look at the remainder in the O-scope.  Once you get rid of the 1kHz tone, whatever remains is the stuff you don't want -- characterize that, add the 1kHz tone back in on paper, and you've characterized the net effect of your DAC and your notch filter.
Note that method 2 requires some significant analog chops.  If your go-to solution involves white proto-board, junkbox capacitors and LM741 amps, you probably need to up your game.


references for techniques to reduce recursive filter startup transients

For a given fequency response magnitude, a minimum phase filter implementation will have the shortest possible transient duration.

A minimum phase filter will have most of its energy packed at the beginning. Hence it's also called as a maximum energy filter. So it will produce steady state outputs as soon as possible (among the set of filters having the same frequency response magnitude)

However, a minimum phase filter cannot be linear phase simultaneously. Because linear phase filter would have a symmetric impulse response but that would also imply transfer function zeros outside of unit circle. Consequently it cannot be a minimum phase filter at the same time as for a minimum phase filter all the poles and zeros of transfer function are inside the unit circle.

Alternatively, if you have some reliable estimation of the initial samples, before $n=0$,  of the input signal, then you can initialize the filter according to those values instead of a full blind zero initialization.


Compensating effects of a system with a known transfer function

That is pretty broad question, as it depends on the system.

However, your problem is extremely common: Every communication system works because a receiver is able to estimate the original events (symbols sent by the transmitter) after the signal has gone through a system with transfer function (the channel).

So, receivers need to do exactly that: reverse the effect of the channel as good as possible. The thing that does that is usually called an equalizer; the equalizer you pick, and how much info about the system you need to make that equalizer work, depends on your "main system".

So, if you need to solve that: go ahead; define your main system as well as possible. Find stochastic descriptions for the parts of the system that you don't know with certainty. Then, find a channel model that fits that model, and look for equalizers used with such a channel.


When upsampling the number video frames to allow playing a video in slow motion, rather than for a frame rate increase for smoothness, you will likely need to modify the audio using a time-pitch stretching/shifting algorithm to stretch the audio out (increase the number of samples) for a longer play duration (to match the increase in slow motion video playback time), and to do so without distorting vocal or musical pitches.

If you just interpolate or resample audio without time-pitch correction, you will mess up all the pitch frequencies in the audio.

See this Wikipedia article for a description of several possible time-pitch stretching/scaling algorithms. There are multiple commercial and open source applications and libraries that can accomplish time-pitch modification/correction with varying degrees of quality and artifacts.


Determine if an impulse is minimum phase and compute a zero-phase equivalent

To find the zero-phase impulse response, simply take the Fourier transform of your signal, find magnitude and find the inverse Fourier transform of that. Mathematically speaking,

$$\mathscr{F}^{-1}\{\text{abs}\left(\mathscr{F}\{x[n]\}\right)\}$$

This essentially says that the entire Fourier transform of your new sequence is real-valued. All real-valued numbers have zero phase.


I figured out the answer. It is an "ionosondes".


Why is incoherence important for compressive sensing?

It is necessary when reconstruction is considered. Simply imagine the case when  $A = \Phi \Psi$  has a high coherence, e.g. all columns are exactly the same and indistinguishable, then there is no way to tell which columns of $A$
correspond to  non-zero elements of $a$  and together were contributed to vector $Y$. If you look at equation 2 below, you can say vector $Y$ is summution of columns of $\Phi \Psi$ multiplied in corresponding elements from vector $a$. When the signal is sparse (first assumption), but we do not know our measurement vector resemble which columns of $\Phi \Psi$ (hence using it to decode), our sparse reconstruction fails.
$$ Y = \Phi X ~~~~~~~~~~~~~~~~~(1)$$
$$ Y = \Phi \Psi a ~~~~~~~~~~~~~~~(2)$$
$Y$ being measurement vector and $a$ being sparse representation of signal of $X$, so $a = \Psi^{-1}X$. Basically, the problem of CS reconstruction is finding our which non-zero elements from the sparse matrix contributed to the measurement and how much was their contribution. If you do not have a unique enough corresponding column for $A$ then you cannot tell which elements were the contributor and how much was the contribution.


Impulse response from difference equation without partial fractions

First of all you can apply partial fraction expansion method to get the inverse DTFT to find the impulse resposne h[n] for this system. You have to assume casual ot anti-causal system.

$$H(e^{jw}) = \frac{3-e^{-jw}}{1 - 0.3e^{-jw} + 0.1e^{-j2w}} = \frac{3-e^{-jw}}{ (1 - \alpha e^{-jw})(1 - \alpha^* e^{-jw})}$$

$$H(e^{jw})= \frac{A}{1 - \alpha e^{-jw}} + \frac{B}{1 - \alpha^* e^{-jw}}$$

Where $\alpha = 0.1500 + j0.2784$ for this particular case. Hence assuming causality and using inverse DTFT look-up tables, deduce as:

$$\boxed{ h[n] = A \alpha^n u[n] + B (\alpha^*)^n u[n] }$$

Time domain recursion is also possible, requiring only two auxiliary values.
For example assuming causality one can use recursion from the LCCDE to obtain $y[n]$ for all $n \geq 0$ from the given LCCDE
$$y[n] - 0.3y[n-1] + 0.1 y[n-2] = 3x[n] - x[n -1]$$

by first rearranging it for $y[n]$ as
$$y[n] = 0.3y[n-1] - 0.1 y[n-2] + 3x[n] - x[n -1]$$

Then the recursion follows for $n \geq 0$ assuming that $x[n]=\delta[n]$
$$
\begin{align}
y[0] &= 0.3y[-1] - 0.1 y[-2] + 3x[0] - x[-1] \\
y[1] &= 0.3y[0] - 0.1 y[-1] + 3x[1] - x[0]  \\
y[2] &= 0.3y[1] - 0.1 y[0] + 3x[2] - x[1] \\
y[n] &= 0.3y[n-1] - 0.1 y[n-2] \\
\end{align}
$$

Note that in computing $y[0]$ you need $y[-1]$ and $y[-2]$ as auxiliary (initial) conditions. Also note that $x[-1]=0$ and $x[0]=1$ and $x[n]=0$ for all remaining $n>0$.

A similar analysis fallows for the anti-causal case.


Emulating guitar pickup position with delays

i won't divulge some proprietary information i have about how to do this, but i'll tell you this hint:

the vibrations on the guitar strings are standing waves.
standing waves as well as any waves on a string are made of of sinusoidal components with frequencies at harmonic intervals, that is each harmonic sinusoid has a frequency that is an integer multiple of the fundamental which is the frequency we will normally associate with the pitch of the note.
with standing waves, there are nodal points for each harmonic at distance that is a submultiple of the length of the string between the bridge on the right and the fret that is terminating the string on the left.
even-numbered harmonics will have a nodal point at the middle of the string and odd-numbered harmonics will not have a nodal point there.
there is a relationship between ration of the even-harmonic energy to the odd-harmonic energy and the position of the pickup and with the position of the pick.

have phun figgering it out.


Computing the derivative from a set of real-time data?

Depending on the level of noise you expect to encounter on your signal, you might actually be able to use the finite differencing you suggest, originally. If you want to know the general trend of your data upward or downward, sum together the, in your case, 24, differences. So, lets say that you have some metric, $D$, for your signal. Lets say that for a given signal, $x[n]$, the first order difference is denoted as $x\prime[n]$ (and is, of course, one entry shorter than $x[n]$). Define $D$ as

$$ D = \sum_{n=0}^{N-2} x\prime[n] .$$

Then, classify your signal as uptrending or downtrending via

$$ D > 0 ~\therefore~\uparrow$$
$$ D < 0 ~\therefore~\downarrow$$

I made some basic assumptions about the nature of the signals in question. I just did the most basic analysis and thought, well, maybe just use a ramp, $r[n]$, corrupted by some noise, $z[n] \sim \mathcal{N}(0,\sigma^2)$. This may or may not be the situation you're looking at, so you'll need to adjust accordingly and see what works for you. But my final model was just the simple
$$x[n] = r[n] + z[n]$$
for the length you gave, $N = 25$.

I went ahead and ran some experiments in MATLAB. I used the following code:

experiments = 1000;
N = 25;

slope_min = -2;
slope_max = 2;
slope_range = slope_max - slope_min;
noise_level = 2;

for i=1:experiments

    slope(i) = slope_range*rand(1)+slope_min;
    r = slope(i)*(1:N);
    z = noise_level.*randn(1,N);

    x = r + z;

    diffsum(i) = sum(diff(x));
end

First, let me show what some of the signals look like...

So, I varied the variance of the noise signal, $z[n]$ and recorded the following correlations between the true slope (up or down trend) of $r[n]$ in comparison with the metric we defined above, $D$. If $D$ does a good job determining the trend of the signal, we should see a strong correlation between downtrends (negative slopes) and negative $D$ and uptrends (positive slopes) and positive $D$.

We can see from the above that $D$ does seem to do a decent job. As the noise level increases (or the signal moves further away from our model), the correlations are not as strong. Try it out on your signals, see if does anything for you.


How to remove constant accelerometer velocity offsetting?

I guess you are sure the device is definitely at rest (v=0) after each event, and not just moving slowly?
If so a couple of possible causes to consider:
NUMERICAL: You might get this type of thing if the integration is inaccurate.
Consider reducing sample interval, and use higher order integration.  However, I would
have expected numerical drift to be sometimes positive, sometimes negative. If it is always negative, numerical seems less likely to me.
INSTRUMENTAL: Does that accelerometer need to be "level"? Is it possible to repeat
the experiment with the accelerometer mounted upside down to see
if the effect reverses?
An immediate suggestion to better understand the problem (which you have probably
already done) is compare the acceleration and velocity plots. As I recall once
the system is underway the two curves should be 90 degrees out of phase, so there
should be time alignment between peaks and zero crossing as in this sketch.

You could check the time alignment at those points. It might indicate if the drift
is arising at a particular part of the cycle (as you suggest).
Pragmatic Adjustment:
A perfect correction would require identification of the exact cause. However, depending on what you are using the velocity curve for it might be acceptable to implement a "pragmatic" fix in processing. As one idea, in an ideal world the areas under the two positive peaks (in accel plot) should balance the area above the negative peak.  The observed drift (dv in picture) should equate to this difference. Hence for example, you may be able to introduce a percentage boost
to the positive acceleration values. With a little thought this should be able
to be automated.
However first I'd suggest looking at the time alignment and see if that reveals anything useful.


Building a Signal Model for Samples from a Sensor

A simple exponential average has the following difference equation:

y(n) = (1-alpha) * x(n) + alpha * y(n-1)

To get a different rise and fall rate, you can dynamically choose alpha depending on whether x(n) is larger than y(n-1). For example, you could build have alpha = 0 when x(n) > y(n-1), and alpha = 0.99 otherwise. This could accomplish the "peak hold" goal, and the second value of alpha would let you control the decay after each peak.


The 3-digit number describes the subsampling of the
chroma (U and V) channels.  A detailed explanation is at

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

In particular, YUV420 means that the U and V channels have half the resolution of the luma (Y) channel, in the horizontal and vertical directions. The sampling method usually present in low to medium quality JPEG images is YCrCb420 (high quality JPEGs are usually YCrCb444 (no subsampling).


You are facing the practical problem of representing a noncausal signal in a programs such as Matlab, Octave, or Python. These programs, in their default, assume that all represented sequences begin at time-index $n=0$.
Therefore, a Hermitian symmetric sequence ( $x[n] = x^*[-n]$ for all $n$ ) ,  which begins at $n=-4$, and ends at $n=4$, will lose this symmetry when interpreted in those programs. The effect is a shift of the signal to the right by $4$ for your particular example.
Hence, what you have compouted is not the FFT of the original symmetric sequence $x[n] = e^{j n}$, but a shifted, non-symmetric, sequence $w[n] = x[n-4]$. Using the DTFT / DFT time-shift property, we see what is actually computed as :
$$ W(e^{j\omega}) = e^{-j4\omega} X(e^{j\omega}) \tag{1} $$ and
$$ W[k] = e^{-j4 \frac{2\pi}{N} k} X[k] ~~~~~ , ~~\text{ for } ~~~k=0,1,...,N-1 \tag{2}$$
Where $N$ is the length of the sequences $x[n]$ and its DFT $X[k]$.
Therefore, eventhough the DTFT / DFT of the desired sequence $x[n]$ has zero imaginary part, as an Hermitian symmetric signal would have, the practically computed $W[k]$ has nonzero imaginary part.
If you wish, you can obtain the non-zero phase samples of the desired DFT $X[k]$ of $x[n]$, from the computed DFT $W[k]$ of $w[n]$, with the following:
$$ X[k] = e^{j 4 \frac{2\pi}{N} k} W[k] ~~~~ ,~~\text{ for } k=0,1,...,N-1 \tag{3}$$
The following MATLAB code demonstrates the issue:
% Demonstrates the Hermitian symmetry and its handling in Matlab/Octave/Python
% environments...

a = -4;
b = +4;
n = [a:b];
x = exp(1j*n);      % Samples of Hermitian symmetrix input: x[k] = x[-k]*

figure,stem( imag( exp(-1j*2*pi*a*[0:(b-a)]/(b-a+1)).*fft(x) ));

The output is:


Random selection of sub-block in an image with matlab


What is the BER for multilevel modulated signal

This is not BPSK! BPSK is binary phase shift keying with binary meaning 2 levels. What the OP is describing is Pulse Amplitude Modulation where multiple bits are mapped to multiple levels. In this case there are 4 levels so it would be 4-PAM. (As Dilip points out in the comments that the levels in 4-PAM would typically be equally spaced while here the inner levels are spaced further than the outer levels, complicating the BER computation slightly as detailed below but still PAM regardless).

The demodulation scheme would similar to BPSK demodulation in that you would multiply the received signal with the recovered carrier and low pass filter. Decision would be made on the output by setting decision thresholds half way between each of the levels (assuming equiprobable data) and then the derivation for BER follows similar to BPSK (solving for the tail probabilities under Gaussian distributions if the noise is assumed to be AWGN as is typically assumed, noting the different tail contributions involved since the levels aren't equally distributed as Dilip has pointed out - that said the thresholds at ±1.5 will likely dominate setting the error rate with significantly smaller error contribution from the threshold at 0).

With equiprobable distribution of symbols and AWGN with standard deviation $\sigma$ this would result in the following Symbol Error Rate:

$$\text{SER} = \text{erfc}(0.5/(\sqrt{2}\sigma) + \frac{1}{2}\text{erfc}(1/(\sqrt{2}\sigma) $$

For smaller $\sigma$ values the second term would likely be insignificant. For larger $\sigma$ there would be additional tail contributions from the next adjacent threshold levels which would also significantly complicate the BER computation as detailed below. Thankfully under all practical applications these tails would certainly be insignificant and can be safely ignored.

To go from Symbol Error Rate to Bit Error Rate, consider how many bit errors occur when a symbol error occurs. Because the mapping done by the OP was properly Gray-coded (only one bit change between adjacent symbols), only one bit error would occur with each symbol error.


If you have two DFTs $A[k]$ and $B[k]$ (note the correct representation of a sinusoid at DFT bin number $1$)

A = [0,-j,0,0,j]; B = [1,1,1,1,1];

with the corresponding time-domain sequences $a[n]$ and $b[n]$

a = ifft(A);    % [0, 0.38042, 0.23511, -0.23511, -0.38042];
b = ifft(B);    % [1,0,0,0,0];

then the multiplication of the time-domain sequences $c[n]=a[n]b[n]$ corresponds to the cyclic (or circular) convolution of the DFTs $A[k]$ and $B[k]$:

$$\text{DFT}\{c[n]\}=C[k]=\frac{1}{N}\sum_{n=0}^{N-1}A[n]B[k-n]_{\text{mod} N},\quad k=0,1,\ldots,N-1\tag{1}$$

where the indices are taken modulo $N$ ($N=5$ in this case), such that they remain in the range $[0,N-1]$:

$$\begin{align}C[0]&=\left(A[0]B[0]+A[1]B[N-1]+\ldots+A[N-1]B[1]\right)\cdot\frac{1}{N}\\
C[1]&=\left(A[0]B[1]+A[1]B[0]+\ldots+A[N-1]B[2]\right)\cdot\frac{1}{N}\\&\vdots\\
C[N-1]&=\left(A[0]B[N-1]+A[1]B[N-2]+\ldots+A[N-1]B[0]\right)\cdot\frac{1}{N}
\end{align}$$

It's easy to see that the cyclic convolution of $A[k]$ and $B[k]$ as given above results in the zero vector, because due to $B[k]=1$ for all values of $k$, each element of the result is the sum of all elements of $A[k]$, which is always zero.


Why do we use the HSV colour space so often in vision and image processing?

The simple answer is that unlike RGB, HSV separates luma, or the image intensity, from chroma or the color information.  This is very useful in many applications.  For example, if you want to do histogram equalization of a color image, you probably want to do that only on the intensity component, and leave the color components alone.  Otherwise you will get very strange colors.

In computer vision you often want to separate color components from intensity for various reasons, such as robustness to lighting changes, or removing shadows.

Note, however, that HSV is one of many color spaces that separate color from intensity (See YCbCr, Lab, etc.).  HSV is often used simply because the code for converting between RGB and HSV is widely available and can also be easily implemented.  For example, the Image Processing Toolbox for MATLAB includes functions rgb2hsv and hsv2rgb.


As we know $E = Pt$ where $E$ is the energy and $t$ is the time. In the continuous domain $t$ is never ending and its beginning isn't known either.

Now, signals in the continuous domain are voltages and as you know from Ohm's Law $V = IR$ and $P = I^2R$.

This basically tells you that for a given resistance $R$, the energy of a signal over a time interval is given by the integration of the power. So,

$$ E = \int\limits_{t_1}^{t_2} P(t) ~dt$$

which in this case becomes

\begin{align}
E &= \int\limits_{-\infty}^{\infty} I^2(t)R ~dt\\
&= R\int\limits_{-\infty}^{\infty} I^2(t) ~dt\\
&= \int\limits_{-\infty}^{\infty} I^2(t) ~dt~~\text{(For $R = 1$)}
\end{align}

Usually in SP world we like to represent the signal as $x(t)$ giving us

$$E = \int\limits_{-\infty}^{\infty} \left|x(t)\right|^2 ~dt$$

Hope this helps!


Deconvolve a symmetric sequence into minimum-phase and maximum-phase sequences

For approximating a minimum phase sequence from a finite length linear phase sequence (such as a windowed Sinc), here are some existing stackoverflow and dsp.stackexchange answers regarding a matlab minimum phase approximation script which uses cepstrum deconvolution:

https://stackoverflow.com/questions/16753508/cepstrum-deconvolution-matlab-code

Derive minimum phase from magnitude

And here's an earlier post on the same matlab script, with more comments, from the music-dsp mailing list:

http://www.music.columbia.edu/pipermail/music-dsp/2004-february/059372.html

And the code:

x = [x; zeros(size(x))];
m = size(x,2);
n = size(x,1);
wn = [ones(1,m); 2*ones(n/2-1,m); ones(1,m);  zeros(n/2-1,m)];
y = real(ifft(exp(fft(wn.*real(ifft(log(abs(fft(x)))))))));
y = y(1:n/2;);

I've found that a greater amount of zero=padding of the original linear phase sequence may help improve the approximation.


Remove AC component from the mixed signal using an elliptic filter

I agree, the 1 Hz signal should be attenuated given the filter response. I suggest sweeping the input frequency to create the frequency response manually as I assume the filter was somehow scaled, possibly by not including a factor of $T$ in the mapping from the analog filter. Try a much higher frequency to see where it is nulled and then narrow in on the actual cut-off.

Further notice that the frequency response is not completely nulled at 1 Hz, but would have about 65 dB attenuation and that it also contributes significant phase distortion to the passband.

Instead of using an elliptic lowpass filter, consider using a 2nd order notch filter for eliminating AC signals. This implementation is detailed further at this post:

Transfer function of second order notch filter


Just to make things clear - this property is not fundamental but important. It is the fundamental difference when it comes to using DCT instead of DFT for spectrum calculation.

Why do we do Cepstral Mean Normalisation

In speaker recognition we want to remove any channel effects (impulse response of vocal tract, audio path, room, etc.). Providing that input signal is $x[n]$ and channel impulse response is given by $h[n]$, the recorded signal is linear convolution of both:

$$y[n] = x[n] \star h[n]$$

By taking the Fourier Transform we get:

$$Y[f] = X[f]\cdot H[f] $$

due to convolution-multiplication equivalence property of FT - that is why it's so important property of FFT at this step.

Next step in calculation of cepstrum is taking the logarithm of spectrum:

$$Y[q] = \log Y[f] = \log \left( X[f] \cdot H[f]\right) = X[q] + H[q]$$

because: $\log(ab) = \log a +\log b $. Obviously, $q$ is the quefrency. As one might notice, by taking the cepstrum of convolution in time domain we end up with the addition in cepstral (quefrency) domain.

What is the Cepstral Mean Normalisation?

Now we know that in cepstral domain any convolutional distortions are represented by addition. Let's assume that all of them are stationary (which is a strong assumption as a vocal tract and channel response are not changing) and the stationary part of speech is negligible. We can observe that for every i-th frame true is:

$$Y_i[q] = H[q] + X_i[q] $$

By taking the average over all frames we get

$$\dfrac{1}{N}\sum_{i} Y_i[q] = H[q] + \dfrac{1}{N}\sum_{i} X_i[q]$$

Defining the difference:

$$\begin{array}
&R_i[q] &= Y_i[q] - \dfrac{1}{N}\sum_{j} Y_j[q]\\
& = H[q] + X_i[q] - \left(H[q] +  \dfrac{1}{N}\sum_{j} X_j[q]\right) \\
& = X_i[q] - \dfrac{1}{N}\sum_{j} X_j[q]\\
\end{array}$$

We ending up with our signal with channel distortions removed. Putting all above equations into simple English:

Calculate cepstrum
Subtract the average from each coefficient
Optionally divide by variance to perform Cepstral Mean Normalisation as opposed to Subtraction.

Is Cepstral Mean Normalisation necessary?

It's not mandatory, especially when you are trying to recognise one speaker in a single environment. In fact, it can even deteriorate your results, as it's prone to errors due to additive noise:

$$y[n] = x[n] \star h[n] + w[n] $$

$$Y[f] = X[f]\cdot H[f] + W[f] $$

$$\log Y[f] = \log \left[X[f]\left(H[f]+\dfrac{W[f]}{X[f]} \right) \right] = \log X[f] +\log \left(H[f]+\color{red}{\dfrac{W[f]}{X[f]}} \right)$$

In poor SNR conditions marked term can overtake the estimation.

Although when CMS is performed, you can usually gain few extra percent. If you add to that performance gain from derivatives of coefficients then you get a real boost of your recognition rate. The final decision is up to you, especially that there are plenty of other methods used for the improvement of speech recognition systems.


You are right, DC means zero frequency. The term probably originates from electrical engineers' use of the terms DC and AC.


$n\textrm{ Hz}$ waveform sampled at $m\textrm{ Hz}$ per second

The basic square wave function signal.square() has a period of 2*pi. Multiplying each time by 2*pi*f is a way of stretching the wave along the time axis to get a wave with a period of 1/f.
To understand what the sample rate refers to, look at the first line:

t=np.linspace(0,1,500,endpoint=False)

If you look at the numpy documentation for this function, you'll see this creates a list of 500 evenly spaced numbers between 0 and 1. This is an array of times. When you then call the signal.square() function with this array, it returns an array of points on the square wave for each of those times. Thus you have sampled the square wave at 500Hz for 1 second.


Image getting completely white on downsampling

When you open an image using imread, you get the image as a uint8 matrix with grayscale (or color) values of 0-255. However, filtering is no possible for this type. Therefore, filter2 converts your image to double.
The function imshow accepts both uint8 and double type. When the input is uint8 it also assumes color values of 0-255. But, when the input is double, it assumes the values are normalised to the range of 0-1. Therefore you should tell imshow that although your input is of type double it has color values of 0-255 by using

imshow(down1,[0 255]);


MATLAB: Compute pole-zero diagram of FIR using only impulse response

Since your response is finite, your system corresponds to a FIR filter, where the coefficients are given by the values of the impulse response. Hence, you can perform Z-Transform of the FIR and find its zeros. Finally, since you would just need to calculate the zeros of the nominator polynomial. See this MATLAB code, following your impulse response example:

h = [1 0 -2 0 0 3];
b = h; % the FIR coefficients of a filter

roots(b)

zplane(b,1);

Output:

ans =

  -1.6379
   1.2438 + 0.4498i
   1.2438 - 0.4498i
  -0.4248 + 0.9309i
  -0.4248 - 0.9309i

So, the zeros are directly given by the roots. For the poles, you have (N-1) poles at the origin.

Also, have a look at e.g. this answer.


Sampling fast sinusoidal signal with low sampling rate

In practice it sounds difficult to measure the amplitude of a 1 MHz signal (or some other very high frequency) with a 1 kHz sampler; in general, the 1 kHz analog-to-digital converter will probably rely on the input voltage being stable for some small fraction of the 1 ms sampling period; if the signal is oscillating 1000 times within that period, it's likely that the sample-and-hold will effect some smoothing/low-pass filtering, or that the conversion may be completely wrong because the signal changed so much during the conversion process.

Assuming you really do have some very fast sample-and-hold so you can record the level of your high-frequency signal "instantaneously", then your best strategy might be to sample at random times (preferably without any quantization of the sampling instant).  If you sample at a fixed rate, you run the risk of sampling at an exact multiple of the high-frequency signal's period.  So if you were taking a snapshot of a 1 MHz signal every 1 ms (i.e., every 1000 cycles) -- so your signal periods were exactly 1000:1 -- you'd see the same part of the cycle every time, and you'd have no idea what the actual waveform looked like.

However, this is statistically unlikely.  In general, when undersampling a sinusoid, you see what looks like a sinusoid with frequency $\Omega_{alias} = |\Omega_{true} - r.\Omega_{sampling}|$ rad/sec, where $r$ is the integer that makes the result of the expression less than or equal to the nyquist rate, $\Omega_{sampling}/2$ .  Except in the case where some value of $r$ makes this evaluate to zero (aliased signal appears as a constant value), you'll sample a sinusoid whose amplitude follows the true amplitude of the underlying high-frequency signal.  But of course you have no idea what the value of $r$ actually is, so you don't know the true frequency.

In general, if you want to know the peak amplitude, you're going to need to observe for a complete cycle of $\Omega_{alias}$, which could be arbitrarily small (or zero, as discussed).  But assuming it's somewhat uniformly distributed, the chance of getting an $\Omega_{alias}$ smaller than 1 Hz with 1 kHz sampling is only 1 in 500 (because we have to include -1..1 Hz, i.e. 1/500th of the 1 kHz interval), so for almost every case we'd be fine just taking the largest value over 1 sec.  You could reduce the chance of error by increasing the observation period.

Or you could have multiple, irrationally-related sampling periods, on the assumption that even if one of them happens to be close to an exact subdivision of the high-frequency period, the others will not.  There will always be a frequency that is the product of each of the individual sampling frequencies that they are all exact submultiples of, but you can make that pretty large if you have more than a few sampling systems.


High pass filter not attenuating signal at start of input

Any filter with such a low cutoff (as compared to the sample rate) will lots of time domaing ringing.
The best way of adressing this is to properly initialize the state of the filter using known properties of your system & signal, but that's fairly tricky.
In your case, there may be a shit cut though: I'm guessing that you have a large bias or DC offset in your original signal. Just subtract the mean (or maybe the mean over the first samples or so, if the bias is drifting). This should subtantially reduce the ringing.


Smoothing power spectrum by convolution with boxcar function

The spectrum is periodic so you need circular convoluioint, not linear convolution
You are smoothing the  complex linear spectrum. Due to the phase of the complex ammplitude you will get significant interference in your smoothing window. Smoothing is typically in the energy domain and the phase is discarded. Complex smoothing is fairly difficult to do.


Your question is broadly answered here: How to extract traffic signs from a photograph?

Your specific case is a little more complicated because the signs are not particularly distinguishing, however the same suggested techniques of cross correlating the images with clean images of original signs (you will probably want a library of images containing all signs for best results) are still applicable.

One suggestion above was to use an infrared camera which will be very effective if the road signs are micro prismatic or glass bead designs.


Why is the total signal response response of zero input + zero state, $y(t) = y_0(t) + h(t)\star x(t)$, not an LTI?

The response of the system is given by

$$y(t)=h(t)\star x(t)+y_0(t)\tag{1}$$

where $x(t)$ is the input signal, $h(t)$ is some well-behaved function, $y_0(t)$ is a function which is independent of $x(t)$, and $\star$ denotes convolution. The model $(1)$ can be used for describing systems with non-zero initial conditions that would be LTI were it not for the non-zero initial conditions (such as an initially charged capacitor).

For $y_0(t)=0$, $h(t)$ would be the impulse response of the corresponding LTI system. For $y_0(t)\neq 0$, $h(t)$ is not the impulse response, neither can the system be characterized by an impulse response.

For $y_0(t)\neq 0$, the system described by $(1)$ cannot be linear, which is easily shown by considering a response $y_1(t)$ to an input signal $x_1(t)$. For a linear system, the response to a scaled version of $x_1(t)$ must be a scaled version of $y_1(t)$, with the same scaling factor. So the response to $ax_1(t)$ must be $ay_1(t)$ if the system is linear. But from $(1)$ we have

$$y_1(t)=h(t)\star x_1(t)+y_0(t)$$

and for the scaled input we get

$$y(t)=h(t)\star ax_1(t)+y_0(t)=ay_1(t)+(1-a)y_0(t)\neq ay_1(t)$$

for $y_0(t)\neq 0$ and $a\neq 1$. Consequently, the system described by $(1)$ is not linear.

In a similar way it can be shown that the system is also not time-invariant. The reason is again the function $y_0(t)$, which is independent of the input signal $x(t)$. If $x(t)$ is shifted, the output doesn't shift accordingly, because $y_0(t)$ doesn't shift.


How can I detect when a signal deviates from a flatline?

This type of problem is part of Detection Theory.  Do you have an estimate of the background noise levels?  How about the statistics of the noise?  If the noise is AWGN, things get a little easier.  You can set a detection threshold based off your noise statistics to achieve a desired probability of detection / probability of false alarm.

Be careful about using the mean of your data to set the threshold.  The "spikes" could bias your mean estimate.  If the noise is not stationary (i.e. it changes levels), using the mean could result in false detections (alarms).  One simple option is to compute a running average and create a threshold off this running average.  If the sample of interest is above the threshold, declare a detection and do not include the sample in the running average calculation.  There're a ton of variations to this simple concept depending on the situation.

If your signal is truly a flat line without noise, you could apply a highpass filter by subtract the previous sample from the current sample.  The result would be an array of zeros with non-zero sample values at the spikes.

For a more formal treatment of detection theory, I highly recommend Steven Kay's book "Fundamentals of Statistical Signal Processing, Volume II: Detection Theory".


Bandpass filters with python for low frequencies

I would bet this is just numerical error in the transfer function.  Try using butter_sos = butter(..., output='sos') instead of ba format, and sosfreqz(butter_sos, ...) instead of freqz.  Does that fix it?

https://docs.scipy.org/doc/scipy-0.19.0/reference/generated/scipy.signal.sosfreqz.html

Edit: I just tried it and it does.  :D

import scipy
import numpy as np
from scipy import signal
from matplotlib import pyplot as plt
def butter_bandpass(lowcut, highcut, fs, order=5, label=None):
    nyq = 0.5 * fs
    low = lowcut / nyq
    high = highcut / nyq
    sos = signal.butter(order, [low, high], btype='band', output='sos')
    w, h = signal.sosfreqz(sos,worN=20000)
    plt.semilogx((fs * 0.5 / np.pi) * w, abs(h), label=label)
    return sos

center_freqs = np.array([0.5* 2 ** (n * 1 / 3.) for n in range(0, 29)])
center_freqs.sort()
lower_freqs = 2. ** (-1 / 6.) * center_freqs
for lf in lower_freqs:
    butter_bandpass(lf, lf*2**(1/3.), fs=1000, order=5)

plt.show()


Does 1 kHz sine tone means $\sin(2(1000)\pi t)$ or $\sin(2(500)\pi t)$?

The trigonometric functions "do not know" what a Hertz is and they do not care either. The only thing they know is that a full circle is $2 \pi$ radians. Whether this circle concludes in days, hours, picoseconds or a slice of it represents the angle a force is applied to some lever, is immaterial.

$2 \pi \omega$ expressed in Hertz, denotes a rate. A rate of going around a circle at the time span of a second. $y = \cos(2 \pi 1 t)$ where $t$ is in seconds, would have concluded 1 circle, composed of $2 \pi$ radians, by the time $t$ ticks to 1.

To make it conclude the circle faster, we multiply the "passing of time" (denoted by $t$) by some number $f$.

Therefore, a 1kHz tone is $2 \pi 1000$ radians per second.

Hope this helps.


Maximum Likelihood Through a Noisy Channel

I'm inclined to think this is true, but so far I've only gotten a counterexample:

Consider the channel:

\begin{bmatrix}0.5 & 0.5 & 0 & 0\\ 0 & 0.5 & 0.5 & 0 \\ 0 & 0 & 0.5 & 0.5 \\ 0.5 & 0 & 0 & 0.5\end{bmatrix}

This a typewriter channel. This matrix is not invertible (it has zero determinant, and rank 2). So the ML estimate of $\Theta$ cannot be gotten from this channel.

However, for $p \neq 0.5$ I think it could well be true.

EDIT: In fact I think the noisy channel aspect of this may well be a red-herring. The ML estimator of $\Theta$ is as you describe, simply because you have a linear model. You don't need any conditions of $C$ other than it should be invertible (i.e. all rows are linearly independent).


Laplacian Filter Implementation in MATLAB

You should use conv2() or imfilter() with your filter.
Unless you want to implement it by hand,


Why Runge-Kutta for Quaternion integration in Kalman filter?

I think the confusion comes from the authors not parameterizing things clearly.  Furthermore, by switching to geometric algebra rather than quaternions, some additional confusion can be cleared up.

The main difference between normal vector algebra and geometric algebra is that we can multiply vectors. So of $e_x$, $e_y$, and $e_z$ are our (orthonormal) basis vectors, we also have $e_x e_y$ which is not a vector but a bivector which can be thought of as an oriented plane element.  (By "pure vector quaternion" they mean a bivector; you are correct that that means $0$ $w$ part and it's reasonable to say it has a $0$ "real" part.)  A key property is that the basis vectors anti-commute with each other, i.e. $e_x e_y = -e_y e_x$.  This leads to $(e_x e_y)^2 = -1$.  In general, any unit bivector squares to $-1$, which means we can apply Euler's formula: $$R = e^{\theta B} = \cos\theta + B\sin\theta$$ where $B$ is a unit bivector and we call $R$ a rotor.  A unit quaternion is just a 3D rotor.  (Complex numbers are just the even sub-algebra of the 2D geometric algebra.  Here we are looking at quaternions, the even sub-algebra of 3D GA.)

If we want a rotation by $\theta$ in the plane of a unit bivector $B$ we use the rotor $R = e^{\frac{\theta}{2}B}$.  (The half comes in because we rotate a vector $v$ via $Rv\widetilde{R}$; see the link above.)  Of course, we want to allow both the plane and the angle to vary so we define the (non-unit) bivector $\Theta(t) = \theta(t)B(t)$ so if $$R(t) = e^{\frac{1}{2}\Theta(t)}$$ then $$\dot{R}(t) = \frac{1}{2}\Omega(t)R(t)$$ where $\Omega = \dot{\Theta}$.  This is their $f = q\omega/2$; their $q$ is our $R$ and their $\omega$ is our $\Omega$.  Now if we have $R(t_0)$ and we want $R(t_0 + \Delta{}t)$ we need to integrate $\dot{R}$ from $t_0$ to $t_0 + \Delta{}t$.

The scheme you describe is roughly analogous to doing forward Euler integration, which essentially assumes their $q$ and $\omega$ our $R$ and $\Omega$ is constant over $\Delta{}t$.  RK4 is just a better integration method.  I doubt there is any special reason they chose to use RK4 as opposed to other integration methods.  It's just the default choice typically.


"Perfect waveform" for cross-correlations?

In general, the resolution of the cross correlation is proportional to the bandwidth of the signal. This is way pseudo-random signals are common for that purpose.
So if you are looking for a signal that will produce sharp and distinct cross-correlation, you should use such signals.

However, if you are looking to improve the time-delay estimation for a given signal, you should pre-filter the signal with high-pass filter as the lower frequencies tend to produce more smooth cross correlation.


How to detect start and finish of temperature control in temperature time series

One option: It's pretty clear that in the uncontrolled range the temperature follows a  simple exponential heating/cooling curve. Something like
$$T_{n+1} = T_n + (T_{Ambient} - T_n) \cdot e^{-\Delta t / t_{Building}} $$
where $\Delta t$ is the time difference between samples and $ t_{Building}$ the thermal time constant of your building (which may be slightly different at different sensor locations).
You could simply test the difference of each new sample to the model prediction and set an threshold of the allowable error. The only two model parameters are ambient temperature and thermal time constant of the building.
The thermal time constant is a function of the thermal resistance of the building and it really doesn't change unless someone opens/closes windows or there is some actual structural work happening. Looking at your data that should be easy to estimate. Ambient you can try to get directly either with a sensor or looking at a reliable weather report. You can also just use the data from the first few sensor to estimate ambient and then test whether the other sensors agree or not.


Ideal low-pass filter in the Frequency Domain

The rule should be something like this i suppose $$D(|U - U_o/2|,|V - V_o/2|)≤ D0$$ where $U_o$ is the length of image in horizontal dimension and $V_o$ is the length of image in vertical dimension


The scheme is not universal.

One reason to do it is to allow the system to be modeled around a set-point or bias level ($u_0$) without having that set-point affect the state space model.

Subtracting the set-point, doing the state space equations, and then adding back on the corrected set-point ($y_0$) lets you do this.


simulator for dsp kit tms320c6713?without using actual hardware kit tms320c6713?

Usually, even before someone builds a processor, they write a specification of the instruction set, and then they build a simulator for that.
The same applies to TI's C67xx series: https://processors.wiki.ti.com/index.php/List_of_Simulator
The wiki article strangely says that TI is moving away from simulation in Code Composer Studio 6, which would honestly be a terrible, terrible, terrible idea and put them back in the 1980s with respect to software engineering practices (no chance to do unit tests on any machine).
The only reason for doing so I could think of is lack of manpower to maintain that simulator, and/or a plan to deprecate the C67xx series altogether.


Why Is The Maximum Likelihood Estimation (MLE) Method Taken as the Benchmark for Comparing with Other Methods?

There are 3 major reasons in my opinion which makes the Maximum Likelihood Estimator so popular:

Intuition
It is very clear what's the logic behind this method and what you maximize. It makes sense even if you have little knowledge in probability and statistics.
Properties
It has great properties (Some are only asymptotically, namely with many samples) such as consistency and efficiency. This means it always makes some sense to use it.
Easy to Derive
It is usually pretty easy to derive the optimization problem defined by the MLE and in some cases it is also simple to solve it.

Now, many parameter estimation methods are compared to the properties of the MLE - Efficiency and Consistency.
If they exhibit similar performance than they are as optimal as the MLE.
Otherwise, they are not, as expected.


Can you quarter the processing time for real, symmetric FFTs?

I'm not sure about quartering the process time but this is what I can think of for now. N point DFT of x.

$X_k = \sum_{n=0}^{N-1}x_nW_N^{kn}$

By the even symmetry assumption

$X_k = \sum_{n=0}^{N/2-1}x_n(W_N^{kn}+W_N^{k(N-1-n)})$

Then we simplify

$X_k = W_N^{k(N-1)/2}\sum_{n=0}^{N/2-1}x_n(W_N^{kn}W_N^{-k(N-1)/2}+W_N^{-kn}W_N^{k(N-1)/2})$

$X_k = 2W_N^{k(N-1)/2}\sum_{n=0}^{N/2-1}x_n\cos(2\pi k\frac{n-(N-1)/2}{N})$

There are probably some errors but if simplifications are valid the 'FFT' size is reduced by a factor of 2 and the twiddle factors are real (this maybe corresponds to your quartering). As in normal FFT the symmetry and periodicity can be exploited to decompose the 'symmetric FFT' into radix-2 steps for instnace.


What helped you in reading and implementing papers and coming up with original ideas in computer vision/graphics?

People working in computer vision/graphics may have very different backgrounds: computer programming, high-performance computing, signal processing, image analysis, optimization, machine learning, etc.

Novel contributions can be wide: from a very fast GPU/CPU algorithm to render textures (I had a reference in mind which I can't recall now) to semantic segmentation. Origins of  ideas may range from huge expertise, trials on billion of images to profound thoughts, through hardware algorithmic wisdom. Your question is far-reaching, and can dig into the depth of the psychology of scientific invention, well-described by some authors in mathematics and physics (I have Raymond Poincaré and Jacques Hadamard in mind). Here are a couple of thoughts.
The ability to (non exclusive):

be and remain curious: listen a lot, don't pretend too much, don't trust doxa (everybody knows - even yourself - that..., nobody does this), look at novelties in other fields
test and probe:hardware, data, hypothesis, models... to the limit (add too much noise, down-sample a lot, reduce learning set, cut the power),
explore unconventional ideas: transfer knowledge from one domain to the other, ask "what if I do..." (I got from Isaac Asimov that most of  science is not like Archimedes' "Eureka (I found it)",but often "Ooh, that funny"), learn XXX that strange theory, go to seminars you are not truly interesting in, or above your level
have people you trust to talk to: especially those who disagree, or dare to say "you're going nowhere", hire good colleagues, better students
take notes, talk to yourself, don't be discouraged: nice works can take decades
teach: basic student questions can drive you nuts and shake your confidence
corollary: have some solid knowledge somewhere (law of excessive learning: you can teach at one level only if you know the level above)
breathe and sleep: get away from the machine, take some fresh air, wander in woods
be patient, but remain aware of your progress: one gain confiance progressively. The first "original" ideas I got were already discovered 30 years before. I'm grateful now went I read a 1-year paper that implements what I am thinking of.
have principles: it is OK to think "This should work!" even if the evidence is shallow. Think about the "AI winter".

Those can be good incentives. On my side, I am lucky to (think that I) have a relatively broad "generic" scientific background (to connect far-fetched domains) and one or two technical skills that I dug deep into. Thus, I am not afraid of switch domains. One main risk I see is siloing: when collaborating with other, it is tempting to (and thus, try to avoid):

remain in our confidence zone,
trust others in their domain (as we would like to be treated, with respect).

One of my few scientific prides (not deep, but quite intriguing) was to work with real-time computing people (my programming skills are very limited), and discover that some "dogma" in the domain was false. My contribution was initially very small, on standard extrapolation, which I though I mastered. Some colleagues trusted me and implemented the stuff. But one colleague who did not understand the derivation asked me again, and again, and again, to understand every step that I thought was "evident". So I reduced the formulas down to the bone. With a couple of image compression tricks, we ended up with a scheme that shook the initial dogma in cyber-physical simulation. Not theoretically, but enough to feel we had "something new".


How to calculate Long Term Average Spectrum of an audio file?

Sounds like a conventional spectral analysis, just do an FFT over your entire  audio file and you obtain an averaged spectrum of the file. A Fourier transform can be viewed as a set of band-pass filters with integration.


How to evaluate expressions without explicitly computing a DFT?

Parts of the answer.

I) is the average value (DC coefficent) of the signal, 2.2

III) is, by Parzeval's theorem, equal to the square sum of the signal, 80


What is difference between outputs of Fourier transform and Fourier series of a periodic square waveform?

The Fourier transform $X(\omega)$ of a $T$-periodic function $x(t)$

$$ x(t+T) = x(t) \quad \quad \forall t $$

having complex Fourier coefficients $c_n$

$$ x(t) = \sum_{n=-\infty}^{\infty} c_n e^{j 2 \pi n t/T} $$

$$c_n = \frac{1}{T} \int_{0}^T x(t) e^{-j 2 \pi n t/T} \ dt \tag{1}$$

can be expressed as a weighted sum of Dirac impulses, where the weights are given by the complex Fourier coefficients:

$$X(\omega) \triangleq \mathcal{F}\left\{ x(t) \right\} = 2 \pi \sum_{n=-\infty}^{\infty} c_n \delta\left(\omega - \frac{2\pi n}{T}\right) \tag{2}$$


Signal analysis: selection of oscillatory parts

By looking at that graph, it seems that simple math (or rather statistics) can be used to solve this problem. I'd divide (segment) your data into equal time frames and then compute the mean and possibly variance (or standard deviation) for each time frame.

The idea here is that the parts you're interested in will have a greater mean (and a greater variance) from the rest of the data. You'll still need a threshold, yes, but this approach is much more robust than just trying to compare the individual data points against the threshold.

For best results, you'll probably need to experiment a bit with different time frame lengths as well as different threshold values.


Doubt on behaviour of frequency responce of a system

It appears to have -40 dB/decade roll-off in magnitude starting at 1 rad/sec and transitioning to -120 dB/decade roll-off in magnitude at approximately 300 rad/sec (this suggests you have 2 dominant poles that begin to be dominant when f is at 1 rad/sec and 4 additional poles that begin to contribute at f> 300 rad/sec since each pole would contribute -20 dB/decade to frequency roll-off.

In comparison, red is rolling off at -80 dB/decade consistent with 4 poles.  The poles would actually be complex given the frequency peaking in all cases.


Local Image Features Extraction for Road Signs

Do you know what signs you are looking for? If yes, maybe you could do template matching (e.g. a normalized cross-correlation, available in matlab). It won't work great when signals are getting closer, since the perspective projection will change their appearance, but it should work for mid-range detection.

You can limit your template-matching search to the right side of the road.

There are also template matching algorithms for affine transformations that might be a next step.

Christos


The (linear or aperiodic) convolution of two vectors $\mathbf x = (x[0], x[1], \ldots, x[N-1])$ and $\mathbf y = (y[0], y[1], \ldots, y[N-1])$ is a vector $$\mathbf z = {\mathbf x}\star {\mathbf y} = (z[0], z[1], \ldots, z[2N-1]).$$
On the other hand, their cross-correlation is a vector
$$\mathbf w = {\mathbf x}\otimes {\mathbf y} = (w[-(N-1)], w[-(N-2)], \ldots,
w[-1], w[0], w[1], \ldots, w[N-2], w[N-1]).$$
Now, the FFT method for computing convolutions computes the cyclic (a.k.a. periodic
or circular) convolution of the two vectors.  To compute the linear or
aperiodic convolution via FFTs, you zero-pad the
vectors to length $2N-1$ so that the cyclic convolution of the zero-padded
vectors is exactly the linear convolution you are looking for.
Thus, after FFTing, term-by-term multiplying, and inverse-FFTing,
what you get is $\mathbf z$. This is exactly what you have done as per
your question. However, when you use the FFT method to compute the
cross-correlation, what you get is the cyclic cross-correlation of
the two zero-padded vectors, and this cyclic cross-correlation is not
$\mathbf w$ but instead
$$\hat{\mathbf w} = (w[0], w[1], \ldots, w[N-1], w[-(N-1)], w[-(N-2)], \ldots,
w[-2], w[-1])$$ so that to get the $\mathbf w$ that you want, you need
to circularly rotate $\hat{\mathbf w}$ so that $w[0]$ appear in the middle
of the vector.


Phaser effect with feedback control IIR filter calculation

Solution can be found in the well-known paper of Vadim Zavalishin. In sections (4.1) and (4.2) there's an explanation on how to implement zero feedback loop around cascade of LP filters. He even mentions later in the paper that exactly the same approach can be used to implement phaser (6.1)


Z-transform of an impulse signal in discrete time

There's a ROC (region of convergence) problem in your derivation.

Consider a decomposition of the signal as $x[n] = x_1[n] + x_2[n]$. Then in order for the $Z$ transform of $x[n]$ to exist, the corresponding $Z$-transforms of both $x_1[n]$ and $x_2[n]$ must exist for a set of $z$ which is the intersection of the individual ROCs of the components $x_1[n]$ and $x_2[n]$. If the intersection set is empty then the Z transform does not exist. Applying this to your signal yields.

$$x[n] = \sum_{k=-\infty}^{\infty} \delta[n-k] = u[n] + u[-n-1] = 1 ~~~~,~~~\text{for all n}$$

Where the individual $Z$-transform and their corresponding ROCs are
$$ u[n] \longleftrightarrow \frac{1}{1-z^{-1}} ~~~,~~~ ROC_1 ~~~ |z| > 1 $$
$$ u[-n-1] \longleftrightarrow \frac{-1}{1-z^{-1}} ~~~,~~~ ROC_2 ~~~ |z| < 1 $$

As you can see the intersection of $ROC_1$ and $ROC_2$ is empty; i.e; there is no set of values $z$ for which both of the components has convergent $Z$ transforms. Hence we conclude that the $Z$ transform $X(z)$ of the signal $x[n]=1$ for all $n$ does not exist.


FloatFFT has a funny array layout. Real and Imaginary components of each datapoint/bin have to be interleaved in the data array you provide.
The FFT plot you show it is not a common FFT plot: is the real and imaginary part of the FFT interleaved: $[ Re\{X[0]\}, Im\{X[0]\}, Re\{X[1]\}, Im\{X[1]\}, ... ]$.

Handle two consecutive points everywhere accordingly. For example: $bin[0] = abs( sqrt( a[0]*a[0] + a[1]*a[1] )$ and so on... Also, follow that layout when you perform the inverse.

Documentation is not the best and you will have to read it carefully.

Bonus: the reason you have to provide a $2 \cdot n$ array to some functions is not that the FFT has to double the data set length, it is because each bin is stored in two consecutive array positions.


Removing noise from audio using Fourier transform in Matlab

From my experience, your best bet is spectral subtraction. You can find a lot of literature about it on the web. Here just a quick overview of what you need to do. You take overlapping windowed blocks of your time domain signal and transform them to the frequency domain using an FFT. Then for each frequency bin you need to estimate the signal-to-noise ratio. There are several simple noise tracking algorithms that perform well if the noise is relatively stationary. Based on the estimated SNR per frequency bin you mulitply each bin with a gain constant between 0 (terrible SNR) and 1 (no noise). Note that you only deal with magnitudes, the phase of the FFT is left unchanged. After modifying the FFT bins with the gain factors, you combine the processed bins back to the time domain using the phase of the original noisy signal.

Here is a classic paper on this subject. It goes a bit further than pure spectral subtraction but the overall system is the same. Don't worry too much about the math involved to derive the optimal gain function.


How to describe correlated noise after the signal is oversampled?

Assume that the input to a simple up-sampler (say by $L$) is $x[n]$, hence
$$y[n] = \begin{cases} { x[n/L] ,\qquad n=kL \\ 0 , \qquad \qquad \ \text{  else } }\end{cases}$$
If you compute the auto-correlation sequence, defined as
$$\phi_{xx}(i,j) = \mathbb{E}(x[i]x[i+j])$$
for an iid white Gaussian process, we have that
\begin{align}
\phi_{xx}(i,j) &= \begin{cases} \sigma^2, \qquad & \text{if } j=0 \\
0, \qquad & \text{else } \end{cases} \\
\phi_{yy}(i,j) &= \begin{cases} \sigma^2, \qquad & \text{if } j=0 \text{ and }i = kL \\
0, \qquad & \text{else } \end{cases}
\end{align}
Unlike i.i.d Gaussian, where $\phi_{xx}(i,j)$ is determined by the lag $j-i$; $\phi_{yy}(i,j)$ is not fully determined by the lag. That said, $y[n]$ is not WSS, whereas $x[n] $ is. With respect to correlation, the quantities $\mathbb{E}(x[n]x[m]) = 0$ for $n \neq m$ is due to the uncorrelation assumption. You can say the same thing on $y[n]$, namely $\mathbb{E}(y[n]y[m]) = 0$ for $n \neq m$.


Obtaining power spectrum from ACF, FFT using Matlab and FFTW

I am also confused about the terms power spectrum, power spectral density and energy spectrum. Are they used interchangeably? I would be grateful for an explanation of this phenomenon and normalization procedure.

EDIT with @jojek`s helpful comment:

power spectrum is NOT the same as power spectral density despite the fact that terms are often mixed up. Units of power spectral density (PSD) are /Hz, e.g. /Hz. Power Spectrum (PS) is equal to PSD*(Equivalent Noise Bandwidth of used window funktion) and has units of  . You often normalise it with N, the number of data points (for discrete time series). For physical processes the PSD and PS are symmetric, thus negative frequencies contain no additional information and only positive frequencies are plotted, negative frequencies are discarded.

energy spectrum is a term of spectroscopy: it describes the intensity of a particle beam as a function of particle energy. What you maybe meant is called energy spectral density, which describes the distribution of energy over frequency; it has units of J/Hz

A good overview is here.

I want to add the term amplitude spectrum/amplitude spectral density, which is the square root of the power spectral density.


Interpreting the inverse fourier transform from a graph

If there are no negative frequencies, $x(t)$ would have an imaginary component (see http://en.wikipedia.org/wiki/Analytic_signal), but we can just assume that the graph was truncated for convenience. So I think you have got it right.

Note: previously I wrote in this answer that the coefficients should have been multiples $\frac{1}{2\pi}$ but they are correct if we assume the negative frequencies mirror the positive.


A first order lowpass filter is usually implemented like this:

$$p[n] = \alpha p[n-1] + (1-\alpha) pi[n]$$

Where $p[n]$ is your filtered power estimation, $p[n-1]$ is the previous result, $pi[n]$ is your new measurement (probably the product of instantaneous voltage and current measurements), and $\alpha$ is a positive parameter just less than 1.

The nearer $\alpha$ gets to 1, the larger the time constant (lower cutoff frequency) of your filter. But beware, especially in embedded systems with limited precision, that getting too near to 1 can make your filter unstable, or at least have problems due to numerical precision.

The cutoff frequency for that filter is around $f_s \frac{1-\alpha}{2\pi\alpha}$, where $f_s$ is your sampling frequency.


Convolution Offset and Factor in Frequency Domain

Those lines that keep the range of the R,G,B channels within the valid [0,255] range work on individual pixels. There is a comparison between the pixel value and 0 / 255 levels. Unfortunately, this operation requires the exact pixel value to be known. However the Fourier transform mixes all the individual pixel values into all frequency bins and unless explicitly (and inefficiently) inverted back for a given pixel, this value is not accessible.

So it's not possible to limit the individual pixel values in the frequency domain.


simulink discrete integrator block problem

All of a sudden I realized I made a silly mistake.
$$
\int \cos(\omega\cdot t)dt=\frac{1}{\omega}\sin(\omega\cdot t)
$$
which means the "integrator: Forward Euler" choice produced the correct result.


Signal approximation using linear combination of functions

The method you chose is a standard least squares approximation, i.e. you minimize the sum of squared errors on the chosen grid. You obtain your coefficients by solving an overdetermined system of linear equations in a least squares sense. I think that's a very sane approach.


What Are The Semantics Of Wav-File Sample Values?

I don't know much about this semantics? of WAV files but their numerical format is the following. (assuming mono)

Given a recording with 8-bit per sample precision, then those samples are unsigned integers taking values between $0$ and $255$. Due to being unsigned, to represent negative values,  there is a bias of $128$, and the sample values are actually interpreted to be between $-128$ and $127$. In this case a value of $128$ represents silence, maximum value $255$ corresponds to an input full scale positive audio waveform and $0$ represents neagitve full scale input waveform.
For $16$-bit WAV files, you have $16$-bit signed integers that range from $-32768$ to $32767$ and $0$ represents silence. Again, $32767$ represents full scale audio input without scaling.

Note that here I refer to a microphone recording through an ADC system which converts analog audio input in the standard range of about $1$ V rms into WAV samples. Please check this 1-V rms for being something higher, about 1.5 Volts or else depending on particular implementation or some specific audio standard. In any case, the analog circuitry is responsible for transmitting a full scale audio to the ADC whan a 1V (or 1.5V) rms wave is present at the input and then the conversion is performed according to the above setting when WAV recording is performed.


Order of applying two Butterworth filters (distributivity)

The short answer: it depends.

The longer version: From MIT OpenCourseWare notes:

  Specifically, if we have several LTI systems cascaded together, the
  output generated by an input to the overall cascade combination does not depend
  on the order in which the systems are cascaded.

In your specific case, you have a cascade of two blocks: one which is a 2nd order Butterworth low-pass filter (which can be shown to be a linear-time-invariant (LTI) system), and another which is your linear transformation function $f(a,b)$. The question then boils down to whether $f(a,b)$, in addition to being a linear transformation, also happens to be a time-invariant transformation.

If $f(a,b)$ happens to be an LTI transformation, then the cascade $\mbox{butter(f(a,b))}$ is equivalent to $f(\mbox{butter}(a),\mbox{butter}(b))$. Thus, in that case one could say that $\mbox{butter}$ is distributive over $f(a,b)$.

If $f(a,b)$ isn't time invariant, then we can't make that claim as to whether the result depends on the order of operations.


How to calculate time of arrival using room impulse response

If you do not have a very controlled room, there might be something attenuating the direct path and a reflection might then be the highest peak. In that case you should accept also a somewhat attenuated peak if it comes earlier than the strongest peak.

Coupling to structures is not that much of an issue. For example, when a sound wave reflects at the interface of wood and air, only a fraction of about 0.0007 (see Ref. 1) of the energy of the wave will be passed to the other medium. The wave will then spend a lot of time in the first medium until much of its energy is passed to the second medium, so you can imagine that the impulse response of this coupling is not very concentrated in time and doesn't disturb the peak detection much.

References:

Acoustics of wood, course in non-destructive testing of wood, ETSI Montes, ETS Arquitectura – Universidad Politécnica de Madrid, June 2005


Show that any continuous-time signal $x(t)$ can be represented as $x(t)= x_e(t) + x_o(t)$

Substitute the last two equations into the first and it should be clear.

iI.e., if you have $x_e(t), x_o(t)$  given to you. Substitute these into the equation for $x(t)$ and you should see that both sides are equal to each other.

(Note that this proves the equality, but doesn't prove that $x_e(t)$ is even or that $x_o(t)$ is odd.)


Consider your carrier signal as $E_c\cos W_{ct}$ and your message signal is $E_m \cos W_{mt}$.

The double side band plus carrier signal is,

$ E_c\cos(W_{ct})+(E_c\cos(W_{ct})E_m\cos(W_{mt})) $

$ E_c\cos(W_{ct})+\frac{1}{2}(E_cE_m[\cos(W_{ct}+W_{mt})+\cos(W_{ct}-W_{mt})]) $

$ Ec*Cos(W_{ct})+\frac{1}{2}E_cE_m\cos(W_{ct}+W_{mt})+\frac{1}{2}*Ec*Em*Cos(W_{ct}-W_{mt}) $

Carrier being $E_c\cos(W_{ct})$

Upper Side band being $\frac{1}{2}E_cE_m\cos(W_{ct}+W_{mt})$ and

Lower Side band being $\frac{1}{2}E_cE_m\cos(W_{ct}-W_{mt})$.

For a signal $A\cos W_t$, the average power is $(\frac{A^2}{R}T)*(\int^T_0 \cos^2 W_t) = \frac{{A_c}^2}{2R}.$

Hence,

Average Power of Carrier is $\frac{{E_c}^2}{2R}$

Average Power of Upper Side band is $\frac{{E_cE_m}^2}{8R}$

Average Power of Upper Side band is $\frac{{Ec*Em}^2}{8R}$


Increasing the number of bins in DFT with limited points and high samplerate

If all you want is a finer bin spacing, then the most straightforward thing to do is to increase $N$ by zero-padding. Note, however, that this will not increase your frequency resolution. Frequency resolution is limited by the total length of time that you observe the waveform. So, if you want to actually be able to resolve frequencies more finely, you need to observe the signal for a longer period of time.

Your first option (essentially interpolating the signal to a higher sample rate) won't do what you want at all; $N$ will be larger, but so will the sample rate of the interpolated signal. If you double the size of the dataset by interpolating by 2, then you'll have to do a DFT of twice the length in order to get the same bin spacing.

Downsampling the signal to a lower sample rate as you suggested in your third option can make sense if, as you said, you only care about frequencies below 300 Hz and your initial sample rate is 44.1 kHz. However, what I said above applies: this won't improve your frequency resolution at all; you'll need to either live with zero-padding to interpolate the DFT spectrum, or observe the signal for a longer period of time instead.


Reducing standing waves with phase randomization

Assuming that you are not changing the phase too much, you are basically turning sinusoids to narrow-band noise. If there is a room frequency response node that strongly attenuates a sinusoid, the total power of narrow-band noise centered at the same frequency will not be attenuated as much. It will survive with a dip in its spectrum. Room modes or frequency response peaks are usually not as severe and narrow as frequency response nodes.

An example speaker and room magnitude frequency response (CC-BY-SA 3.0 by Gregory Maxwell)

However, I wonder if your application affords the degradation of audio quality due to the modulation. Also, the room frequency response is not affected so if the test audio is noise, you won't get any sort of improvement.


What is the maximum dynamic range of an N bit signal?

From taking a look at the presentation, it seems that they might be using a slightly different definition of dynamic range than is typical. Usually, as in the Wikipedia article on the topic, dynamic range is defined as:

$$
\text{dynamic range} = \frac{\text{largest possible representable value}}{\text{smallest possible representable value}}
$$

For an $N$-bit (unsigned) signal, this is equal to:

$$
\text{dynamic range} = \frac{2^N-1}{1} = 2^N-1
$$

However, their discussion of dynamic range is interspersed with discussion of quantization noise. Therefore, I posit that they instead define dynamic range as:

$$
\text{dynamic range} = \frac{\text{largest possible representable value}}{\text{largest possible quantization error}}
$$

For a uniformly-quantized quantity such as this, the maximum quantization error is equal to half of one bit. That leads to a dynamic range of:

$$
\text{dynamic range} = \frac{2^N-1}{0.5} = 2(2^N-1)
$$

The extra factor of 2 gives you an approximate increase of 1 bit in "dynamic range" when measured this way.  I assume that's what is meant by an $N$-bit digitized signal providing $N+1$ bits of dynamic range.


What should I do with FFT result array from microphone input?

Musical tone or pitch is different from FFT magnitude peaks.  A musical waveform is rarely a pure sinusoid, integer periodic in aperture.

Please, instead research pitch detection/estimation algorithms (several listed in wikipedia) that do a lot more than simply look at FFT results to determine the tone or pitch of the human singing voice.


A question about Gibbs phenomenon example in Richard Lyons book

I don't have the book, but from what I understand from your question, you take a discrete desired frequency response with $32$ elements, apply an IDFT to get a $32$-tap filter, and then you remove taps to see the effects of truncation in the time domain. Well, this is not exactly what the window design method is about, and it is no surprise that if you don't remove filter taps in the time domain and apply a DFT that you get back your desired $32$-tap frequency response (because, as you've figured out yourself, IDFT and DFT are inverse operations of each other).

However, your actual frequency response does not have $32$ values in the frequency domain but it is continuous in frequency. It is given by the discrete-time Fourier transform (DTFT):

$$H(e^{j\omega})=\sum_{n=0}^{N-1}h[n]e^{-jn\omega}\tag{1}$$

where $N$ is the number of time domain taps (filter length), and $h[n]$ are the filter coefficients. You can compute a sampled version of $(1)$ by using the DFT, and for example, if $N=32$, you can (and should) compute more than $32$ frequency domain samples by using zero-padding. In Matlab/Octave, you can use the commands fft or freqz to compute sampled versions of $(1)$ with an arbitrary number of frequency domain samples.

To help you improve your understanding of the window method, and of DSP in general, I would like to recommend to you a really good (and free) book: Introduction to Signal Processing by S.J. Orfanidis.

Also take a look at this answer, this answer and this answer to related questions.


I will try to tackle questions 1 and 4.

  1) Entropy- It is desired that the entropy of the system be maximized. Maximizing entropy means no symbol is better than the others or we do not know what the next symbol / outcome would be. However, the formula states a negative sign before the summation of the probability logarithms. Thus, it means we are maximizing a negative value!!

No, the values of the logarithms themselves are negative, so the negative sign makes them positive.  All probabilities are a real number from 0 to 1, inclusive.  The log of 1 is zero, and the log of anything less than 1 is negative.  This may seem problematic since the log of 0 is $-\infty$, but we are really trying to maximize the expected value of these logs, so when we multiply by the probability itself the entire value approaches 0, not $\infty$.  Entropy peaks when the probability is $1/2$.

  So,why do we want to maximise entropy since it would mean that we are maximizing the uncertainty of the next symbol whereas we want to be certain about waht the next occurence of the symbol would be.

No, in communicating information we absolutely do NOT want to be certain about what the next symbol will be.  If we are certain, what information is gained by receiving it?  None.  It is only through uncertainty about what the transmitter will send that we can receive any information.

  4) What information does mutual information between two vectors convey?

When there is mutual information between two vectors knowing something about one tells you something about the other.  Mathematically this equates to the following- knowledge of one vector affects the probabilities of the other vector.  If they were independent, this would not be the case.

An example of mutual information is digital walkie-talkies.  One vector is the bit stream that the first walkie-talkie sends.  The second vector is the signal that the second walkie-talkie receives.  The two are obviously related, but due to noise and unknown channel conditions the second walkie-talkie cannot know for certain what the first one sent.  It can make some really good guesses based on the signal, but it cannot know for sure.


PyAudio recording is taking much longer than planned

First off, a small thing:

The sampling rate usually has to be one of the standard sampling rates.
Closest to your 1024*40 would be 44100.

The hardware has to support the sampling rate. Some drivers will resample to the given rate, but that will be somewhat slow as it is in software. Try using a proper sampling rate.

The real bottle neck is the plotting.

Your are calculating the FFT for every block, and plotting it.

The FFT is blazing fast, no worries there.

Plotting on the other hand is horribly slow.  The plotting routine is not optimized for speed, and you are also rescaling and formatting the plot on every go around.

Comment out the lines that do the plotting.  Leave the FFT in.  You'll see a tremendous change in speed.

Edit:
Just noticed.  Besides plot being slower than molasses, you also have a pause in there.  It is only 1 millisecond, but I don't know how precise it really is.  Many timers under Windows use 15 millisecond steps, so that pause could be zero, or 1 millisecond, or 15 milliseconds.

In any case, it is not doing you any good.  Take it out.


A 100x100 pixel image is always a 100x100 pixel image.

You can print it with 100 DPI printer, and it will produce 1x1 inch picture.
A 50 DPI printer will produce 2x2 inch picture because it makes dots that are twice as big.

Same thing when viewing 100x100 image on 15" laptop with 1280x720 screen, it is larger because the pixels are larger than when viewing on a 15" laptop with 1920x1080 screen.

Or same thing when viewing the 100x100 image on 24" 1920x1080 screen, the image is smaller than when viewing the same image on 32" 1920x1080 screen.


What if the Down sampling factor is not an integer?

Resampling (by rational or even irrational ratios) can be done by low-pass filtering in conjunction with high quality interpolation of all the samples needed for the new rate, directly.  No two-step up-sampling following by downsampling is required (although that is one possible implementation for simple ratios).

A Sinc kernel, being a reconstruction formula, is the ideal interpolator for band-limited sampled waveforms.  In practice, a windowed Sinc interpolator is usually as good or better than the FIR or IIR filters used inside other resampling methods.

Also, information may be lost in going to any lower sample rate (even in your original 4:1 resampler) if the original sample data contains higher frequency spectrum that has to be removed by the anti-alias low pass filter (required to prevent aliasing at the lower sample rate).


Identify which of the three signals is closest to a sinusoidal curve

Is it sufficient to identify the «most sinuoidal» of those 3, or would you also want linear projections of those (consistent with a IMU sensor tilted vs the plane of motion)?
A simple solution might be to do a windowed fft and pick the direction where the «crest factor» of the fft magnitude was largest (best explained by a single sinoid).
Edit:
It appears thah your «sine» is somewhat variable. If the above approach fails, then perhaps removing the DC component («highpass filter») followed by removing the high frequency edges and impulses (lowpass filter) for an effective tuned bandpass filter, followed by squaring and summing the samples could be used.


Is Unsharp Mask (USM) Equivalent to Applying Laplacian of Gaussian Filter Directly on the Image?

Unsharp Mask is a sharpening filter.
Intuitively, you apply high pass filter on an image and add the scaled result to the original image.
So the equation you posted is accurate:
$$ o = f + \alpha (h \ast f) $$
Where $ h $ is an High Pass Filter.
If we implement our high pass filter by $ e - g $ where $ e $ is the unit impulse and $ g $ is a low pass filter implemented by a Gaussian Filter you'd get:
$$ o = f + \alpha (h \ast f) = f + \alpha ((e - g) \ast f) = f + \alpha ( f - g \ast f) $$
So basically the result is add to the image the scaled difference between the image and a low pass filtered version of the image.
This is exactly (With some quantization steps) what Photohsop is doing (See Example 001, and Example 002).
Regarding your question, Laplacian of Gaussian (LoG) is an high Pass Filter. So it can replace $ h $ from above.
As you can see, you can't only use it directly but scale the result and add it to the original image.
Difference of Gaussians
As can be seen in the Difference of Gaussians page at Wikipedia, there is a connection between difference of gaussians and LoG. It is explained in Tony Lindeberg - Image Matching Using Generalized Scale Space Interest Points - Appendix A:

Intuitively, we can approximate a Scaled Unit Impulse by a Gaussian Kernel with very small standard deviation.
Now, the difference between Unit Impulse (Or its approximation) and LPF gives us High Pass / Band Pass. It is easy to see in Frequency Domain:

So the logic is: Unit Impulse - Wide Gaussian (Low Pass Filter) ~= Narrow Gaussian - Wide Gaussian = Dog ~= Log. Where DoG and LoG are basically High Pass Filter based on the Gaussian Kernel.
The Equation
From :

The equation says that: $ (1 + \alpha) e - \alpha H $ is the sharpening filter which is correct. Let's rewrite it:
$$ (1 + \alpha) e - \alpha H = e + \alpha e - \alpha H = e + \alpha (e - H) $$
So $ e $ being the Unit Impulse, hence $ e - H $ where $ H $ is a low pass filter (Specifically one could use Gaussian Kernel) gives us an High Pass Filter. Scaling it and adding it wo the neutral item with respect to convolution (The Unit Impulse) gives the sharpening filter. Applying on the image:
$$ f \ast (e + \alpha (e - H)) = f + \alpha (f \ast (e - H)) $$
As written above, Unsharp Mask, which is a sharpening filter, is adding to the image the scaled convolution of the image with an High Pass filter.


As noted in the following link, baud rate is estimated from "the spectrum of the signal after nonlinear processing."  In their case, as noted in the first three sentences of the section entitled "Proposed Spectral Approach for Bit Time Estimation" - p. 916), the nonlinear processing is the first derivative followed by an absolute value operation:

http://www.iasj.net/iasj?func=fulltext&aId=51991

Others have proposed/analyzed estimation using different nonlinear operations, such as squaring of the signal or using higher derivatives.  Obtaining the references for the above paper and doing a thorough search would be a good place to start to learn more.


library for computing differential signatures from 2D planar curves

Recently I have worked on a projective invariant line signature:

  Projective splines and estimators for planar curves
  Thomas Lewiner and Marcos Craizer International Journal of Computer Vision, 2000
  http://thomas.lewiner.org/pdfs/projective_splines_jmiv.pdf

That paper proposes a projective length and a projective curvature estimators for plane curves, when the curves are represented by points together with their tangent directions. The goal of the paper is an invariant spiral construction in the canonical projective frame.

If you like to implement this paper, the paper is pretty self contained. And I posted the details which are kind of hidden in my blog along with some code.

If you don't need projective invariance then I would suggest Fourier Descriptors or Elliptical Fourier Descriptors. There are many MATLAB implementations available for that.

And this work might be quite nice to take a look, if you are concerned about numerics and invariance.


Frequency filtering with a DFT and meaning of removing complex conjugates

Here is a different way to look at it:

The violation of complex conjugate symmetry produces and "error". That error is, that the impulse response becomes complex. By applying the "real()" operator on the inverse FFT, the error is eliminated again.


Understanding index transformation in derivation of Fourier transform for sampling rate reduction

It's not derived, it's just chosen in a smart way such that the relationship between the decimated and the original sequences becomes obvious.
It's just a rearrangement of the terms of the sum. As a simple example, take an infinite sum of numbers $a_r$:
$$S=\sum_{r=-\infty}^{\infty}a_r\tag{1}$$
Under certain conditions that we don't need to bother with now we can rearrange the sum and write it like
$$\begin{align}S&=\ldots+a_{-2M}+a_{-M}+a_0+a_M+a_{2M}+\ldots\\&\ldots+a_{-2M+1}+a_{-M+1}+a_1+a_{M+1}+a_{2M+1}+\ldots\\&\vdots\\&\ldots+a_{-2M+(M-1)}+a_{-M+(M-1)}+a_{0+(M-1)}+a_{M+(M-1)}+a_{2M+(M-1)}+\ldots\tag{2}\end{align}$$
with some arbitrarily chosen integer $M$. So we just start with element $a_0$ and add every $M^{th}$ element, then we move to element $a_1$ and add again every $M^{th}$ element, etc. If we do this $M$ times, we've added all elements, just like in the originial sum $(1)$.
Eq. $(2)$ can be written as
$$S=\sum_{i=0}^{M-1}\sum_{k=-\infty}^{\infty}a_{i+kM}\tag{3}$$
which means that we expressed the index $r$ as $r=i+kM$.