What are the parameters of a calibrated spherical camera? I'm aware if we have a pinhole camera model several parameters describe the specific camera (such as aspect-ratio, focal length, principal point, distortion parameters etc).

There're standard procedures to estimate such parameters that one can follow.

I was wondering however if I replaced the pinhole camera with a spherical camera (such as 360 camera for example) what parameters describe such camera?

I would suppose radius is an example of parameter, are there any other relevant parameters in such case?

Is there any reference you can point out?

...if we have a pinhole camera model several parameters describe the specific camera (such as aspect-ratio, focal length, principal point, distortion parameters etc).

"Distortion parameters" does not sound like a typical pinhole camera model. A pinhole camera does not have a "finite apperture" or way of focusing light other than a tiny little opening which all rays of light (that will go on to produce the image) go through. You can take lens distortions into account by adding them to the pinhole model as additional operations but that then is not the standard pinhole camera model.

The pinhole camera model is a very simple transformation of a point from a global 3D coordinate system to a local 2D coordinate system. It does nothing more than that. It maps the angle at which a point appears within its field of view to a point on its focusing plane. Another useful way to be thinking about lenses is in terms of mapping points in front of them to points behind them.

The pinhole camera model is the simplest of those transformations and requires:

The position of the camera in global coordinates (Where is it?).
The orientation of the camera (Where does it point to?)
The "Focal Length"
The sensor size
The principal point

The focal length in particular is exactly the "thing" that creates this simple map. The shorter the focal length, the wider is the field of view within which points get mapped to a plane behind the pinhole and vice versa.

But notice here, rays in front of the pinhole enter it at some angle and leave exactly behind the pihnole at the same angle.

  I was wondering however if I replaced the pinhole camera with a spherical camera (such as 360 camera for example) what parameters describe such camera?

A pinhole 360 camera is a pinhole camera with 360 degrees field of view. Exactly the same set of parameters but now practically points can enter the pinhole from anywhere around it and depart at the same angle behind it. Notice here that in this view the "focal length" is still valid as a "ratio" of lengths.

But in reality, there are no 360 degree "pinholes". Instead, the wide field of view is acquired in some other way. That might be a mirror for example or multiple cameras with various lenses.

In this case, it is even more useful to be thinking in terms of rays entering the lens and rays departing from the  lens (or optical system in general). The general model in this case is established around a transfer matrix that maps "entry points" to "exit points" and basically can model any sort of configuration you like.

In terms of references, Geyer and Daniilidis:A unifying theory for central panoramic
systems and practical implications is probably still relevant and since then, Daniilidis has also put together a volume of relevant work, the contents of which you might find useful, even if it is just for more relevant terminology to get you closer to what you are dealing with.

Hope this helps.


Acoustic features distinguishing the male and female voice? I would like to use MFCC analysis to distinguish between a male voice and a female voice.

Are there any particular quantifiable characteristics consistent with, and exclusive to, either gender? (E.g., energy, fundamental frequency, excitation, the presence of certain waveforms, etc.)

Main difference is the frequency of the fundamental, which is about an octave higher for women with the split point around 160 Hz or so. A fundamental lower than 160 Hz is most likely a male and a fundamental higher is most likely a female. A good overview over a number of studies on the topic can be found at http://www2.ling.su.se/staff/hartmut/f0_m%26f.pdf


gaussian kernel and bionomial coefficients The following image is from robert collins' ppt from http://www.cse.psu.edu/~rcollins/CSE486/

How to calculate sigma ?

It shows [1 4 6 4 1] is corresponding to sigma=1.

I learned that usually 3*sigma is best to represent. So kernel size is 7*7 is suited. But its size is 5*5.

I don't know why [1 4 6 4 1] is matched to sigma=1.

The well-known twentieth-century American columnist Dorothy Parker once reviewed a book
saying

  "This is not a book to be put down lightly; it should be thrown with great force."

Be that as it may, Professor Collins is suggesting that a bar plot of the binomial
coefficients on the $n$-th row of Pacal's triangle, after suitable normalization,
will look from far away like a Gaussian distribution with mean $n/2$ and standard
deviation $\sqrt{n}/2$.  He then further obfuscates this by using the fact that for
even values of $n$, we can write the standard deviation as $\sqrt{k/2}$
where $k = n/2$
while for odd $n$, it is simply $\sqrt{n}/2$, inviting the reader to arrive at
the mistaken conclusion that there is a fundamental difference between the odd
and even rows.

To get a better understanding of where the $n/2$ and $\sqrt{n}/2$
comes from, look for reading material on the DeMoivre-LaPlace approximation
which is a special case of the Central Limit Theorem.


Convolution property Is the convolution of two equal signals the same signal?

I have this system:

and i have to find $h_1[n]$ (which i've been struggling for a while), given

$$ h_2[n] = u[n] - u[n-2] $$

and $h[n] = h_1[n] \star h_2[n] \star h_2[n]$ given by:

I can find a solution if $h_2[n] \star h_2[n] = h_2[n]$, but I don't know if that holds.

First, we define $h_3[n] = h_2[n] \star h_2[n] = [1, 2, 1]$. From this result, we know that $h_1[n]$ must have 5 elements (so that $h[n]$ ends up with 7 elements).

Let's define $h_1[n] = [g_1, g_2, g_3, g_4, g_5]$. We can find these as follows, using the definition of discrete convolution.

First, we know that $g_1 h_3[0] = h[0] = 1$, so $g_1 = 1$.

Then, the next value of the convolution is $h[1] =  g_2 h_3[0] + g_1 h_3[1] = 5$, from which we deduce that $g_2 = 3$.

Then, $h[2] = g_3 h_3[0] + g_2 h_3[1] + g_1 h_3[2] = 9$, or $g_3 = 2$.

You can stop here, since (by symmetry, in this particular example) $g_4 = g_2$ and $g_5 = g_1$.

In general, the problem can be seen as solving a system of linear equations. In this example, we need to solve the system:

$$\begin{bmatrix}
h_3[0] & 0 & 0 & 0 & 0 \\
h_3[1] & h_3[0] & 0 & 0 & 0 \\
h_3[2] & h_3[1] & h_3[0] & 0 & 0 \\
0 & h_3[2] & h_3[1] & h_3[0] & 0 \\
0 & 0 & h_3[2] & h_3[1] & h_3[0] \\
0 & 0 & 0 & h_3[2] & h_3[1] \\
0 & 0 & 0 & 0 & h_3[2]
\end{bmatrix}
\begin{bmatrix}
g_1 \\ g_2 \\ g_3 \\ g_4 \\ g_5
\end{bmatrix}
=
\begin{bmatrix}
h[0] \\ h[1] \\ h[2] \\ h[3] \\ h[4]  \\ h[5] \\ h[6]
\end{bmatrix}$$

Thanks to @Hilmar for pointing out a boneheaded mistake in my original answer!


1st order filter output Two 1st order filters :

--[ b0, b1, a0, a1]
1 [ -1.40374673978609e+008 1.40374675978609e+008 1 1.000000235752411 ]
2 [ 2.22044577169459e-016 2.22044632680604e-016 1 -1 ]

gives these responses in Octave:

When using those same coefficients (in parallel order or even just one set of them) in another software I get different results. Example:

Magnitude response isn't a problem (unless not 'flat' in range 20Hz-20kHz) but the phase response is 180° even when implement only one of those filters. Shouldn't phase be 90° as it is in Octave plot.

Q1: is it just because of the transfer function (1st order (Octave) vs biquad (other software)) or is there something else beind this difference?
Q2: is it possible to plot these two filters as parallel system in Octave/MatLab

Matlab results.

The filters coefficients have a big dynamic variation. So you have those results. For second filter, coefficients b are about zeros.

The first filter

b = [-1.40374673978609e+008 1.40374675978609e+008];
a = [1 1.000000235752411];
freqz(b, a)

The second filter

b = [2.22044577169459e-016 2.22044632680604e-016];
a = [1 -1];
freqz(b, a)


Resize line art image I want to print some line drawings on A3, but they aren't very large - around 600x800px.

I tried obvious functions in Gimp, like resize then enhancing with wavelet sharpener or unsharp mask, with selective gausian blur before it, but print results are no good.

They are not easy for vectorization either, but perhaps tweaking some tracing tool can give me better results then resizing bitmap. Before I start exploring this option, I thought to ask for advice:

Can anyone share a recipe, either by using programming algorithms (numpy/scipy or matlab) or by using tools available in Gimp or Photoshop, to accomplish quality resize of line drawing?

Take a look at pixel art scaling algorithms, particularly the hqx algorithms.  They're typically limited to a integer scaling factor (i.e. 2x, 3x, 4x) so they aren't completely general purpose but they should do quite well on line drawings.


normalize fft produced by fftw_plan_dft_r2c_1d fftw_plan_dft_r2c_1d takes an array for real numbers and produces an array of complex numbers 1/2 the size of input array + 1.

per fftw's website:

  As a result of this symmetry, half of the output Y is redundant (being the complex conjugate of the other half), and so the 1d r2c transforms only output elements 0...n/2 of Y (n/2+1 complex numbers), where the division by 2 is rounded down.

When I normalize do I normalize by the size of the fft or the resulting fft_size / 2

    // normalize the fft
    for (unsigned int x = 0; x <(fftSize/2); x++)
        array[x] = array[x] / complex<double>((fftSize/2), 0);   <<<<<

or

    // normalize the fft
    for (unsigned int x = 0; x <(fftSize/2); x++)
        array[x] = array[x] / complex<double>(fftSize, 0);       <<<<<

Reading the manual, even though only half the data is returned, the full length FFT is calculated:

So you should normalize by whatever the full length is normalized to: fft_size.


Audio Signal Separation - Identifying Interest points I asked this question here: Audio Separation of .wav signal but it wasn't clear, so, here is my second attempt:

First off, assume that I have a .wav file containing a sentence as follows:

"My name is Michael" I would like to extract, from this, the following:

"My" -> Phoneme (1)

"Name" -> Phoneme (2)

"Is" -> Phoneme (3)

"Michael" -> Phoneme (4)

This means that I have taken my 1D signal, and split it into a 2D signal (vector) that contains these particular words/phonemes which I can then analyse and identify. I would therefore, like to compute this in the time domain and not the frequency domain. Just to clarify again:

I take in a 1D signal containing a sentence, split this sentence up into different parts which contain this data: vect[0], vect[1], .... vect[4] Let's say in matlab I did the following command wavwrite(vect[0], ....) then it would output the word "My" and putting all the blocks together would give me the full sentence back.

Here is my "real-world problem" instead of Phonemes, I have bat calls, the length of each bat call is unknown at this stage. But here is a typical sample of a bat call: Here and for each of these bat calls, these need to be separated from the inputted signal and stored inside a vector (Just like the example above), this, then allows me to identify each of the bats and perform analysis on them.

This, like, the sample would give me a 2D vector containing each of the bat calls: "Bat1", "Bat2" ..... "Bat[n]" it is unknown the amount of time the bats have been recorded for, or, what the length of each of the bat call therefore is.

What I have done so far:

I have obtained the bat signal, processed it and I am given the following (Which is plotted):

I have also Emphasised the signal using the following formula:

rawSignal[i] = rawSignal[i] - (0.95 * rawSignal[i-1]);

And then I have Compressed the signal using the following:

    float param = 1.0;
for(unsigned i=0; (i < rawSignal.size()); i++)
{
    int sign = getSignOf(rawSignal[i]);
    float norm = normalise_abs_value(rawSignal[i]);
    norm = 1.0 - pow(1.0 - norm, param);
    rawSignal[i] = denormalize_value(norm, sign);
}

Which then gives me an output of the following:

I'm unclear to where I should go from here in identifying single elements ("calls") from this signal. Since, if I use zero-crossing and/or calculating the total energy of the signal and thus using a threshold then it will just remove the noise and I'm left with a compressed version of the signal.

Speaking to someone, they suggested that I should try and use the Cochleagram domain however, I'm not familiar with this and there is very little research on this available.

If anyone has any suggestions, or the algorithms that I could use then please suggest them.

(a follow-up to my suggestion on the previous question), you can use the spectrogram and ICA to help:

A similar shorter sound file:

import wave, struct, numpy as np, matplotlib.mlab as mlab, pylab as pl
def wavToArr(wavefile):
    w = wave.open(wavefile,"rb")
    p = w.getparams()
    s = w.readframes(p[3])
    w.close()
    sd = np.fromstring(s, np.int16)
    return sd,p

def wavToSpec(wavefile,log=False,norm=False):
    wavArr,wavParams = wavToArr(wavefile)
    print wavParams
    return  mlab.specgram(wavArr, NFFT=256,Fs=wavParams[2],detrend=detrend_mean,window=window_hanning,noverlap=128,sides='onesided',scale_by_freq=True)

wavArr,wavParams = wavToArr("bat_speech.wav")
hf = pl.figure(); ax=hf.add_subplot(1,1,1)
ax.plot(wavArr)

Now take a look at the spectrogram:

Pxx, freqs, bins = wavToSpec("bat_speech.wav")
Pxx += 0.0001
freqs += (len(wavArr) / wavParams[2]) / 2.
hf=pl.figure(figsize=(12,12));
ax = hf.add_subplot(2,1,1);
#plot spectrogram as decibals
hm = ax.imshow(10*np.log10(Pxx),interpolation='nearest',origin='lower',aspect='auto')
hf.colorbar(hm)
ylcnt = len(ax.get_yticklabels())
ycnt = len(freqs)
ylstep = int(ycnt / ylcnt)
ax.set_yticklabels([ int(freqs[f]) for f in xrange(0,ycnt,ylstep) ])

We can clip this at 8000Hz or so it looks like, or don't bother cleaning it up.

Now you have frequencies which can proxy for multiple sources in BSS, so you can play with PCA, ICA, normalization, etc.  For example, see if you have some components you can isolate:

from sklearn.decomposition import PCA, FastICA
ncomps = 7
# reduce dimensionality with PCA
pca = PCA(n_components=ncomps)
y = Pxx.copy().T
pc = pca.fit(y).transform(y)
# run ICA
ica = FastICA(n_components=ncomps,random_state=42)
z = ica.fit(pc).transform(pc).T
hf = pl.figure()
for p in xrange(ncomps):
    ax = hf.add_subplot(ncomps,1,p+1)
    ax.plot(z[p])

or see if the spectrogram is enough to let you do your segmentation:

hf = pl.figure()
ax = hf.add_subplot(1,1,1)
ax.plot(np.sum(Pxx,axis=0))

EDIT:  Just realized you said you didn't want to use the frequency domain, but it may help you isolate your phonemes which you can extract.  Anyway, running on your SeroWeb.wav I get this:


Relation between Covariance matrix & Energy of a random signal Let's say I have the below random signal:
$ Y[n] = [y(n), y(n-1), y(n-2), \ldots, y(1)] $

I have two random variables now: The first one $X_1 $ which express the maximum eigenvalue of the covariance matrix of $Y$. The second one $X_2$ which express the energy of the random signal.

Now my question is: Are the two random variables independent or dependent, when whether signal samples $y(n),y(n-1)\ldots$ are IID or correlated with each other?

By intuition the two random variables should be correlated!? Isn't it right!?

I'm not totally sure I understand your question: could you give a concrete example? For instance, what's the randomness? Anyway, here are some relations between the quantities you're interested in.

The energy of the signal is $$\mid \mid Y \mid \mid_2^2 = \sum_{i=0}^n y_i^2$$

The Auto-correllation matrix $R$ is formed as $YY^T$ (outer product between the signal and it's transpose).

Since this matrix is symmetric, it's matrix 2-norm is the spectral radius: $$||YY^T||_2^2 = \mathrm{max}|\lambda_j|$$

Another way to write the 2-norm would be as $$||YY^T||_F^2 = Tr(YY^T)$$ It's fairly stragtforward to prove that $$||* ||_F^2 \geq ||*||_2^2$$ for any matrix, and they co-incide when $YY^T$ has rank 1. It's also true that $$|| * ||_F^2 \leq r|| *|| _2^2 $$


Filtering of complex exponential functions In my problem the signals $x(t), h(t), w(t), y(t)$ and $y_w(t)$  are defined as follows:

\begin{align}
x(t)&=\sum_{i=1}^M a_i\delta(t-\tau_i) \\
h(t)&=e^t\exp\left(-e^t\right) \\
w(t)&=\begin{cases} 1, & \text{when}\quad |t|\le\frac{\tau_0}{2}\\
0, \quad &\text{elsewhere} \end{cases}\\
y(t)&=x(t)*h(t) \\
y_w(t)&=y(t)w(t)
\end{align}

Can any one recommend me a way how to get the $x(t)$ for only measuring $y_w(t)$ and knowing the $\tau_0$ by using the FFT.

Hints:
\begin{align}
X(\omega)&=\sum_{i=1}^M a_ie^{j\tau_i\omega}\\
H(\omega)&=\Gamma\left(1-j\omega\right)\\
W(\omega)&=\tau_0\textrm{sinc}\left(\frac{\tau_0\omega}{2\pi}\right)
\end{align}

HINT:

If the given $X(\omega)$ is the Fourier transform of $x(t)$, what does $x(t)$ look like? It will turn out that removing exponential terms with certain values $\tau_i$ is equivalent to suppressing certain time intervals in the time domain.


How does spatially white noise differ from diffuse noise In noise reduction, having read that directivity index is used to reduce diffuse noise and white noise gain reduces spatially white noise, I would like to know the difference between both the noise types. Is it that diffuse noise is identically distributed but not white? Hope some one could help :)

In a diffuse noise field the noise energy flows in all directions with equal probability. This does not mean that the noise at one point is uncorrelated with the noise measured at another point. On the other hand, for spatially white noise, also called incoherent noise, the noise signals at two different points in space are always uncorrelated.

An example for a diffuse noise field is noise in a car or in an office. An example for spatially white noise would be amplifier noise in a microphone array, where the noise in each channel is independent of the noise in other channels.


Damped sinusoidal form FFT of signal I'm doing an assignment for the course Signal analysis where I have to analyse a signal. I've tried quite some things now but it's still bothering me that the FFT is looking weird, and not looks like the 'normal look' FFT's we learned in class.

FFT (absolute values of complex values):

FFT (no absolute values):

The FFT seen in the image above is zoomed in on the frequency range 0-30Hz. The rest of the frequency range does not show a lot of (high) peaks, which probably are caused by noise.

The signal is created during a method of welding, using an oscilloscope with a sampling frequency of 1000Hz. I've filtered the signal to remove noise, and after that the signal is converted to the frequency spectrum using the fft function of MATLAB.

Signal before and after filtering:

My general question is, can the shown FFT be valid or did I make a mistake? I estimated the ground frequency to be around 5.5Hz, can I say this when I take one period of the big sinusoidal wave? I also noticed there are about 64 little sinusoidal waves inside one (ground??) period, is this an high harmonic wave form?.

If my theory is right, what causes the fft to be a damped sinusoidal form?

The code I use is basically the following. I leave the part of the noise filtering out because I don't think it's necessary for this question. The dataset is an matrix of 40100 rows.

fs = 1000;
cleanSignaal = data(:,4);
fftSignal = fft(cleanSignaal)/lenght(cleanSignaal);
f = fs/(2*length(fftSignal)):fs/length(fftSignal):fs;
plot(f,abs(fftSignal));
xlim([0 fs(m)/2]);
title('Fast Fourier Transform')
xlabel('Frequentie (Hz)')
ylabel('Magnitude')

What you have is a very sharp rectangular boxcar shaped function.  It is approximately centered at the mid point of your data set.  Your magnitude spectrum looks correct, in that it is a really low fundamentally frequency with a ton of harmonics.  Your FFT without ‘absolute values’ is presumably plotting just either the real or imaginary values.  The oscillation in that plot occurs because of where that pulse is centered in the data set.  It’s unlikely that plot has any useful information.  I have no idea what you are trying to find from this endevour, but if you were interested in the period of the pulse you could zero of the end of the data set, providing more granular frequency data that you could use to measure the fundamental.  That or you could just measure it directly, which is probably a lot easier.


High-pass filter in Python (Scipy) This code is taken from a pitch detection algorithm. It is called before pitch detection to remove low-frequency noises.

def highpass_filter(y, sr):
  filter_stop_freq = 70  # Hz
  filter_pass_freq = 100  # Hz
  filter_order = 1001

  # High-pass filter
  nyquist_rate = sr / 2.
  desired = (0, 0, 1, 1)
  bands = (0, filter_stop_freq, filter_pass_freq, nyquist_rate)
  filter_coefs = signal.firls(filter_order, bands, desired, nyq=nyquist_rate)

  # Apply high-pass filter
  filtered_audio = signal.filtfilt(filter_coefs, [1], y)
  return filtered_audio

I have trouble understanding what the parameters filter_stop_freq, filter_pass_freq and filter_order mean. I have read this but it did not help me.

Filters generally can't just have a single cut-off frequency (otherwise their complexity will be huge).  So a transition region is usually needed.

In this case, the filter_stop_freq is that frequency below which the filter MUST act like a stop filter and filter_pass_freq is that frequency above which the filter MUST act like a pass filter.

The frequencies between filter_stop_freq and filter_pass_freq are the transition region or band.

Because this is designing an FIR filter, filter_order is just the length (number of taps) of the filter. The higher the order, usually, the better the fit (and the smaller the transition region can be).


How to use a FIR on a signal with a drifting central frequency? I have been working on this subject for a while and I couldn't get a proper answer.

I am trying to estimate the phase of the input signal :

$x(t)=Ae^{j(at^2+bt+c)}+n(t)+i(t)$

Where $n(t)$ is a AWGN, $i(t)$ is an interfiring signal and $a,b,c$ are constant values.
As you can see, the central frequency of the useful signal changes over time ($f_c(t)=2at+b$).

First, I have to filter out the interfering signal so I implemented a low-pass FIR filter H and I compute the output of the filter as follow $y[n]=\sum{h[k].x[n-k]}$. However because of the moving frequency my signal will not stay in the filter passband very long.

One solution I had was to estimate the central frequency of the signal every 1 second and shift the input signal by this estimate (which means the I add to the FIFO of my FIR the samples $x[k]e^{-2\pi f_ck/F_e}$) but when I estimate the new central frequency I will have samples with different central frequencies in my FIFO and the filter doesn't like that (I can see phase jump at the time of frequency change). So can anybody tell me how I am supposed to tackle this issue ?

Thank you !

EDIT : Thank you very much for your answers, I realized while reading your answers that I simplified my problem too much. The actual evolution of the phase of my interest signal is much more complicated than a second order polynomial. My signal is more like
$x(t)=Ae^{jd(t)}+n(t)+i(t)$ with $d(t)$ a very complex function taking into account the distance and the speed of the object emitting the signal (which is a pure carrier). Also, the $i(t)$ is a function representing other object emitting other pure carrier at different speed and position (so my signal is more like $x(t)=\sum{A_ke^{jd_k(t)}}+n(t)$). I want to get the phase of each signal separately so first of all I did a small accumulation of signal to compute the spectrum of the signal and get the approximate frequency of each signal. Then I apply several FIR filters, one for each signal hoping I would get each $A_ke^{jd_k(t)}$ separately. The problem I have is that the frequency of each signal moves a lot and I have to "adapt" my FIR filter (and I don't know how to do that...)
I talked about a second order phase polynomial because it is a good model of the function $d(t)$ on a short term observation.

Thank you once again for the time and effort you put in your answers.

Without understanding more detailed results on what i(t) is, one solution is to estimate either the phase or frequency first, and then pass that result through a filter.

An approach to do the frequency estimation digitally directly on the I (real) and Q (imaginary) samples of x(t) is the cross product frequency discriminator. After hard limiting x(t) to remove the AM components (which would result in samples on the unit circle in the complex plane), the frequency can be estimated using:

$$F \propto  I_1Q_2-I_2Q_1$$

Where $I_1$,$Q_1$ and $I_2$,$Q_2$ represent successive samples. This is because the cross product is proportional to the phase between two vectors, so in this case the phase is the phase versus time, and therefore frequency. This is shown in the math and diagram below.

Frequency is the change in phase versus time ($d\phi/dt$). To determine the phase between two hard limited successive samples in time $V_1$ and $V_2$ given by:

$V_1 = e^{j\phi_1} = I_1+jQ_1$

$V_2 = e^{j\phi_2} = I_2+jQ_2$

For small angles, the phase $\phi_2-\phi_1$ can be reasonably approximated by the imaginary component of the conjugate product $(V^*_1) (V_2)$:

$V^*_1V_2 = e^{-j\phi_1}e^{j\phi_2} =e^{j(\phi_2-\phi_1)}$

For small angles $cos(\phi) \approx \phi$,

So therefore for small angles

$\phi_2-\phi_1 \approx Im(e^{-j\phi_1}e^{j\phi_2})$

$= Im((I_1-jQ_1)(I_2+jQ_2))$

$=I_1Q_2-I_2Q_1$

The above is extremely simple with low processing requirements useful for hardware implementation but depends on the small angle criteria for accurate angle estimation. This will work well as long as the sampling rate is high enough that the phase between samples doesn't deviate significantly. Alternatively if further processing is feasible, an ATAN2 operation can be instead be done on the two successive I,Q samples for full 4 quadrant accurate phase estimation (and therefore frequency given the phase difference between two successive samples over a much wider frequency range).

A simple delay and multiply frequency discriminator could also be used, followed by a  filtering of the resulting frequency versus time waveform. I have further detailed an example of using a simple delay and multiply frequency discriminator in this recent answer: Demodulation of FSK signal

Both approaches require hard limiting your waveform first to remove any amplitude variation of your signal. A differentiator as done through the differencing in these approaches is a high pass filter, the hard limiter removes all AM noise, and the post detection filtering will remove the noise induced "jumps" that you see in the result.


Physical significance of resolution in image processsing What is the physical interpretation/significance of resolution?

Technically, it is the number of pixels used to denote an image but what is it that this resolution tries to convey to us?

So the question really is, "what is a pixel?"

A pixel when you access the data/sample of the image is just a number. The capturing device has generated this data represents the visible area on a particular surface and the color of the pixel represents the light that is falling on that surface.

Hence, when resolution of the capturing devices is high, you are able to capture equivalent of finder details from the incoming source of light and its' variations. when the same image is brought down to a resolution of say 1/4 which means the resultant new image represents the average color/intensity of that 4x4 window around that point. Hence, the lower resolution more "pixelating effect" do you see.

There is one more aspect to see here. The details in image also represents the size of the object subject to the scene geometry. So for example, if you capture a shot in camera which has one flower and a leaf, we identify that if the leaf is roughly taking twice as more pixels than flower, then leaf in real life is twice as bigger (from that perspective) than flower. However, when zoomed out, and pointing to the macro scene, the same number of pixels could be occupied by a full tree. So we can find the relative dimension of the object based on their area in the image because image is nothing but a projection of that scene on the canvas. It is possible to know, however, the exact dimension of the object as well if we can apply some known distances on the canvas or if we know the distance of the object from the camera film.

So in essence, pixel represent the light intensity or the color of the object at a fixed grid on the scene - the finer the resolution - more pixel reprsent the same object in finer details and more information. And, distance between the pixels represents the relative distance between two points on the scene.


How can a changing signal be Fourier transformed? Fourier transformation is used to split a periodic signal into frequencies, and calculate the phase and shift for each frequency. But what if we would like to record an audio, which is changing, for example, because someone starts to play with a new instrument?

I know that if the recording is 44100 Hz, then we record 44100 samples in each second.  But how is it Fourier transformed? The only way I can think of is that one second is divided into 25 segments, and we do the transformation to get the phase and shift for each frequency in each segment.

You are correct that it is often broken up into chunks. There are multiple ways to combine those chunks (Google: fft overlap save).

It is often important to window the data (again Google: hamming window fft). One critical reason to do this that gets to the heart of your question: the fft assumes the chunk you fed it continues infinitely forward and backward in time. Obviously an audio file does not have this property.

For your example you might want to chunk your audio file into groups of powers of 2 like, 1024,2048, or 4096.


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

A full-width rectangular window (a.k.a no window) would best emphasize spectral content near the edges of an FFT aperture, at the cost of artifacts caused by convolution of the spectum with a Sinc function, also called "spectral leakage".  If you can ignore or compensate for these spectral artifacts then there's little reason to window.

But you might be asking the wrong question.  If you want the fastest possible response to some specific frequency band, you should be using a minimum phase filter (or a bank of them, as needed) or a close approximation to such, not a faster block DFT.  Then you can update with each new sample instead of waiting to process a block, and tune the time-frequency response of the filter(s) for each specific band of interest.


Conjugation in Fourier Transform I have a very simple question.

In Oppenheim book, it says that:

  If CT Fourier transform of $x(t)$ is $X(j\omega)$ then, CT Fourier
  transform of $x^*(t)$ is $X^*(-j\omega)$.

What I can't understand is what is $X^*(-j\omega)$? Or what is $X^*(j\omega)$? Can you tell me what I should understand when I see $X^*(-j\omega)$ or $X^*(j\omega)$?

For example $x(t) = (a + bj)t$ and $X(j\omega) = (c + dj)\omega$ what should be $X^*(-j\omega)$?

Thanks

$X^*(j\omega)$ is the complex conjugate of $X(j\omega)$. So if

$$X(j\omega)=X_R(\omega)+jX_I(\omega)$$

then

$$X^*(j\omega)=X_R(\omega)-jX_I(\omega)$$

and

$$X^*(-j\omega)=X_R(-\omega)-jX_I(-\omega)$$

To gain some more understanding, consider the definition of the Fourier transform:

$$X(j\omega)=\int_{-\infty}^{\infty}x(t)e^{-j\omega t}dt$$

So we have

$$X(-j\omega)=\int_{-\infty}^{\infty}x(t)e^{j\omega t}dt$$

and

$$X^*(-j\omega)=\int_{-\infty}^{\infty}x^*(t)e^{-j\omega t}dt=\mathcal{F}\{x^*(t)\}$$

as claimed in the quote of your question.


Transform matrix creation in compressed sensing In Compressed sensing ECG reconstruction I would like to use biorthogonal wavelet basis as transform/sparsifying domain. I found out that wfilters command in MATLAB will give 4 type of filters, but which filter I would use to make the transform matrix?

Either decomposition lowpass/highpass or reconstruction low pass/high pass or is there any other better ways to make the transform matrix??

I recommend not forming the transform matrix directly. That will cost more memory and take longer to compute the transform. Instead, use the idwt/dwt or waverec/wavedec functions in place of multiplying by A/A.T.
Many reconstruction solvers and convex optimisation toolboxes support supplying function handles for this instead of a matrix. If you are implementing your own solver it should also be easy enough to integrate these functions instead of the matrix-vector products.


FFT frequency resolution I'm having some problems understanding the FFT. Is the frequency resolution in the spectrum calculated as

$\frac{\textrm{sampling rate}}{\textrm{number of FFT points}}$ or $\frac{\textrm{sampling rate}}{\textrm{0.5 * number of FFT points}}$?

Asking this because the spectrum is symmetric for real-valued inputs. So, say I have $f_s = 1000$ Hz and $N$ = 1024, where $N$ is the number of FFT points. Now, is the frequency resolution $\frac{1000 \textrm{ Hz}}{1024} = 0.9766$ Hz or $\frac{1000 \textrm{ Hz}}{0.5 * 1024} = 1.9531$ Hz?

Suppose you have a signal $x[n]$, with $n \in {0, 1, ... N - 1} $. The same-size DFT is defined by:

$$
X[k] = \sum_{n = 0}^{N-1} x[n] e^{-j\frac{\Large{2 \ \pi n k}}{N}}
$$

The frequency resolution is going to be how many Hz each DFT bin represents. This is, as you have noted, given by $\frac{f_s}{N}$.

If on the other hand you had zero-padded your signal, such that $N_{zp}$ is greater than  $N$, then a more apt term of frequency granularity is given by $\frac{f_s}{N_{zp}}$

  Asking this because the spectrum is symmetric for real-valued inputs.

That is irrelevant. The frequency resolutions/granularities are given by the above.

  o, say I have fs=1000 Hz and N = 1024, where N is the number of FFT
  points. Now, is the frequency resolution 1000 Hz1024=0.9766 Hz or 1000
  Hz0.5∗1024=1.9531 Hz?

If your sampling frequency $f_s = 1000$ Hz, and you are taking an $N=1024$ (same-sized) FFT, then your frequency resolution is $\frac{1000}{1024}$, which is equal to 0.9766 Hz/bin. If your $N_{zp} = 1024$ (FFT length after zero-padding), then your frequency granularity is 0.9766 Hz/bin.


How does one know how many times to downsample in Harmonic Product Spectrum? I have been working on a Harmonic Product Spectrum Algorithm. Now, all the literature I've read about the subject tells me to downsample an N amount of times. How does one determine what this N value should be? Here is my implementation of Harmonic Product Spectrum so far. Feel free to tell me if I have gone wrong somewhere.

private int HarmonicProductSpectrum(Complex[] fftData, int n){
    Complex[][] data = new Complex[n][fftData.length/n];
    for(int i = 0; i<n; i++){
        for(int j = 0; j<data[0].length; j++){
            data[i][j] = fftData[j*(i+1)];
        }
    }
    Complex[] result = new Complex[fftData.length/n];//Combines the arrays
    for(int i = 0; i<result.length; i++){
        Complex tmp = new Complex(1,0);
        for(int j = 0; j<n; j++){//multiplies arrays together
            tmp = tmp.times(data[j][i]);
        }
        result[i] = tmp;
    }
    //Calculates Maximum Magnitude of the array
    double max = Double.MIN_VALUE;
    int index = -1;
    for(int i = 0; i<result.length; i++){
        Complex c = result[i];
        double tmp = c.getMagnitude();
        if(tmp>max){
            max = tmp;;
            index = i;
        }
    }
    return index*getFFTBinSize(fftData.length);
}

Take a look at this paper on page 4. It has a great graphical representation of the algorithm. I'll try to write up some pseudocode here that you can adapt to your language (looks like Java to me) as necessary.

/* Harmonic Product spectrum PSEUDOCODE
 * @param fftData Discrete Fourier transform of data
 * @param N Number of times we downsample the spectrum to get HPS
 */
int HPS( Complex[] fftData, int N )
{
    // Find magnitude of the FFT
    Real fullSpectrum[] = absOfComplex(fftData);

    // Keep only the positive frequencies (DC to Nyquist)
    Real spectrum[] = discardNegativeFrequencies(fullSpectrum);

    // Make a new array to store HPS
    Real hps[] = copyOf(spectrum);

    // Perfrom HPS:
    // Go through each downsampling factor
    for (int downsamplingFactor = 1; downsamplingFactor <= N; downsamplingFactor++)
    {
        // Go through samples of the downsampled signal and compute HPS at this iteration
        for(int idx = 0; idx < spectrum.length()/downsamplingFactor; idx++)
        {
            hps[idx] *= spectrum[idx * downsamplingFactor];
        }
    }

    return findIndexOfMax(hps);
}

This isn't bug-proof but it should get you started.


Simulation at Baseband or Passband A naive question I guess. Would it make a difference in the results when a system is simulated at Baseband or at Passband. If I have a plot of a BER vs SNR curve that has been plotted at Baseband, will the same change if I perform the simulation at Passband with a carrier frequency of 2.4GHz? AFAIK, it should not make a difference, but I would like to know if what I am thinking is correct.

Another question is that if I have two signals from two transmitters and add these signals at the receiver, does it make a difference if these two signals are passband or baseband. Do the results vary? Since Passband simulation consumes a lot of time, would it be correct if I simulate the behavior at Baseband and assume that the result would be almost the same at the Passband?

The only thing I could see, with the information you gave, is that the way you get from passband to baseband and vice-versa could affect your performances. It should not, but it could.

Let's say you have to down sample to get to baseband, and up sample to do the contrary and you use filters to do so. If your filters have variations of 1dB in their passband and a cutoff frequency that cuts some of your baseband signal, then your performances could be deteriorated. In an ideal situation it should not, but some considerations (restriction on the number of coefficients of your filters, for example) may imply some performance degradation.


Training Fractionally Spaced Channel Estimator and Equalizer If you are attempting least squares channel estimation with a fractionally spaced channel estimator, do you want the training sequence to also be fractionally spaced or symbol spaced? It looks like you could do either.

The training sequence should be at the same spacing as the equalizer when considering its sampling at the input to the equalizer. Adaptive algorithms converge to the least square solution based on the error between the received sequence and the transmitted sequence (known when a training sequence is used). Further, the equalizer can only determine a solution where channel energy is present (for example, if you transmitted a single tone, you would only find the single phase/magnitude weight for that frequency)

If you use only one sample per symbol to derive the channel, the solution will suffer from "band edge aliasing" as I depict in the figures below, and the channel solution will suffer from the folded energy at the band edges.

For further details see this post where I have derived the operation of the LMS Equalizer and the solution using the Wiener-Hopf equations (for when you have the luxury of block post-processing, but also helps illuminate the operation since this is the same solution that the adaptive algorithms would converge to).

Compensating Loudspeaker frequency response in an audio signal


Estimate Taylor series coefficients from samples of a function Say I have measurements of a function $y = y(x)$, sampled at $x_i$ with some noise, that could be approximated by a Taylor series expansion. Is there an accepted way of estimating the coefficients for that expansion from my measurements?

I could fit the data to a polynomial, but that's not quite right, because for a Taylor series the approximation should get better the closer you are to a central point, say x = 0. Just fitting a polynomial treats every point equally.

I could also estimate the various orders of derivatives at my point of expansion, but then I need to make decisions about what differentiating filters to use and how many filter coefficients for each. Would the filters for the different derivatives need to fit together somehow?

So does anybody know of established methods for this? Explanations or references to papers would be appreciated.

CLARIFICATION

In response to the comment below, my sampling is a rectangular window from an infinite function, that isn't necessarily band-limited but does not have strong high-frequency components. To be more specific, I'm measuring the variance of an estimator (measuring displacement in a medical ultrasound signal) as a function of a parameter of the estimator (the level of deformation or strain of the underlying tissue). I have a theoretical Taylor series for the variance as a function of deformation, and would like to compare it to what I get from simulation.

A similar toy example might be: say you have a function like ln(x), sampled at intervals in x with some noise added. You don't know what function it really is and you want to estimate its Taylor series around x=5. So the function is smooth and slowly varying for a region around the point you're interested in (say 2 < x < 8), but is not necessarily nice outside of the region.

The answers have been helpful, and some kind of least-squares polynomial fit is probably the route to take. What would make an estimated Taylor series different from a normal polynomial fit, though, is that you should be able to shave off higher-order terms, and have the polynomial still approximate the original function, just within a smaller range about your initial point.

So maybe the approach would be to do a linear polynomial fit using only data close to the initial point, followed by a quadratic fit with a little more data, cubic using a little more than that, etc.

Instead of exact polynomial fitting, you could use a least-squares fit, which will find the polynomial of specified order that minimizes the total squared error between the fit and the measured $(x_i,y_i)$ pairs. This can help to mitigate the effects of noise on the fit.

Given measurements $y_i$ of a function $y = f(x)$ at domain values $x_i$ ($i = 0, 1, \ldots , N$), choose a polynomial order $M \le N$ (if $M = N$, then you're down to exact polynomial fitting, as the $N$ points uniquely determine an $M$th order polynomial). Then, set up a system of equations that are linear in the desired polynomial coefficients $p_k$:

$$
y_i = p_Mx_i^M + p_{M-1}x_i^{M-1} + \ldots + p_1x_i + p_0, \;\;\;\;\;i = 0, 1, \ldots , N
$$

The least-squares problem can be solved by arranging the measurements into matrix-vector form:

$$
\mathbf{A} =  \left [ {\begin{array}{cc}
 x_0^M & x_0^{M-1} & \cdots & x_0 & 1  \\
 x_1^M & x_1^{M-1} & \cdots & x_1 & 1  \\
 \vdots & \vdots & \vdots & \vdots & \vdots \\
 x_N^M & x_N^{M-1} & \cdots & x_N & 1 \end{array} } \right ],
\;\;\mathbf{y}= \left [ {\begin{array}{cc}
 y_0  \\
 y_1 \\
 \vdots  \\
 y_N  \end{array} } \right ]
$$

The least-squares solution generates the vector of polynomial coefficients $\left[p_M, p_{M-1}, \ldots , p_0\right]$ that minimizes the total squared error in the above linear system. The solution can be computed as:

$$
\mathbf{\tilde p} = (\mathbf{A^TA})^{-1} \mathbf{A^Ty}
$$

It is worth noting that the matrix $(\mathbf{A^TA})^{-1} \mathbf{A^T}$ is also known as the pseudoinverse of the matrix $\mathbf{A}$. You can then use the least-squares polynomial coefficient vector $\mathbf{\tilde p}$ to evaluate the polynomial at any other $x$ values that you wish.


Why do we need to 'estimate' the PSD I've read about the PSD of one OFDM symbol for instance being an estimate of the PSD, but how can it be an 'estimate', when it literally is the PSD of that OFDM symbol.
What is actually being estimated? I get that the OFDM symbol is supposed to be a single record of a WSS random process, and therefore its an estimate of the PSD of the random process because you can't calculate the expected value. But why would you want to calculate the PSD of 'the random process' and not just a single OFDM symbol – surely that's the granularity of focus.

why would you want to calculate the PSD of 'the random process' and not just a single OFDM symbol – surely that's the granularity of focus.

Because you usually don't care about the single symbol as much as you care about making statements on e.g. the channel in general. So, you estimate the PSD.


Hilbert-Huang Transform Implementation I am looking for a implementation of Hilbert-Huang Transformation. Ideally, in Java, but any language with similar syntax would prove useful. Could someone point me to a specific implementation, preferably in Java?

Hilbert-Huang starts with empirical mode decomposition (EMD). I know one HHT code is available on Matlab central but I personally find it not very robust and extremely sensitive to edge effects. A more robust solution may be found here : http://perso.ens-lyon.fr/patrick.flandrin/emd.html

Once the EMD is performed, the missing step to complete HHT is the Hilbert transform which is fairly straight forward.

More details and links are mentioned in this post.


how to get the phasor(frequency domain) representation of a digitalized signal? I know that we can get the phasor representation(frequency domain) of a signal through Fourier transform, but the picture above gives a method that can get the phasor representation of a signal in a different way, in this method according to the writer, it numerically "demodulate" the signal to create a phasor representation of the signal, however, i can't understand that, anyone can help me? Thanks in advance!

This is simply about representing a sinusoidal signal $x(t)$ in the following way:

$$x(t) \ = \ \cos(\omega_0t+\phi) \ = \ a \cos(\omega_0 t) \ - \ b \sin(\omega_0 t)$$

The 'method' described in your question is about identifying the constants $a$ and $b$. The phasor of $x(t)$ is then $X=a+jb$, because

$$x(t)=\text{Re}\{Xe^{j\omega_0t}\}$$

Using basic trigonometric identities you get

$$a=\cos\phi\quad\text{and}\quad b=\sin\phi$$

But you can also compute these values by correlation (as shown in your question), which is equivalent to computing the Fourier coefficients of $x(t)$ (apart from a sign difference for the coefficient $b$). Of course, the only non-zero coefficients are the coefficients corresponding to the fundamental:

$$a=\frac{2}{T}\int_0^Tx(t)\cos(\omega_0 t) \ dt,\quad  T=\frac{2\pi}{\omega_0}\\
b=-\frac{2}{T}\int_0^Tx(t)\sin(\omega_0 t) \ dt$$

These last two equations correspond to the first two steps described in your question.


Image FFT radial component reduction I'm using MATLAB's fft2 to generate the fourier power spectrum of an image. Naturally the spectrum is two dimensional as well, showing the direction of the direction and the frequency of the spectral components. I'd like to get rid of the directional component and just look at the frequency, i.e. reduce the dimension of the spectrum by one.

The way I see it, one can either:

sum over a ring containing all pixels with
$$
r_1 \leq \sqrt{x^2 + y^2} \leq r_2
$$
(which would make it tricky to select the correct boundary radii to really count the pixels that correspond to a certain frequency)
or one can convert the image to a polar coordinate system and sum over the axis corresponding to the radial component of the matrix.

I wonder which method would be preferrable and for what reasons and what normalization I would have to make (of course, for the first method, there would be much less pixels contained in the ring element for lower frequencies than for higher frequencies and I guess some similar scaling problem will occur for the second method)

The two methods are equivalent over discrete images and in both of them it is "...tricky to select the correct boundary radii to really count the pixels that correspond to a certain frequency".

Whether "selecting pixels by distance" or converting from rectangular to polar coordinates, there will always be the problem of what do you do when you have to "land" between pixels. That is, when you have to "land" at an image site, within the acquired data range but between known values.

The answer here is of course "interpolation". Whether you are trying to select all harmonics that are $\approx 14.48765\ldots$ pixels from the centre of the image or you are converting $(x,y)$ pixel $(15,5)$ to its $(\rho, \theta)$ corresponding location of $(\approx 14.41479\ldots , \approx 81.36745\ldots)$, you still have the problem of "landing" between pixels.

You could of course $round(\cdot)$ everything, which would result in a nearest neighbour interpolation but that doesn't look very attractive. You are likely to have to use more complex interpolation schemes.

By the way, by convention, in digital image processing, the two dimensional discrete fourier transform's output is shifted by $\frac{N}{2}$ in both dimensions in order to bring the low frequencies to the centre of the "image". Therefore, you are likely to have to use fftshift too.

As far as the last part of your question is concerned, there should be no need for "scaling" if float or double data types are used, both to calculate the spectrum and its sum over different frequencies. If you are asking what sort of scaling you should do in order to be able to display the image (?), then, please note that you may have to use a logarithmic transformation on the spectrum values before you convert the image to something like unit8 to prepare it for display. Happy to clarify this later part more if more information emerges.

Hope this helps.


Is digitization of discrete valued signal necessary for it to be processed? While reading the book DSP by Proakis and Manolakis i came across the concept of digitization of signals. The Author states that Digitization(rounding  it off to nearest integer) of signal is necessary in order for it to be processed by computer. If i write the following code :-

float signal[5]={ 1.5 , 1 , 4.4 , 3 , 2 , 3 };
for(int i=0;i<5;i++){
    signal[i]*=2;
}

However, In this above code I am able to process the signal without rounding it off to the nearest integer. Please help me clear my doubt.

Your 1.5 and 4.4 are scaled integers (scaled to the nearest 10th).  C float variables are also a form of scaled integer (exponent indicating the scale factor).  So your data is still rounded to the nearest of smaller scaled integers.

In binary digital computers everything gets converted and broken down into 1-bit integers (or sets of them) before processing by digital logic gates.

Analog computers work differently.


How to reverse frequency spectrum of a 2D image? Reversing spectrum of an audio signal is said to be often used in DSP – i.e. applying a low-pass filter is quite a simple operation, compared to a high-pass. So they shift the specrtum by multiplying each other sample by -1, apply the Low-pass filter and then reverse again.

I am wondering if the same processing technique could be applied to images: reverse frequencies, apply edits/processing, reverse again.

Please advice how can a spectrum of a 2D signal be reversed?

Is it something as simple as each other sample inversion as for a 1D audio signal? I'd prefer some code snippets / visual examples to math equations.

I also wonder if FFT of an image, reversing the frequencies, IFFT, FFT, another reversing, IFFT can be considered lossless to the source image?

Upd. I found an amazing github repo "javascript-labs" containing FFT / IFFT implementations in JavaScript. It helped me understand the spectrum inversion better. They just swap quadrants of the FFT image to reverse the spectrum.

Also I discovered to myself that a 2D FT can be parallelized into 1D FT's for vertical and horizontal pixel sets.

My research goal is to find a reversible transform of a source image "S" to another image "T", such as that T carries no visual clues to the content of the original image S. A transformation that is reversible, and also not dependent on precise pixel-to-pixel correspondence. I.e. it should be possible to recover original image S after taking a photo of a printed poster with transformed image T – of course with adequate losses / distortions / color regression. Hypothesis that reversing spectrum of an image would scramble it original content has failed.

I think that you should see it as doing it in one dimension and then the other, meaning you'd invert the values of every other pixel in every row, then once it's done do it again for every column, which I suppose means you can do it directly in 2D by inverting pixels in a checkered pattern like this:

1 -1
-1 1

As for your FFT/IFFT question yes it should be lossless, save for small computation errors, and I'm not sure but it might only work with even sizes so the Nyquist and DC components can switch places. However in both cases you'll have to deal with negative pixel values, so that kind of limits you to your own algorithms, unless you do $f(v) = 0.5+0.5v$ on all your values and then undo it.


Notch filtering line noise outside of frequency band of interest Practically speaking, if one is interested in a frequency band well-separated from the line noise (say, for example, the 10-20 Hz band, with 60 Hz line noise), would it be advisable to notch out this line noise before bandpass filtering to the desired frequency range? It seems that, given of course an acceptable filter rolloff, it would be better not to perform an additional filtering step if it were not necessary; however, I'm wondering whether this may lead to some unforeseen problem.

A general rule of thumb is to filter as little as necessary, because every filter distorts your desired signal, even if just a little bit. Also: a notch filter is advised in a situation where narrow-band interference is present within your signal's band. The notch filter will reject the interference while having minimal effect (ideally) in your signal.

In the problem you describe, my preferred approach would be to go with a low-pass filter only, and increase its order to achieve sufficient rejection of the 60 Hz noise.


What's the correct graphical interpretation of a series of MFCC vectors? I'm studying speech-recognition, in particular the use of MFCC for feature extraction. All examples I've found online tend to graph a series of MFCC extracted from a particular utterance as follows (graph generated by me from the software I'm writing):

As you can see in the graph above:

the x-axis is used for each of the MFC coefficients (from 1 to 12 in this example)
the y-axis is used for the values of the coefficients (ranging aprox from -12 to 42 in this example)
you have as many lines as frames or feature vectors you have extracted (140 in this example).

Now, this doesn't make too much sense to me, because what we are seeing here, is the superposition of all the feature vectors at once, losing any time information. I'm having a hard time to understand how this representation is useful.

In my mind, I would represent the extracted vectors as follows (again, graph generated by me):

In the graph above:

the x-axis is the frame or vector number (1 to 140)
the y-axis are the coefficient values (again, from -12 to 42 aprox)
you have one line for each feature (12).

To me, this representation should be more useful because you can see the evolution in time of each particular feature, and in my mind that should have a stronger impact on how to apply comparison algorithms on spoken words.

Maybe the two representations are equally valid and useful for different purposes, very much like when you need to study a signal in the time domain or in the frequency domain, but in the case of speech recognition I would expect the evolution in time of each individual feature to be more meaningful than the density of values for each feature (and perhaps I'm completely wrong :P).

So, two questions in fact:

Why is the first representation the one than seems to be widely used and not the second one?
When you want to compare two sets of extracted MFCCs, for example by using Dynamic Time Warping - DTW, and related to this topic, do you compare the feature vectors (i.e. 140 vectors of 12 features) or the frames (12 vectors of 140 frames)? (in other words, MxN or NxM?)

Thanks!

The way MFCC's are always used is by feeding them into the classifier. This can be done on a frame-by-frame basis (12x1 vector), or by concatenating (12xN) - same as a spectrogram.

Thus for DTW, you must perform the classification by calculating the distance between 12D vectors. It's a Dynamic Time Warping, so the difference must be calculated between frames in time domain, not between samples in a "single feature domain".

Here is a graphical representation of that:

Just an addendum on MFCC plotting. There is no point in looking at the line plot of all MFCC vectors stacked together (first image), and it's not really useful to look at the time evolution of all MFCC's (second image). No concept of right graphical representation exists for MFCC's unless you want to achieve something specific.

Firstly, a distribution (histogram) of each MFCC can be plotted to see what are the typical values, etc. Especially if you have a GMM model trained for a signal, you could plot the distribution of test sample and overlay with a PDF of GMM in that dimension.

The second approach would be to plot the MFCC's similarly to a spectrogram (as shown in the picture above). One might want to use it for "nice representation" or seeing the overall trends and patterns. But it's almost impossible to tell what kind signal is hiding behind such representation.

For example, here is an MFCC plot of a very simple signal. What signal is that? If you are ready then please hover the yellow space.

   Exponential sine sweep!

Was it easy to tell? Definitely, by looking at the spectrogram it was immediately obvious. Maybe there are few people who guessed since they played with such examples before. Now something more complicated:

Ready?

   Birds chirping - you can listen here

I hardly think so that this MFCC representation is of any use, except filling the space. Now few people might guess the signal based on a spectrogram, but it's impossible to do so by looking at MFCC's.

So if someone is telling you that he can tell what is happening with the signal based on "de-correlated cosine fits of log energies within triangular filter banks in frequency domain", then he is either a liar (with a prior of $99.9\%$) or a computer.  Very rarely there are occasions when looking at raw MFCC values is useful - like debugging and sanity checking (are coefficients constant for stationary signal).


Does fading memory mean impulse response with finite support? Given a continuous-time impulse response $h(t)$, bandlimited to $B$.
The discrete-time $h[n]=h(n/(2B))$ is then a unique and perfect representation of $h(t)$ and a discrete-time system $h[n]$ is then a perfect representation of $h(t)$ (impulse variance).
However, $h[n]$ generally requires infinitely many terms. I am looking for a name of systems $h(t)$ where $h[n]$ has finite terms (that means $h(t)$ has finite support) or can arbitrarily well approximated by finite terms.
Obviously the name for $h[n]$ (the discrete-time equivalent) is a FIR filter.
Is a system with fading memory the right term?
If not, why not and what is the right term?
PS: FIR seems to be a terms almost exclusively used for discrete-time systems. Why does it not exist for for continuous-time systems?

What you want is theoretically impossible: if $h(t) \neq 0$ only on a finite interval, then it is not bandlimited. Fortunately, all real-world signals and systems have the property that $H(f) \rightarrow 0$ as $|f| \rightarrow \infty$. This allows you to approximate $h(t)$ very closely with only a finite number of samples. For engineering purposes, this is often more than enough.
Note that the property above means that $h(t)$ has finite energy and is absolutely integrable (it is in $L_1$). For the math to work out, it is often also required that $h(t)$ is in $L_2$, or that its square is absolutely integrable. These are the terms that are commonly used to describe these kinds of signals and systems.


MTI RADAR, DSP algorithm question Im new to RADAR signal processing.
I using this doc for ref: https://www.mathworks.com/help/phased/examples/ground-clutter-mitigation-with-moving-target-indication-mti-radar.html
I trying to rebuild the code there, I dont need all the definitions of the antennas and signals because I already have file with I/Q pulses example and I need to process them.
I stuck on the part where I should do [Matched Filter] -> [3-Pulse Canceller]
The way I do it, which i told is wrong, but the result just look better for me is:

calculating Matched Filter for every pulse (assuming I working with 64 pulses) and I arrange the filtered data on Matrix like this:

take the matched filter result row by row and entering the 3PC filter like this:

and than continue the processing..
the matlab code:
%% Matched Filter
mf_buffer_Q = filter(coeff,1,buffer_Q,[],2);
mf_buffer_I = filter(coeff,1,buffer_I,[],2);

%% Three-Pulse Canceller
h = [1 -2 1];
TriPC_buffer_Q = filter(h,1,mf_buffer_Q,[],2);
TriPC_buffer_I = filter(h,1,mf_buffer_I,[],2);

the plot result:

the other way to do it is what i think is wrong but I told is the right way is:
%% Matched Filter
mf_buffer_Q = filter(coeff,1,buffer_Q,[],1);
mf_buffer_I = filter(coeff,1,buffer_I,[],1);

%% Three-Pulse Canceller
h = [1 -2 1];
TriPC_buffer_Q = filter(h,1,mf_buffer_Q,[],2);
TriPC_buffer_I = filter(h,1,mf_buffer_I,[],2);

that means:

the plot result for this:

what output seems better?
what is the right way to do this 2 filters? one horizontal and the second vertical
or I just should do both filters horizontal.
thanks.

The second way is how it is done.
In the fast-time (or range-bin dimension) you are right to perform the matched filter.
However, Doppler information is gathered from sampling pulse-to-pulse. The MTI filter essentially subtracts pulses so that if they have similar phase, or the phase is not changing at all, then there is little to no Doppler being generated and thus the target is not moving.
In other words, the three-pulse canceller applies a filter across the Doppler dimension for each range bin, as the second pulse matrix suggests.


Why do Dyadic filterbanks downsample the high pass signal portions? I'm currently programming a dyadic filter bank and have a question. I notice in all of the visual representations:

(from here (Dyadic Analysis Filter Bank)), the high pass filtered output for each stage is downsampled. That is, when the signal is split into two subbands, the low pass band is downsampled and so is the high pass band.

The question: why does every representation insist on downsampling the high pass subband?

I understand that you've just halved the frequency content (by half banding the signal), thus we can downsample by two and still satisfying the sampling theorem (sample rate $\ge$ $2\times$bandwidth). However, what confuses me is the fact that we've still aliased the high pass signal such that the frequency content is 'mirrored' around Fs/2. Essentially all of the high frequencies in the highpass subband now correspond to low frequencies in the aliased signal, and all the low frequencies in the highpass subband now correspond high frequencies.

Doesn't this make the analysis more difficult? As in, when you want to apply a highpass filter to the highpass subband, you are now required to apply a lowpass filter to the downsampled aliased version?

[Good question, that made me rethink of stuff I believed natural. I shall incorporate them in future lectures]

Downsampling the highpass (and the lowpass) provide you with a critically sampled filterbank, in other words, not redundant. Critical sampling is a price to pay in some applications like compression, and in the context of finding an orthogonal transform. So yes, this is a limitation: one should take care of aliasing, and of the imperfect character of the designed filters. And yes too, filterbank coefficients are subject to folding/mirroring. If one were to analyse these coefficients as if they were a genuine signal, for instance with linear filters or Fourier transforms, that could raise some concerns.

But remember that the DFT can be interpreted as a filter bank as well, downsampled by the length of the signal (plus zero-padding). Each DFT coefficient becomes, in this interpretation (correct me if wrong), a single "sample", rejected to a baseband.
I understand your concern, but it is very rare to apply a linear filter to highpass subbands. On the contrary, one usually employs nonlinear filters and selection tools.
Would you perform linear filtering or Fourier analysis on it without caution? Of course, if you keep this coefficient, set the others to zero, and inverse the DFT, you get a signal back, a sine, and then you could perform traditional processing on it.

The same applies to FIR (even IIR ones) filter banks. Here, the downsampling is much smaller. And the series of coefficients longer. But in many cases, these coefficients are not "signals" per se. Very often, they are not processed "linearly": they are cancelled, shrunk, thresholded. Most of the work focuses on both amplitudes AND  location, due to the relative locally of the filters (contrarily to DFT).

If you want to keep with the linear interpretation, the best is probably to set the other coefficients to $0$ and reconstruct with the inverse FB. Then, you have a real signal back again, and you can safely Fourier, filter, etc.

But not all wavelet representations do that: "I notice in all of the visual representations" and "why does every representation insist on downsampling the high pass subband?". This is not always the case, and you have a wealth of oversampled, shift-invariant, etc. filter banks, not only dyadic, that do a much better job.

For instance, you can avoid the downsampling on the highpass (and also on the lowpass), and you end up with a more redundant representation, that can be useful for some processing cases. Then, you can process the subbands more normally.

For some (personal) background on (oversampled or not) filter banks, which content I can speak of: in 1D Optimization of Synthesis Oversampled Complex Filter Banks, in 2D A Panorama on Multiscale Geometric Representations.


Frequency response of filter  signal mix This might be a weird question but here's the setup: I have a few biquads that filter a signal $x[n]$ and output the filtered signal $y[n]$. I can calculate the frequency response of those biquads with the z-transform. Now I mix $x[n]$ and the filtered signal $y[n]$ such that the mixed signal $y'[n]= a \cdot x[n] + (1 - a) \cdot y[n]$ where $0\leq a\leq 1$. How would I go about calculating the frequency response of $y'[n]$?

The Z transform is linear. So if
$$ Y(z) = H(z) \cdot X(z) $$ then
$$ Y^{'}(z) = a \cdot X(z) + (1-a) \cdot Y(z) = a \cdot X(z) + (1-a) \cdot H(z) \cdot X(z) = (a + (1-a) \cdot H(z)) \cdot X(z) $$

Hence
$$ H'(z) = a + (1-a) \cdot H(z) $$


Understanding the frequency scale of a spectrogram The graph below was derived from a raw seismogram recorded during an earthquake over a timespan from t=0 to t=1400 seconds (not shown on x-axis).

The original seismogram $s(t)$ is not shown, but it is known that the original signal is the convolution between the source wave $w(t)$, the effects of the medium $g(t)$, and the effects from instrument itself $i(t)$; that is,

$$s(t)=w(t)*g(t)*i(t)$$

The graph above is what remains after removing the instrument response as shown below:

I think this instrument response is telling me that the amplitudes of the original seismogram $s(t)$ at frequencies higher than $0.02Hz$ were amplified at a constant gain of approximately $1E9$, while amplitudes lower than $0.01Hz$ are effectively dampened (by gains significantly smaller than at those frequencies higher than $0.02Hz$).

But what exactly can I effectively say about the shape of the original seismogram from this information?  If these frequencies are 'dampened' in the broadband representation above, what, if anything, can be said about the amplitudes and shapes observed in the original graph?  Can I simply say that longer waves should be present in the original signal but not say anything about when they occur?  Could I somehow convolve these graphs to generate the original signal (or at least closely approximate it)?

I think what you're saying is that you have a linear, time invariant (LTI) system model for both the propagation medium and the response of the measurement instrument itself. Given that you know the (at least approximate) magnitude response of the instrument, you would like to best estimate what the actual signal at the instrument's input looked like.

This is similar to the process of equalization used in communications systems. Since you have some knowledge of the system (in this case characterized by the impulse response $i(t)$ or the frequency response plot that you showed), you can attempt to undo the effects on the signal that are imposed by the measurement instrument. More generally, this type of operation is known as deconvolution: given an output of a convolution and optionally knowledge of one of the inputs to the convolution, one would like to calculate the other input. In this case, one of the convolution inputs is $i(t)$ and the output is the measurements reported by your instrument.

It appears that you don't know the exact impulse response of the measurement device, only its magnitude response. This isn't enough information in itself to execute the deconvolution, so you'll need to make some assumptions about the phase response of the instrument in order to move forward with that approach. This is a very powerful technique, however; as the Wikipedia article notes, the approach can be used to compensate for the effects of the seismic propagation medium also (since it's also modeled as a convolution with an impulse response), potentially allowing you to gain more insight as to the original seismic event.

Edit: And, as to your questions about how you're interpreting the plot, your understanding seems correct. The system has a highpass characteristic with a flat passband starting at approximately $0.02\ \text{Hz}$, with the magnitude response decreasing by a factor of 10 per decade of frequency below $0.02\ \text{Hz}$. This could be modeled well as a single-pole highpass filter.


What are the state of the art methods for person segmentation and person pose identification I know that HOG is pretty state of the art for person detection, but in its original form HOG is only a detection method, and I don't believe it is particularly fast either (there is a significant speedup which takes the sliding window and rejects unlikely possibilities early, but it has been patented).

What is the current state of the art approach for determining the pose of a person (including close poses like a head shot, and hand shot, etc...), and for segmenting the person from the background?

The setting here is still images.

Check out Yang and Ramanan's work with flexible mixture of parts.  It starts with a central part, but pairs of parts that are attached are detected jointly by estimating the probability that the image area being queried supports the co-occurence of the attached parts in the image.  Error is presumed in the attachment and is minimized using a spring-like model.  Each part is detected using HOG-style techniques, and all pairs of parts in the human 'tree' are iterated over to determine the detection.  In my experience it has a very good false alarm rate, but requires hardware acceleration for real-time applications.


main direction in 2D map of data It may be a very simple question. We would appreciate any comments, guide or complete solution. We prefer however incomplete solutions but novel practical ideas.
 - How to find the main direction of variation i.e., main diagonal in the following example figure (north-east,south-west)?
Any algorithm and idea for coding or snippet is welcome.

Updates:
We found that the procedure mentioned in the accepted answer is almost PCA (Principal Component Analysis) in its simplest implementation. PCA worked very well. The resulting direction is perfectly matched with what we could expect, visually.

You could compute the covariance matrix from your elevation data. The eigenvector that belongs to the larger eigenvalue will give you the main direction of variation:

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


If Digital Signal has frequency components then why we don't use it As we know there are two kinds of signals: analog signals and digital signals.

And there are two common techniques for multiplexing the signals: Time Division Multiplexing (TDM) and Frequency Division Multiplexing (FDM).

As I understand it, FDM is commonly used to multiples analog signals, while TDM is commonly used to multiplex digital signals.

What I don't understand is:  Why do we not (commonly, at least) use FDM to multiplex digital signals or TDM to multiplex analog signals?

Analog and digital
First off, you need to understand what an "analog" signal is, and what a "digital" signal is, how they are different, and how they are similar.
The term "analog" comes from the old distinction between "analog" and "digital" computers.  A "digital computer", even a very primitive one of decades ago, has always been more or less what we know now -- a box that manipulates "bits" (or, in very old systems, decimal digits), where the information is represented by numbers encoded in the bits.
An "analog computer", on the other hand, used voltage values in the computer to represent numeric quantities.  Ie, where the digital computer might represent the value "6" as the binary bits 110, the analog computer might represent the value "6" as a voltage of 6 volts.
(I could go on here for several paragraphs explaining more about how an analog computer works, but it suffices to say that, except in very limited applications, analog computers are no longer used -- where needed their function can be simulated with digital computers more cheaply and accurately.)
But the point is that an "analog signal" is called that because the voltage levels in the signal are analogies ("analogs") for some "real-life" phenomenon.  Eg, in an analog audio signal the voltage level is an analog for the sound pressure that would create the sound -- a higher voltage level means more sound pressure, etc.  Further (and it will be important later), the "X-axis" -- the time dimension -- is also "analog".  That is, in most cases (except where, eg, a recording has been sped up or slowed down) one second of time in the analog signal corresponds to one second in "real life".
In the image below, peaks in the signal represent higher air pressure (and valleys represent lower air pressure).

On the other hand, the ups and downs of a digital signal have no direct relationship to "real life".  Rather, they represent binary digits (hence "digital") in a numeric code used to represent whatever "real life" phenomenon the signal corresponds to.  In the diagram below, the bits (each line represents a separate bit stream) might represent characters of ASCII text, or they might represent readings from a heart rate monitor.  Looking at the data you can't even tell if it might be a pleasant sound (like the guitar), or a discordant sound (like the cymbal).
Also note that the X dimension is not an analog for anything -- an entire symphony might be transmitted as a digital signal in ten seconds, or a high-resolution photograph that was snapped in an instant may be transmitted, bit-by-bit, over a period of seconds or minutes.  (The old Mariner space probe would take hours to transmit a single high-resolution image, because the bit rate was so slow at that distance.)

But note that both types of signals can be represented by a fluctuating voltage on a wire, and, in fact, you can hook a digital signal (if the data rate is slow enough) directly up to an ("analog") audio amplifier and listen to it.  (It will come through as a noisy squeal).  So at this fundamental level -- a voltage that fluctuates over time -- analog and digital signals are basically the same.
(It also needs to be noted that, from the standpoint of "information theory", the maximum amount of information that can be encoded into a signal occupying a given bandwidth is the essentially the same, regardless of whether it's digital or analog, or how it may be multiplexed.  However, while one could spend days discussing information theory and entropy I won't go into it further, because it's for the theoreticians, and because it makes my head hurt to think about it.)
Multiplexing
Next we need to understand "multiplexing".  In the US, at least, a "theatre multiplex" is a group of motion picture theatres in a single building.  Similarly, when you "multiplex" signals you combine several signals into one all-encompassing signal.  But just as you wouldn't want to go to a motion picture theatre and see 6 movies all projected onto a single screen simultaneously, when we multiplex signals we do so in a way that the receiving station can separate the combined signals back into the individual ones.
The earliest form of signal multiplexing(*) was common radio.  Several radio stations in the same area can all broadcast their signals into the same "ether", but each uses a different assigned frequency, so that, by tuning your radio receiver to a specific frequency (actually a frequency band) you can receive the "hard rock" station you want to hear, and not have to simultaneously hear the "easy listening" station your parents listen to.
This, in essence, is Frequency Division Multiplexing (FDM) -- the available spectrum is divided into frequency bands and each signal is combined with a "carrier" frequency to occupy the appropriate frequency band.  (The original signal, before multiplexing, is referred to as a "baseband" signal, because it has not been combined with a carrier.  This sort of contrasts with "broadband", where multiple signals have been multiplexed, but unfortunately the term "broadband" has multiple confusing definitions.)
FDM used to be used extensively in the phone system, for long-distance lines.  A single pair of wires (and associated amplifiers), which might have a frequency response from 20Hz to 100KHz, carried a dozen or so audio signals, each at an assigned frequency just like broadcast radio or TV.  Later coaxial cables were used, increasing the available bandwidth to several MHz and allowing several hundred audio signals per cable.
(A significant drawback of this scheme was that many amplifiers were required -- one every ten miles at least, and one every mile for some schemes.)
Now, it should be clear that, given that there's no fundamental difference (electrically) between an analog signal and a digital signal, digital signals can indeed be broadcast using FDM.  In fact, if you have a cable modem, the cable, with over a hundred analog channels, is frequency multiplexed (just as over-the-air TV is).  Your cable modem signal occupies (sharing with your neighbors) one of those "analog" channels.  So your cable modem signal is being transmitted via FDM (though with several added twists).
Now for Time Division Multiplexing (TDM).
Time division multiplexing involves what it name implies -- the time axis of the combined signal is divided into discrete segments, and different individual signals are transmitted in the different segments.
TDM actually goes back to before regular digital data.  The Bell Telephone company in the US experimented in the 50s with chopping up analog signals into lots of thin slices and interleaving the slices from several different analog signals together, allowing the several signals to be sent over a single pair of wires without having to resort to FDM.  So in this case analog signals were transmitted via TDM.
But doing this was quite complicated.  The analog signal needed to be "sampled" rapidly enough that the signal didn't change much between samples.  For an audio signal with a maximum frequency of about 3.5KHz (the range for common telephone lines) sampling must be done at at least twice that frequency (so about 7 thousand "slices" per second) and then those slices must be interleaved with others (maybe five or ten).  (This was known as "pulse amplitude modulation".)  At the other end the slices need to be separated, and the missing space between slices filled in by "smoothing" the signal.  And all this needs to be done without seriously changing the amplitude (sound volume) of individual slices.
Of course, the above scheme was very difficult to get right and didn't really provide much advantage relative to FDM over the same wires, so when true digital capabilities came around they dropped it like a hot rock.
Why digital is different
So far I've been showing how there's really no difference between analog and digital.  Now I'll reverse course -- digital is really different in some very important ways.
First off, digital signals don't degrade like analog signals.  Old analog telephone networks introduced a lot of distortion into the sound as it was transmitted cross-country and reamplified many times (as often as once every mile).  A digital signal does develop distortion in transmission, but since it's digital it can be "regenerated" at the end of each hop.  This alone made TDM much more practical for digital, since there was no longer a need to maintain accurate signal amplitude, and the cables and amplifiers become simpler and less expensive (and farther apart).
But that's not all.
Consider, for instance, if you use Skype to talk to someone on the Internet.  Your voice is first converted to an analog signal by the microphone in your computer.  But it's then immediately digitized, and the digitized bits are "buffered" -- stored up awaiting transmission.  Algorithms in Skype will examine the digitized sound and discard periods of silence (since there's no need to transmit silence), and they'll otherwise process the sound to minimize the number of bits needed to represent it.
Eventually (after a seeming eternity of computer instructions, but in an instant from your ear's point of view) a group of sound bits will be formed into a block, have the appropriate internet headers appended to the front, and be queued for transmission over your Internet connection.  But, since you happen to be simultaneously using your web browser, the sound data may have to wait a few more milliseconds for your web request to be transmitted.
On the receiving end, Skype accepts the packet of transmitted data from the receiving computer's Internet connection and sticks it back in a queue, taking care to insert it in the right order in the queue based on a "timestamp" in the data packet.  (This extra care is necessary because data packets may be delivered out of order in some cases.)  Then another part of the Skype program removes packets from the front of the queue and "schedules" them to play back the recorded sound at the appropriate time based on the timestamps.  So the received data is not played immediately, but instead waits (perhaps a few more milliseconds) until "the time is right" to be played, so that the time dimension of the original data (and hence the fiction of continuous voiced transmission) is maintained.
So what happened here?  The voice data from Skype was time-division multiplexed with your web browser requests, and, once out on the Internet trunk lines, both were time-division multiplexed with thousands of other messages.  But, critically, the multiplexing was not based on a fixed schedule (eg, dividing the time dimension into 10 millisecond slices and "dedicating" the channel to transmitting each signal for 10 milliseconds), but instead an "adaptive" algorithm was used, where each signal only occupied the channel for as long as it needed.  Even though this technique requires that each data packet contain an address to identify it's destination (whereas with a fixed multiplexing scheme the destination would be defined by the slice location), it's still far, far more efficient than a scheme that relies on fixed slices, either in time (primitive TDM) or frequency (FDM).
Summary
Now note the differences between TDM (as used in networking systems) and FDM.  FDM has the advantage of being relatively simple (witness standard AM radio).  It can carry analog or digital signals equally well, but, because it does either equally well, it offers no great advantages to digital.
TDM is only marginally practical for (undigitized) analog signals -- the complexity and potential loss of signal quality are major negatives.  But digital, it seems, was "made for" TDM -- the data is already "sliced", and accurately maintaining amplitude or time reference is not required.  And since "dead air" is represented by a complete absence of data packets, not by simple silence on a phone line, "adaptive" multiplexing, where only "useful" data is transmitted, and not on a rigid schedule, is much more practical.  For these reasons much more digital data can be transmitted in a given bandwidth using (adaptive) TDM vs FDM.
So it's not that digital CAN'T be transmitted using FDM and analog CAN'T be transmitted using TDM, but rather there are practical reasons why doing so would be more expensive, less efficient, and less reliable.
(I need to add that even given all the above FDM is used for digital data in some circumstances.  For instance, the link to/from communications satellites is very "broadband", and no single data stream could practically use the entire bandwidth.  So the overall channel width is sliced into frequency bands, with each occupied by a separate data stream.  In this case, since the digital data stream is already time-division multiplexed and there is very little "dead air", the advantage of further TDM would be slight.)

(*) There were some crude TDM schemes used for telegraph signals as far back as the 1850s, but that's mostly lost in the mists of history.


image type after an ifft reconstruction I reconstruct images from MRI k-space using ifft and root-sum-of-squares method.

imRef = ifftshift(ifftshift(ifft(ifft( ifftshift(ifftshift( kspacedata,1),2),[],1),[],2),1),2);
imSOS    = squeeze(sqrt(sum( abs(angle(imRef)).^2, 3)));
imagesc(abs(imSOS))

k = abs(imSOS);
disp(class(k)); %double

Class of image matrix k is shown as double.

Is the image matrix k an indexed image or intensity image?

The sample values from reconstructed image k is shown below.

4.2753    4.9807    4.5435    5.6548    6.1303    5.0229
3.3805    4.6260    5.1594    5.2692    4.1187    4.4885
5.8990    5.5275    4.3493    5.6182    6.7237    4.3071
6.4885    4.6861    4.4086    3.5034    5.2378    4.7466
6.1099    4.6995    4.1673    4.7408    3.8915    5.7531
5.4006    3.1289    5.6541    5.8782    4.6568    3.8166

By referrring http://in.mathworks.com/help/matlab/creating_plots/image-types.html,
I learned that for indexed images, the data matrix is represented by integers which are indices to color map. And if it were intensity images, I feel the data range is somewhere between 0 and 1. If k represents an indexed image, how can I know the colormap associated with it?

Thanks.

FFT and IFFT are linear operators, and as such, the results only make a lot of sense in a linear intensity space, not if indexed into a non-linearly mapped space.


Why doesn't sampling a periodic continuous-time signal yield a periodic discrete-time signal? I have been studying signals and systems lately and I have came across the following claim:

  The uniform sampling of a periodic continuous-time signal may not be periodic!

Can someone please explain why this statement is true?

If the ratio between your sampling frequency and the frequency of your signal is irrational, you will not have a periodic discrete signal.

Assuming you have a 1-kHz sine wave and you sample at 3000*sqrt(2) Hz. You will have approximately 4.2 samples per period. However you will not be able to sample the sine wave exactly at the same place. Hence your digital signal will not be periodic.

However, if you sampled the same 1-kHz signal at 4 kHz, you would get a periodic discrete signal. The period would be 4 samples.


Hilbert transform Use as a Single Sideband Frequency Translator in Analog communication How Hilbert transform helps us in communication? What is the practical significance of it? Specifically how is the Hilbert Transform used to create single side-band frequency translators?

There are many applications for the analog equivalent of the Hilbert Transform, often implemented with “90 degree power dividers” also referred to as “quadrature splitters” as well as broadband variants implemented with phase tracking networks.

One application is single sideband frequency translators done with IQ mixers; I have another posting that details this further here:

Frequency shifting of a quadrature mixed signal

Another application is phase shifters and vector modulators where the signal is split in quadrature (To I and Q components) then combined in 0 degree phase after weighting the I and Q components according to the desired phase and amplitude. This can provide frequency translation similar to the approach with IQ mixers. The architectural difference is the use of absorptive bi-phase attenuators in each quadrature leg instead of reflective mixers which is more appropriate for static phase and amplitude control.

Other common applications for the Hilbert Transform that are done digitally can also be done in the analog with the purpose of converting a double sideband IF signal to be single sided (positive or negative frequency only) or creating a baseband analytic signal (which must be complex).

The link I gave above that demonstrates using the Hilbert Transform to frequency shift a signal is converting the signal to single side band frequencies only. Since the OP further asked this question I will explain that further in this context.

First consider Double Sideband frequency translator which would use a single mixer and as shown below is demonstrated as an upconverter translating a sinusoidal tone with frequency $\omega_s$ to a carrier of frequency $\omega_c$. Here we see that both sidebands of the baseband signal with frequency $\omega_s$ have been translated up to be centered on $\omega_c$:

In comparison, the implementation for Single Sideband frequency translator uses two 90° power dividers as shown in the diagram below. I have designated one output for each of the power dividers to represent the imaginary axis by putting a label $j$ on those ports just to show similarity to the Hilbert Transform: there is no reason I cannot define one output to represent such even though many may argue that both are real signals-- whether a signal is real or imaginary rests entirely in our definition and consistency with our rules of math. Regardless the same block diagram would similarly work with real signals, and is generally considered as such. If you swap one of the quadrature legs the resulting output would be the other sideband.


What do we call the areas with "1" in a binary image I am looking for possible professional terms. In a binary (1 or 0) image, normally, the areas with "1"s are where the information is. I wonder if there is any professional name for it in image processing. Thank you!

This problem born when you use pixel intensity as physical mass, for example to calculate mass center or inertia moments, in this case the background should be 0 and object should be 1.

Btw your question is not clear because in a binary image {0,1} have the same information content and generally there are not difference.

There is not a name for 1 pixel and 0 pixel, often in common convention 1 is white and 0 is black, but it's not a general rule. Maybe you can call 1 high value pixel and 0 low value pixel but I think that are not common names and they can make confusion in some contexts.


RF for Data Transmit and Receicve for Beginners Is there any good matrial avialable for data specific DSP? I have gone through some great DSP intros for beginners in python but all of them revolve around the characteristics of audio signal

A possible combination is a copy of Digital Communications by Proakis & Salehi along with GNU Radio.

With SDR, you actually create / process the signals that carry your data over the RF medium.

Therefore, to an extent, "freshman level theory" is all you are stuck with if you have to write an end-to-end SDR system.

But, you could keep the book (or, a book) as a reference for points you don't understand when going through something like GNU Radio's PSK Demodulation resources for example.

I am not sure what your expectation is from SDR, but taking the comment (about SS7) as a starting point, it's not like the SDR has an API which you use to simply "send data". Although, I am sure that you could add that layer too once you have your basic DSP dataflow set up.

The beauty of SDR is that you take full control of every building block from transmission to reception. So, there is nothing stopping you from writing a Phase Locked Loop (which you are going to need) from first principles or use one of the existing ones (or here).

But you still have to have some awareness of what is the "loop bandwidth" for example, in the case of the PLL.

Hope this helps


How to find out whether guassian or salt& pepper noise present in an image, Is there any algorithm Is there any algorithm for identifying the noise in an image
I have to identify whether the noise present in a mammogram image is guassian or salt& pepper.

Exploit the fact that salt&pepper noise has only 2 values and median filter does excellent job in getting rid of it.

Filter your image with median filter
Find pixels where median filter did most of its job
Plot histogram of those pixels
Compute fraction of 0s and 1s

threshold = 0.2; % adjust to perform best on your images
im = imread('../lena.bmp');
imd = double(im);
imd = imd./max(max(imd));
imdg = imnoise(imd,'gaussian');
imdsp = imnoise(imd,'salt & pepper');
subplot(2,2,1), imshow(imdg)
subplot(2,2,2), imshow(imdsp)
imdgf = medfilt2(imdg);
imdspf = medfilt2(imdsp);
g = imdg.-imdgf;
sp = imdsp.-imdspf;
gm = abs(g)>threshold;
spm = abs(sp)>threshold;
subplot(2,2,3), hist(imdg(gm))
subplot(2,2,4), hist(imdsp(spm))
waitforbuttonpress


Methods for calculating signal quality metric i'm sure this is a fairly large topic area, but I'm not quite sure where to start. I have a signal in which I can change several parameters, such as the frequency, gain, etc. The modulation is a differential BPSK, and my current metric for determining the signal quality is what percentage of my auto correlation (correlating the current bit signal with the previous bit) is above 0. If it's all above or below zero, a 100% signal quality. if it's half and half, then it's 0%, and then it's linear from there. Is there a better way of doing this? I'm wanting to be able to change my signal parameters in order to increase my signal quality, where a higher signal quality should be a higher reliability of receiving a message correctly.

I'm also looking at doing demodulation in the frequency domain, so then my calculations would probably need to change. Any input and references would be appreciated. Thanks!

From your question it sounds like you are looking for a SNDR metric.  In that case you are measuring a distortion level by comparing the distorted signal to a reference.
For a discrete-time signal $y$, comparable to reference $r$, the SNDR would be:
$$ SNDR = \frac{E[r^2]}{E[(y-r)^2]} $$


human speech noise filter Does anyone know of a filter to attenuate non-speech? I am writing speech recognition software and would like to filter out everything but human speech. This would include background noise, noise produced by a crappy microphone, or even background music. I have already implemented a first order filter that compensates for the 6 dB roll-off of the power spectrum, but I'm still hearing noise (though the speech sounds a lot clearer). I have thought to use a low-pass filter, but I am iffy about doing that for two reasons:

I don't know whether or not a low-pass pre-filter will interfere with the rest of the speech processing. Even though the human ear can only detect sounds lower than 20 kHz, I don't want to risk eliminating any higher order harmonics that might be necessary to process the speech (though I don't know if this is the case or not. But I don't want to take any chances).
I understand that the excitation of certain consonants (such as f, h, and s) are almost entirely white noise. I don't want to implement a noise filter that will eliminate good noise, so to speak.

Ideally, I would like to be left with only the speech of the person speaking into the microphone. If you have any ideas, or there is something that I'm missing, please let me know. Much appreciated!

A speech communication channel as used in telephony typically has a frequency response of 300 Hz to 3 kHz. Although this rejects a lot of the energy in normal speech, intelligibility is still quite good - the main problem seems to be that certain plosive consonants, e.g. "p" and "t", can be a little hard to discriminate without the higher frequency components.

So you're probably looking for a "sweet spot" somewhere between using the full 20 Hz - 20 kHz bandwidth typically found in consumer audio and the most aggressive filtering used for voice comms (see above). I would suggest starting with a bandpass filter from say 50 Hz to 8 kHz. It will probably only improve SNR by a few dB at best, but it may help, particularly if you have a lot of high frequency background noise.