texts: Breaking a convolution into smaller pieces

texts: For a project I need to do convolution and i use gpu for calculations. Sometimes I have to deal with kernel sizes of 50x50 and this size of kernel is sufficiently large that it chokes the gpu. (not enough memory svailable by the gpu) I need to find a way to break the kernel into smaller sizes (8x8 or similar) and do the convolution this way (i.e. piece by piece and later stitch together) so that I can realize gpu enhancements. What is the best way to break these kernels into smaller pieces? (I have no prior knowledge of kernel size so the solution must be something I can deal with in run time. )


texts: Does a linear phase FIR filter shifted in frequency preserve linear phase?

texts:
Let's say I have a symmetric (and therefore linear phase) FIR low pass filter with real coefficients.

If I then shift this filter in some direction in frequency by multiplying its coefficients with a complex exponential (to get a band pass filter with complex coefficients), will it still be linear phase?


texts: Checking convolution property II of z transform

texts: I have two sequences x and y of lengths, say 5 and 10.  I multiply them in time domain element by element.  I get another sequence.  Now this as per the convolution theorem should be equal to the inverse z tranform of the convolution of the two z transforms of the original sequences.  Now, I know how to get the ZTs of the original sequences.  How do I get the convolution of  these ZTs?


texts: Harmonic Product Spectrum algorithm. I dont understand one step

texts:
I was looking some Harmonic Product Spectrum algorithm examples and I came up with this:

//Implement Harmonic Product Spectrum
for(k = 0; k < buffer_size / 8; k++)
{
    sum[k] = magFFT[k] * magFFT[2*k] * magFFT[3*k];
    // find fundamental frequency (maximum value in plot)
    if( sum[k] > max_value2 && k > 0 )
    {
        max_value2 = sum[k];
        fund_freq = k;
    }
 }
 fund_freq1 = fund_freq * 8000 / buffer_size;

What I dont understand is the reason behind the last line, I tested it and it works but I couldn't find an explanation for that multiplication. Can anyone help me please?

Source of code: SOURCE


texts: Fourier transform of cosine to the power of 3

texts:
How can I find the Fourier transform of

$$ f(x) = ( \cos(x) )^3$$

I know that for $ g(x) = \cos(x) $

$$\mathcal F \Big\{ g(x) \Big\} = \mathcal F \Big\{ \cos(x) \Big\} = \pi \Big [ \delta(w-\pi / 2) + \delta(w+\pi / 2) \Big ]$$

But using this pair of Fourier transform how to obtain the $ F \Big\{ f(x) \Big\} $ ?? Is there a direct/simple way to do that?


texts: retrieving original data from phase and magnitude of Fourier transform

texts:
I use this snippet of python code to transform data to Fourier phase and magnitude and then retrieving original data.

from random import randint as RI
import numpy.fft as FFT
import math

w = 4
h = 4
random_range = 255

vals = [[] for i in range(h)]

for i in range(h):
    for j in range(w):
        vals[i].append(RI(0,random_range))

fftc = FFT.fft2(vals)

magnitudes = [[] for i in range(h)]
phases = [[] for i in range(h)]

for i in range(h):
    for j in range(w):
        magnitudes[i].append(math.sqrt(fftc[i][j].real**2+fftc[i][j].imag**2))
        phases[i].append(math.atan(fftc[i][j].imag/fftc[i][j].real))

for_ifft = [[] for i in range(h)]

for i in range(h):
    for j in range(w):
        rl = magnitudes[i][j]*math.cos(phases[i][j])
        im = magnitudes[i][j]*math.sin(phases[i][j])
        for_ifft[i].append(complex(rl,im))

ifftc = FFT.ifft2(for_ifft)

print vals
print ifftc

When i compare outputs, results are different. for example two typical outputs are:

//original data
[[115, 4, 33, 91], [228, 123, 46, 111], [83, 227, 81, 55], [229, 188, 138, 41]]

//retrieved data from phase and magnitude informations
[[ 309,182.25,105,57.25],[ 33,118.75,120.50,116.25],[ 148,128.75,79,94.75],[ 165,112.25,-6.50,29.75]]

I know the theory behind it. for example for calculating phase and magnitude we do this :

Magnitude = Square_Root(real*real+imaginary*imaginary)
Phase = Inverse_Tangent(imaginary/real)

and for calculating complex number from magnitude and phase we do this :

Real = Magnitude * cos(Phase)
Imaginary = Magnitude * sin(Phase)

So, what is the cause of problem. Is it a data scaling issue.


texts: $ {L}_{0} $ Pseudo Norm Minimization in Compressive Sensing

texts: I have recently taken up studying compressive sensing related papers. Some things are not very clear to me or may be I am not able to visualize the scenario as is said. Like how $l_0$ norm minimization is NP hard? Can somebody explain with an example (step by step example may be?) of a matrix and vector for an equation $y=Ax$, how $l_0$ minimization of the same can be NP hard?


texts:
I'm currently testing the performance of an energy detector using two tests. One test involves finding the optimal threshold to achieve a probability of detection(Pd) of 0.9 and a probability of false Alarm (Pfa) of 0.1.

Having completed the first test, the next test I want to do is to see how varying SNR affects both Pd and Pfa, given a constant threshold. If Pd is the probability of the detector detecting a signal given the signal of interest is being transmitted and pfa is the probability of the detector detecting a signal given the signal of interest is not being transmitted how exactly can I vary SNR to get Pfa. I'm sure its something that I'm not quite understanding but if the signal of interest is not being transmitted and all that the receiver is picking up is noise, how is there a SNR in the case of Pfa?


texts: If sampling rate is high does that increase the error or decrease the error?

texts:
Suppose we are tracking a computer mouse on the screen, if we increase the frequency of the sampling rate(to very high values) then would that increase the error in getting the exact position of the mouse or decrease it?

Would we get lots of outliers and we would have to average out?

What I am trying to ask will it get smooth lines or every little flicker in movement will be noted?


texts: Classifying animated GIFs as animation or film

texts:
Is there any known way of differentiating between frames of real footage and animated content?  It seems like there's probably range of colors more likely to appear in a histogram for real footage, or maybe a level of color variation more likely to be in animation.

The case I'm looking at is converting GIF images to h264, so there's a max of 256 colors to work with.  I'm mostly trying to figure out when I should use the animation profile so the output is screwed up as little as possible in the transcode.


texts: How do DAWs align the phase when working with plugins?

texts:
I'm no signal processing guru, so I'm looking for a pretty basic explanation that I'll understand. I've seen this question about Phase Locked Loops, but I have a hard time understanding that, especially in relation to digital audio workstations (DAWs).

What I want to know is how DAWs align the phase of signals when sending them through AU or VST plugins. How do they know how much the plugins will delay a signal they send through it?

I may be mistaken in my assumption that any correction is actually necessary, but then how do DAWs correct for delay of each channel in parallel processing?


texts: Properties of a Kalman filter of a non controllable system

texts:
Let's say I have a standard system
\begin{align}
x(t+1)&=Ax\\
y(t) &=Cx(t)
\end{align}

As you can see $B=0$, so the system is not controllable.

For the steady state Kalman filter I'd say that if $A$ is Hurwitz then it converges.

But what if $A$ is say,  $\mathbf I_n$ (i.e. eye(n)), what can I say about the optimality of the SS.  filter varying $C$?


texts: Resampling and removing high frequency noise?

texts:
I am currently working on a simple sampler that will allow me to load in a wav file and use my MIDI keyboard to play the loaded wav sample at the frequency according to the note played.

Now I need this pitch shifting to be done in real time so I opted to use resampleing to achieve different pitches. I don't care that the sample speeds up or slows down when doing so.

What I am doing is say for example sake I have a wav file of a sine wave at 261.63hz (middle C) and this file has 1000 samples (example). Now I want to resample this sine wave so it plays at 146.83hz. I am dividing the desired frequency of 146.83hz by 440hz which gives me 0.3337. I then divide up the 1000 samples by 0.3337 to get the desired frequency. I will admit I am fairly new to all this but I think I am doing this correctly? This now leaves me with a smaller array of samples that when played gives me a rougher 146.83hz sine wave.

Now the issue is I am getting high frequency noise in the 2k to 20k range which varies on each note played. I think I am suppose to use interpolation to remedy this but I have not figured out how to implement it properly in my code. Do I resample then interpolate the resampled data OR do I resample and add in the interpolated samples between each sample?

I tried linear, cubic and hermite interpolation but neither made a difference in removing the high frequency noise. Though I don't believe I am using them properly. I got the interpolation the code from http://www.musicdsp.org so I'm  pretty sure it's correct.

Any help or suggestions would be great.


texts:
I am having trouble understanding how the antenna picks up signals that utilize FSK.

Because FSK alters the frequency and most antennas pick up signals on a very select frequency, how exactly does it pick up the signal?
Or does the antenna pick up on a range of frequencies. And if so, doesn't that limit the bandwidth?

I am having this question, because when one tunes a radio to a certain signal, it picks up signals of that frequencies only. But I believe radio stations use FSK to transmit their data also.


texts: Fourier Transform of Alternating Periodic Rectangular Pulse

texts:
I'm having trouble determining Fourier transform of signal. I have 2 ideas on how to solve this problem.
Given the signal is periodic I could use formula for Fourier transform of periodic signals:

$$X(j\omega) = \sum\limits_{k=-\infty}^{\infty} C_k \cdot 2\pi \delta(\omega-k\omega_{0})$$
where $\omega_0 = \frac{2 \pi}{4T}$.

Also, I could make signal $g(t)$ non periodic rectangular signal, and $h(t)$ alternating pulse train, and then $g(t)*h(t)$ should give me signal I needed.
Although I have some ideas I'm stuck solving this problem. Please help. Thanks


texts: detection of periodicities in n-dimensional signals

texts:
Generally speaking, what analyses are necessary and sufficient for the detection of periodicities in an n-dimensional signal amounting to a discretely sampled density distribution over n-dimensional space?

Assume that the signal-to-noise ratio and sampling frequency in all dimensions are adequate for detection of any periodicities that are present.

Would an n-dimensional Fourier transform and an n-dimensional autocorrelogram suffice?
Would these be limited to detecting periodicities that are aligned
with one of the basis vectors of the n-D space?

Suggested references
   for further reading would be appreciated.


texts:
What kind of bit error rate do they get from deep space communications (Pioneer, Voyager, et.al.), and what kind of modulation and FEC allows them to recover messages with that microscopic level of received signal power?

Are there more modern modulation methods and encoding schemes for similar channel conditions?


texts:
I am a hamradio hobbyist with little knowledge in SP.  I have to decimate the sampling rate of the received FM signal so I can get to the sampling rate that is supported by the sound card of my computer.

In the manual it says that I can choose the decimation number by dividing the initial sampling rate by bandwidth.Is this a general formula that I can always use? I tried google but could not find this specific formula.
Thank you


texts: When calculating SNR, is noise included in signal?

texts: Given a signal $x(t) = s(t) + n(t)$ where $s(t)$ is the desired signal voltage and $n(t)$ is the noise, should the signal to noise ratio of this signal be 20log(xrms/nrms) or 20log(srms/nrms)? i.e., should I use the desired signal or the total signal for the signal to noise ratio?


texts:
Can someone please explain the illustration (figure 1.19) at the bottom of this image? It looks like there are four sampling points but I don't understand what the different curves represent. Actually except from the dashed curve I can't even distinguish the different curves from each other. Can someone please clarify this so that the illustration makes sense?

(Image source: Digital Signal Processing Principles Algorithms and Applications)


texts: Verifying the result from lowpass filter

texts:
I have a event record from acceleration censer.
I'm trying to remove gravity from the event record by low/high pass filter as following site shows.
https://developer.android.com/reference/android/hardware/SensorEvent#values
As the first step I tried to extract Gravity by using following formula as a low pass filter.
G[t] = 0.1*G[t-1]+(1-0.1)*a[t]

G[t] is extracted gravity as the result.
G[t-1] is previous result .
a[t] is censer value in current row of event record.

I calculated it on Excel sheet for each X,Y,Z axis.
And I'd like to know how to verify the result.
I thought I can verify by summarizing  G value for X,Y,Z from low pass filter will be 1023 based on following picture.

Green is X
Blue is Y
Red is Gravity
G only weighs on Y at the first.
After rotate X-Y,G will weigh on both X and Y.
There will be G on X and G on Y
However, Gx+Gy+Gz is not 1023 during censing so far.
Here is  top 8 line of the data.
time-x  x   time y  y   time-z  z   time-zz zz
0   904 0   -20 0   504 0   -320
0.006   892 0.004   -16 0.006   504
0.008   888 0.015   -4  0.015   488
0.02    888 0.018   8   0.02    492
0.022   892 0.022   16  0.022   492
0.031   876 0.029   20  0.032   500
0.034   876 0.034   4   0.038   496
0.049   884 0.051   4   0.056   488


texts: Covariance matrix of an adaptive filter input

texts:
I run many times in equations containing the trace of covariance matrix of an adaptive filter input. But it is not really clear what it is.

For example in this paper the input covariance matrix is

$$\textbf{R} = E\left(\textbf{u}_i^*\textbf{u}_i \right)$$

where:

$\textbf{u}_i$ is row vector of input (seems to be regression vector in given time)
$^*$ denotes Hermitian conjugation (complex conjugation for scalars)
$E$ is not explained, but I would expect it is the expected value

Questions:

What does it mean?
How does the matrix look like, what is the shape of it?
How do I estimate it from input data $\textbf{U}$ - matrix where one row represents one input vector?


texts: Resolution of Discrete Fourier Transform is 1/T - Mathematical proof?

texts:
In many articles I see that the frequency resolution of the Discrete Fourier Transform (DFT) equals Fs/N where Fs is the sampling rate and N is the total number of samples. Fs/N is equivalent to 1/T where T is the total time the signal was sampled.

So the longer the time of sampling took place the better the resolution of the spectrum (allowing to separate different tones)

I'm looking for the mathematical proof of the 1/T resolution as I do not see it directly from the definition of the DFT.


texts: Calculating cutoff frequency for Butterworth filter

texts:
I have a problem while calculating cutoff frequency, suppose we have these specs.

Firstly, I calculated the order of the filter and got $N=5.8858$ and round it up to get $N=6$.

Now I'm supposed to get $\Omega_c$. Using these equations:

\begin{cases}1+\left(\frac{0.2\pi}{\Omega_c}\right)^{2N} = \left(\frac{1}{0.89125}\right)^2 \quad&\quad (1)\\
1+\left(\frac{0.3\pi}{\Omega_c}\right)^{2N} = \left(\frac{1}{0.17783}\right)^2\quad&\quad (2)
\end{cases}

Now with $N=6$ and $T=1$, substituting in $(1)$

\begin{align}
\left(\frac{0.2\pi}{\Omega_c}\right)^{12} &=\left(\frac{1}{0.89125}\right)^2 - 1 =0.25893\\
\implies (\Omega_c)^{12} &= \frac{{(0.2\pi)}^{12}}{0.25893}\\
 &= 40.29
\end{align}
But in the textbook it says $\Omega_c = 0.7032$, what I did wrong? Any help would be appreciated.


texts: Causal Noise Free 1st Order Derivative in Discrete Domain

texts:
I need to have causal noise free first order derivation (Derivative).
Now I am using a simple finite differences formula:

$$ \frac{ x \left( n \right) - x \left( n - 1 \right) }{ {T}_{s} } $$

The filter above, being High Pass, is super sensitive to noise.
Do we have a noise free or noise reduced first order differentiation which does not need more than three point for calculation?


texts: Resampling with and without replacement for estimating significance of spectral components

texts:
In order to test the significance of spectral components, it seems reasonable to randomly sort the data in order to destroy all the serial correlations / spectral order e.g. 100,000 times, and then estimate the significance or C.I.s by comparing your sample spectrum against the spectra of all the resorted data. Essentially H0 = there is no serial correlation / spectral order, and therefore the signal is white noise.

Here the resorted data would be resampling without replacement. My question is - is this correct or incorrect for spectral estimation? Would it be better to use resampling WITH replacement? And, if so, then why?

Thanks all for your responses! Indeed I am aware it is called bootstrapping - and just to clarify, my intention is not only to resort the time series 10,000 times, but to also create 10,000 spectra from these resorted versions of the time series in order to establish how many times my spectral component of interest is above the resorted time series spectra.

So H0 is that, over 10,000 comparisons, my spectral component is not higher in amplitude than the corresponding component in the resorted spectra more than 5% of the time. My H1 being that the spectral component is greater in amplitude (and therefore not due to chance) in at least 5% or more of the comparisons.

Indeed, I think you can go further than this and also calculate how many times it is BELOW the resorted time series spectra - this would, presumably, give you an indication of its persistence (i.e. it is more persistent than white noise and is therefore unlikely to change in a given period).

The code below calculates: the spectrum for a sine wave X[k] (length 30 time points), an average spectrum for the sine wave resampled 20000 times without replacement, an average spectrum for the sine wave resampled 20000 times with replacement. The original signal spectrum could as well be an average spectrum of, say, 200 signals of interest.

 Fs=1; %Sampling Frequency
 T=30; %Total signal duration in seconds
 t=0:(1./Fs):(T-(1./Fs)); %Time vector
 f1=1; % A frequency figure (in Hz)
 x=cos(f1*t);
 X=fft(x);
 X = abs(X);

 for i = 1:20000
     ps = round(rand(1,30)*29)+1;
     x_p_all(i,:)=(x(ps))';
     X_hat(i,:) = abs(fft(x_p_all(i,:)));
 end

 X_hat_mean = mean(X_hat);

 for i = 1: 20000
     ps2 = randperm(length(x));
     x_p_all2(i,:)=(x(ps2))';
     X_hat2(i,:) = abs(fft(x_p_all2(i,:)));;
 end

 X_hat_mean2 = mean(X_hat2)

 plot(X_hat_mean)
 hold on
 plot(X)
 plot(X_hat_mean2)

It can be seen that the resampling with replacement averages to be the same as the resampling without replacement average (as expected). Obviously, the energy of an individual resampled signal will be different for with and without replacement, but bootstrapping undertakes the process thousands of times.
Now, if we don't permute and take the average, we have the option to compare our spectrum against each individual permutation, where each comparison would be a comparison like this:

Clearly, in this instance, the spectrum of the resampled signal (in this case, without replacement) looks far from white noise because it contains spurious periodicities as a result of the short time duration and limited number of frequency bins. My method compares these resampled signal spectra against the original signal spectrum e.g. 10,000 times to establish significance - a similar method to bootstrapping (which does not, as far as I am aware, require Bonferroni correction).


texts: Shift operators in graph signal processing (GSP)

texts:
A key element of Graph Signal Processing (GSP) is the shift operator, that describes how a signal changes when shifted on the graph.
There exists many different proposals of shift operators. Several rely on various kinds of Laplacians, others on random walks. Some preserve the signal sum, others preserve its energy.
Also, most shift operators apply to vertex signals, but some apply to edge signals.
I am interested in building a comprehensive list of shift operators in GSP literature (for vertex or edge signals), and in understanding the pro and cons of each definition. I would love to classify these shifts operators in meaningful categories for comparison.
What are your favorite GSP shift operators, and why?
The following references provide many examples:

On the Shift Operator, Graph Frequency and Optimal Filtering in Graph Signal Processing by Adnan Gavili, Xiao-Ping Zhang
The Dual Graph Shift Operator: Identifying the Support of the Frequency Domain by Geert Leus, Santiago Segarra, Alejandro Ribeiro, Antonio G. Marques
Unitary Shift Operators on a Graph by Bruno Scalzo Dees, Ljubisa Stankovic, Milos Dakovic, Anthony G. Constantinides, Danilo P. Mandic
Translation on Graphs: An Isometric Shift Operator by Benjamin Girault, Paulo Gonçalves, Eric Fleury
Flow Smoothing and Denoising: Graph Signal Processing in the Edge-Space by Michael T. Schaub, Santiago Segarra

However, shift definitions are scattered over a quite wide and heterogeneous set of papers.


texts: Non Linear Recursive / Online Least Squares

texts:
I am searching for a recursive or online non linear least squares algorithm. I want to spread the computation out as new data is sampled like in the linear Recursive Least Squares or the LMS.

Ideally a recursive Levenberg–Marquardt algorithim would exist as Levenberg–Marquardt works great on my non linear problem but need to reuse all the samples to calculate a new estimate.

The online recursive non linear estimators I have found so far are:

The Extended Kalman Filter
Unscented Kalman Filter
Particle filter

Some articles I have found interesting:

Two step estimator as described in the paper "Optimal recursive iterative algorithm for discrete nonlinear least-squares estimation (1996)" by Gordon T. Haupt, N. J. Kasdin, George M. Keiser, and Bradford W. Parkinson
The Modified Extend Kalman Filter as described in the paper "A recursive algorithm for nonlinear least-squares problems (2007)" by A. Alessandri, M. Cuneo, S. Pagnan, M. Sanguineti

Do other online non linear estimators exist?
What are the popular online non linear estimators? I have minimal experience with non linear estimation.


texts: Z-domain transfer function to difference equation

texts:
So I have a transfer function $ H(Z) = \frac{Y(z)}{X(z)} = \frac{1 + z^{-1}}{2(1-z^{-1})}$. I need to write the difference equation of this transfer function so I can implement the filter in terms of LSI components. I think this is an IIR filter hence why I am struggling because I usually only deal with FIR filters. I have tried to simplify the filter, and I get:

$H(z) = \frac{z+1}{2(z-1)}$

This gives me the gain (K = 0.5) and the poles as +1, and the zeroes as -1, hence filter is stable.

Can anyone help me? I usually divide through by the denominator, and hence get the difference equation, but I can't in this case. Then it's just a  case of looking at the difference equation and implementing the filter with delays, multiples, adders etc.


texts: Digital 16-bit PCM audio signal scaling and STFT spectogram amplitude

texts:
I'm sampling a 16-bit PCM audio signal at 44100 Hz and I'm trying to understand the amplitude level of both the digital audio signal and its corresponding spectrum. The final goal is creating audio fingerprints for visualization and sound comparison. By now, I'm conscious I'm not calibrating my hardware so that dB scale is just a relative measurement to itself

First, for the 16-bit PCM audio signal I do understand the signal (Audio ADC lectures) ranges from -32768 to 32767, scaled to [-1 to 1) sig_scaled = sig_int16 / 32768 so that both -32768 and 32767 correspond to the maximum amplitude in db (0 db):

max_db = 20*np.log10(abs(+-1))= 0 dB

*Here, I do also understand it could be scaled to [0, 1) sig_scaled=(sig_int16-+32768)/65536

Meanwhile, the minimum one (~96db) corresponds to:

min_db = 20*np.log10(1/65536) ~= -96 dB

Which scaling is better [-1 to 1) or [0 to 1)?
How to avoid the "RuntimeWarning: divide by zero encountered in log10" when the ADC lecture or sig_int16 value is 0?

I think [-1,1) makes more sense since a sine wave goes from -1 to 1 but the obtained spectograms for both scales are similar.

Second, regarding the spectrum I've performed an stft on the [-1, 1) scaling, obtaining the following image (stft_db = 20*np.log10(stft/np.max(stft))):

What's the meaning of the -159.93 dB minimum and 0 dB maximum value? / Shouldn't the minimum value be -96 dB or bigger?
Which scaling ([-1, 1), [0, 1) or dB amplitud) should the STFT be performed on?

Thank you folks, let me know if you need more explanations or code!


texts: what is the output of BoW after an image has been trained with SIFT algorithm and k-means

texts: SIFT algorithm provides a 128 dimensional feature vector that is used for image classification.When all the interest points(key points) are taken together and K-means clustering is applied,the image is then clustered into many clusters.My question is how does this clustering/grouping of feature points together help in creating a dictionary(of equivalent image words).Does the output after applying this BoW algorithm yield a single vector for an image.


texts:
I need to find the Fourier transform of the following signal:

$$ne^{-an}u[n]$$

The answers start by using the rule of the basic signal:

$$a^nu[n] \rightarrow \frac{1}{1-ae^{-j\omega}} $$

and then uses the property of differentiation in the frequency domain. However I think there's something wrong with the way the rule is used. The answers assume that:

$$e^{-an}u[n] = a^nu[n] \rightarrow \frac{1}{1-ae^{-j\omega}} $$

Is there something wrong or am I missing something?

Here's the answer for reference:


texts:
Assuming we have time domain signal $$x(t)\,\,\,\,\,\,\, 0\leq t< T$$ and we decide to upsample the signal by factor $L\geq1$.

Then the samples of this oversampled time domain signal are represent as

$$x(kT/L)$$

I have read in paper that setting $L=1$ corresponding to the case of Nyquist rate sampling, does anyone know why this is case

If $L=4$ how many zeros are padded to obtained the upsampled signal?

I provided the text that made me ask the question in first place.


texts: Impulse response measurement : practical questions

texts:
To measure system impulse response I can do the following steps :

system input : it is excited with an excitation signal having a duration t_exc with a bandwidth bw_exc (for example sine linear sweep)
system output : it is measured during t_out with a sampling frequency f_s
transfer function H is obtained applying fft(out)/fft(in)
impulse response h is obtained applying ifft(H)

I have several questions about details of this method. I know that these questions really depends on system studied, but I would like to know the general guidelines :

How to determine the appropriate excitation signal duration?
Do you recommend excitation signal windowing?
How to determine appropriate output signal measurement duration?
before applying fft, how much zero padding do I need to add to excitation signal and output signal?
how to avoid zero division when excitation signal max frequency < f_s/2 ?


texts: What is the pitch salience? I see it a lot when looking for papers about multiple pitch detection but haven't found a proper definition of it.


texts: Augmentation for EEG signal classification using Deep Learning

texts: Augmentation is a technique that we use in deep learning for expanding the training dataset. It includes different ways of modifying an image and adding it to the training dataset. My question is if we are using augmentation for EEG classification and say we generate spectrograms and feed them as images to a CNN can we use Augmentation to increase the training dataset?


texts: Why do we not encode frequency in 3 dimensions in MRI?

texts:
If I understand correctly, to select a particular slice (of an object, during an MRI scan) that is orthogonal to the direction of the main magnetic field applied, we apply a magnetic field that varies linearly in strength in the direction of the main magnetic field (conventionally the z direction). Subsequently, in either the x or y directions within that slice we apply another linearly increasing magnetic field. This means that within that slice, protons (of the same molecules) precess at frequencies linearly varying in a chosen direction. This means that we can approximate the slice as a bunch of rectangles, each of which induces different frequency signals in the RF coil (although in fact no two points along the direction of the gradient have the same frequency). Now, instead of applying a third linearly varying B field in the next dimension (x or y, depending on the previously chosen one), magnetic pulses are applied in this dimension for varying periods of time so as to shuffle the phases of protons precessing at the same frequency due to the previous 2 magnetic fields. I understand that by fourier transforming the sum of all these guys (all the different signals recorded at different phase differences, for each frequency band, by means of tipping the magnetisation into the x-y plane) we can isolate particular voxels' contribution to the signal. But why, instead of doing this complicated phase pulsing, don't we just apply the third magnetic field in the left-over dimension and do business as usual?
To my mind this method sounds a lot simpler, and the only clarification I can find for phase encoding is that I am actually mistaken in my understanding of the entire process. Someone please explain!
I'll include a paragraph of what I'm talking about (I can't copy paste the text and I am definitely not transcribing it):


texts: Estimating scaling and translation between two consecutive temporal images

texts:
I have two camera images taken from a camera mounted over a vehicle at time $t$ and $t+1$. I need to estimate the scaling and translation in  the images in spatial domain itself so that information from two images can be fused to get a better image. I have camera calibration parameters, vehicle dynamics parameters for the calculation.

Is it possible to estimate scaling efficiently? For translation I have used cross-correlation peak based approach, although it is not correct all the times.


texts: MRI K-space to image: why track frequencies in two dimensions?

texts:
While trying to better understand how an MRI goes from k-space to an image, I came across this wonderful website that explains how you would represent an image as a collection of rows of pixels  where each row of pixels is represented as a sum of waves of greyscale intensity.
First you select a single row from an image.

Then you apply the Fourier transform to represent the row as a collection of waves.

What I don't understand is, if we can get the frequencies for each row and accurately reconstruct what a row of the image looks like via frequency encoding, and if we already know where each row is in relation to the other rows (position of the blue line), why do we need to do the same basic thing along a different dimension (the y-axis using phase encoding)? Isn't it sufficient to simply stack the rows to get the full image?


texts:
I can create a Linear Quadratic Gaussian Integral(LQGI) controller very easy by using GNU Octave. LQGI is in the area of Optimal Control Theory.

But there is something called Robust Control Theory. And every time I'm searching for that I got these in my  search field:

$H_2$
$H_{\infty}$

Questions:

What are they?
How do I apply them to control the processes?
Are they good for control systems?
Are they better that PID and LQGI?
Are they easy?
Are they for nonlinear system?
Can they have a Extended Kalman Filter?
What's the difference between them both?
Can you show me an easy example how to apply them?
Other?


texts: IIR filter relationship between number of poles, order and second order sections

texts:
I have ported Iowa Hills Code Kit to a mcu to calculate IIR filter coefficients. They are used to initialize the CMSIS Biquad Cascade IIR Filters Using a Direct Form II Transposed Structure.
Trying to compare the output of my mcu program with a similar python script written using SciPy I realized that there is a discrepancy between the two, as explained below.
SciPy
Calculating a 4th order bandpass IIR filter using the following in SciPy,
sos = scipy.signal.iirfilter(4,
                            [1500, 2500],
                            btype='bandpass',
                            analog=False,
                            ftype='bessel',
                            output='sos',
                            fs=50000)

results in 8 zeros and 8 poles which is then converted to 4 Second-Order-Sections. This gives

My understanding however, is that the code above is producing an 8th order filter and that a 4th order filter should be represented as 2 second order sections and should have 4 poles and zeros. I might be missing something very obvious, but I'm not sure why signal.iirfilter produces this result!
IowaHills
To calculate a 4th order bandpass IIR fitler using IowaHills however, I would have to input numPoles=2 to get 2 second order sections as shown below,
TIIRCoeff coeffs = {0};              // filter coefficients struct (IowaHills)
TIIRFilterParams fparams;            // filter parameters struct (IowaHills)
fparams.IIRPassType = iirBPF;        // bandpass filter
fparams.ProtoType = BESSEL;          // proto bessel
fparams.BW = 0.10;                   // filter width
fparams.dBGain = 1.0;                // gain
fparams.Gamma = 0.0;                 // transition bandwith
// NumPoles = [order / 2] for BandPass and BandStop
// NumPoles = [order]     for LowPass and HighPass
fparams.NumPoles = 2;                // 4th order
fparams.OmegaC = 0.41;               // centre frequency for bandpass
fparams.Ripple = 0.0;                // only for chebyshev/elliptical
fparams.StopBanddB = 0.0;            // only for chebyshev/elliptical
coeffs = CalcIIRFilterCoeff(fparams);

this gives

Can anyone explain why the relationship between number of poles and filter order is different in the two programs? I'm assuming both of them are correct so what is the right way to interpret the outputs?


texts:
I am trying to understand the basics of IQ modulation. The amplitude of a demodulated IQ signal is given as $\sqrt{I^2 + Q^2}$ which is easy to understand from the phasor diagram.

Yet, I don't understand the IQ power: ${I^2 + Q^2}$ . Where is the load impedance in this expression ? I'd expect something like

$Z_{load} = 50 \Omega $

$P = \frac{I^2 + Q^2}{Z_{load}}$

What is missing in my understanding ?


texts: Training own classifiers using trainCascadeObjectDetector for face detection

texts:
I want to train a cascade classifier to face detect. I read this page. But I wont use trainingImageLabeler app, because I already have face and nonface database sized 24x24. How I use trainCascadeObjectDetector function with my own database?

Is there any example except this you know?


texts:
I am new to OpenCV and DSP. I want to code a filter which act like a real optical filter. When you filter the image, the output looks like shotted with a $550nm$ optical filter lens. Like a Bandpass filter for image. I wanna filter all lights below $550nm $ wavelength.

Is it possible with OpenCV? How can Ii do that? Filter2D or Gabor Filter?

Sorry for my English and sorry if the question is irrevelent.

Thanks


texts:
here is my scenario: one speaker is sending a chirp, and a microphone is recording it as well as all following echoes. The function of the chirp is already know.

My goal: subtract the original chirp from the recorded sound, and leave only the echoes.

I think this questions is hard for me in these ways:

To match the original chirp and recorded sounds, I don't know which
approach to use: cross-correlation, or minimize sum(abs(difference
between two vectors)). The difference here means element-wise difference.
The echo may be mixed with the chirp itself, so I'd better only do
cross-correlation on the first small portion of the chirp.
The amplitude of recorded sound depends on the microphone.

Thank you very much!


texts: impulse response cascaded with time reversed of itself

texts:
Consider a filter with real-valued impulse response $h[n]$. The filter is cascaded with another filter whose impulse response is $h'[n] = h[-n]$, i.e. whose impulse response is the time-reversed version of $h[n]$.

$h[n]$ cascaded with $h[-n]$

i think we have to perform convolution

if so what is phase of new cascaded impulse response


texts: DFT of Transient vs. Analytical Solution

texts:
I'm trying to compute the Fourier Transform of a lightning impulse waveform.  The python code I'm using is below.  The problem that I am running into is that the magnitude of the DFT and the magnitude of the analytically derived Fourier Transform are off by a factor equal to 1 / time period, where the time period is the length of time that I am using to model the signal for the DFT.

from math import *
import numpy as np
import matplotlib.pyplot as plt

def Impulse(t, V_0=1.04, alpha=0.0146e6, beta=2.467e6, t_chop=float("inf")):
    if t<=t_chop:
        return V_0 * (exp(-alpha*t)-exp(-beta*t))
    else:
        return 0

npImpulse = np.vectorize(Impulse, excluded=['V_0', 'alpha', 'beta', 't_chop'])

# Sampling rate and time vector
total_time = 500e-6 # seconds
n_samples = 2**24
sampling_rate = n_samples/total_time # Hz (samples/s)

t_values = np.linspace(0, total_time, n_samples)

# Time constants
tau1 = 68.2e-6
tau2 = 0.405e-6

y_values = npImpulse(t_values, V_0=1.0, alpha=1/tau1, beta=1/tau2)

# plot with various axes scales
plt.figure(1)
plt.plot(t_values, y_values)
plt.grid(True)
plt.show()

# FFT calculation
sp = np.fft.rfft(y_values) # real FFT, + freq only
freq = np.fft.rfftfreq(y_values.size, t_values[1]-t_values[0])

# Normalization
# This is where things are not right.  Should be dividing by n_samples, not samples/sec
sp = sp / n_samples

# Calculation of Fourier Transform analytically
sp_analytic_real = tau1/(1+(2*pi*freq*tau1)**2)-tau2/(1+(2*pi*freq*tau2)**2)
sp_analytic_imag = -(2*pi*freq*(tau1**2)/(1+(2*pi*freq*tau1)**2)-2*pi*freq*(tau2**2)/(1+(2*pi*freq*tau2)**2))
sp_analytic_mag = np.sqrt(sp_analytic_real**2 + sp_analytic_imag**2)

plt.figure(2)
plt.plot(freq, np.abs(sp), label='FFT')
plt.plot(freq, sp_analytic_mag, label='Analytic')
plt.xscale("log")
plt.yscale("log")
plt.legend()
plt.grid(True)
plt.show()

If I set my calculation time period to 1s, the two curves match.  If I divide the magnitude of the DFT return by rfft by the sampling rate (or alternately multiply by the 1 / time period), they match.  I can't for the life of me figure out why that should be.  Can someone clue me in to the obvious that I am missing?


texts: I have the scatter plots for a Hamm and Tukey(0.25) windowed data. The complex data plotted from Hamm windowed forms a pattern, why?

texts: The Hamm windowed complex data forms some sort of a quantized pattern whereas in the same code if tukeywin(0.25) is applied normal noise is observed. Any idea why ?


texts: Implementation of Block Orthogonal Matching Pursuit (BOMP) Algorithm - Fix Given Code

texts:
This is my implementation which doesn't work:

function [result,gamma,T]=block_OMP4(phi,signal,N,d,maxIter)
[m n]=size(phi);
r=signal;
gamma=[];
T=[];
i=1;
tol=1e-4;
err=1;
x=zeros(N*d,1);
val=[];
figure,
while  (i<=maxIter && err>tol)

    maxv=[];
    for j=1:N
        startx=(j-1)*d+1;
        endx=startx+d-1;
        block=phi(:,startx:endx);
        t=block'*r;
        maxv(j)=norm(t,2);
%
    end
    [maxvl,maxvindx]=max(maxv);
    gamma=union(gamma,maxvindx);
    x2=lsqr(phi,r);
    xb=zeros(N*d,1);
    for i=1:length(gamma)
        stratx=(gamma(i)-1)*d;
        endx=stratx+d-1;
        xb(stratx:endx)=x2(stratx:endx);
        x(stratx:endx)=x(stratx:endx)+xb(stratx:endx);
    end

    y=phi*x;
    r=signal-y;

    val(i)=norm(r);
    err=norm(r)
    i=i+1;

end
result=x;
end

Anyone has a working implementation?


texts: GNU Radio - General Filtering Questions - Any feedback is appreciated!

texts:
This is a fairly general question I feel like.

So, the scenario I am presenting is this:

I am generating random data using the random source module and feeding that data into the Constellation Modulator, thus producing a Basebase modulated signal. This takes us to my first question, why is the signal from (-f to +f) and not just (0 to +f)? Attached is a picture.

Next question is this, if I implement a low pass filter from (0 to f) f being a frequency relatively close to 0 Hz, why do I not see the rejection following with the roll off? Also, why is the negative frequency component still there?

I am teaching myself quite a bit of dsp and gnuradio, so any constructive feedback would be great!

If you look at the picture that actually has some decent pulse shaping, it seems to me that the cutoff frequency is far to high as well, but it's only at this range that I get decent pulse shaping.

Thanks!

-DasBoy


texts: Why do we get ripples in the frequency spectrum of a chirp, and how do I remove it?

texts:
The ripples in the amplitude spectrum of a linear spectrum is undesirable. I desire to obtain a flat spectrum response. Hence, how can I remove these ripples? Also, which factors (like time bandwidth product, chirp rate) are affecting the ripples, and how?

I designed a chirp with a frequency from 10 to 90 Hz in MATLAB. However, I obtained Fresnel ripples (which normally occurs) in the passband, along with overshoots at the band edges. I need to study the application of chirps on a material and use it for defect analysis. For this, I need to optimize the whole energy applied to the chirp. Hence, I want to get rid of these ripples. So, what parameters like chirp rate and time of sweep will affect the system?


texts: Synthesis discrete time signal from fourier coefficients

texts:
Following information is given about a signal $x[n]$

$x[n]$ is real and even signal
$x[n]$ has a period $N=10$ and Fourier coefficients $a_k$
$a_{11}=5$
$\frac1 {10}\sum_{n=0}^9 |x[n]|^2=50$

How can we obtain $x[n]$ from these information?
I know the first info makes the Fourier coefficients real and even, third info makes the Fourier coefficient $a_1=a_{11}=5$(since Fourier coeffs are also periodic with N=10) and fourth info is Parsevals theorem.


texts:
I have a question:
What is the difference between
$ I \cos(t) - Q \sin(t) $    and   $ I \cos(t) + Q \sin(t) $

except that they are phase shifted 90 degrees ?
On the internet I can find both modulators with adding quadrature data and substracting quadrature data.

http://www.ni.com/tutorial/4805/en/  at the end
http://zone.ni.com/reference/en-XX/help/372058G-01/nirfsa/iq_modulation/

thank you


texts: How can a stable system have a RHP pole?

texts:
I am working with a control system which has an unstable pole in the process. So I have the following transfer functions for the controller and for the process:

G = 10/((s+10)*(s-1));
K1 = K*(s+1)/s;

where G is the process and K1 is the controller. In the controller I have to choose a gain K such that the system is stable, and using the Routh Criterion I have found that it is stable for $K>1.125$, so I have chosen has value 4 because it gives me good performances. Now, I have that if I do the root locus I have:

L1 = K1*G;
figure;
rlocus(L1)

so I have a pole in the RHP. But I have that with an appropriate value of the gain of the controller we have a stable system. But in this case the pole should be in the LHP, not in the RHP.

So, what I don't understand is how it is possible that the system is stable even if it is present a pole in the RHP?

Can soebody help me?


texts:
I have a signal generated from an accelerometer, what I want to do is to receive the incoming signal and check if it is in a certain interval, if it is then I can apply interpolation on it(but interpolation is not the main issue). The basic problem is how to divide a signal into chunks Online and by online I mean that I don't already have an upper time limit of the signal.

A graph to describe further:

Green Signals are the actual samples I am receiving. Let suppose I have received first 2 Green Bars then I'll wait for  3rd green point and use interpolation to find a point at 10ms after that the 4th point will arrive and I want to use it to improve my interpolation by including 3rd and also 4th point until I reach 20ms where I will use the same method used on 10 ms point.

Until now what I have done is divided whole time sequence into chunks and performed interpolation but that is not very practical as in practice I won't have complete sequence already. I'm using python but it is not important. Any suggestion in this regard is highly appreciated.

If this question is in the wrong place please let me know which platform should I use.


texts:
I’m a newbie and wondering if I calculate the power spectrum from the peak amplitude or from the amplitude rms ?

This is my matlab script so far:

%-----------------------
% Signal time domain
%-----------------------

Fs = 2000; % Sampling frequency Hz
dT = 1/Fs; % Sampling period
L = 1024; % Length of signal [points]
t = (0:L-1)*dT; % Time vector

% Form a signal containing a 50 Hz sinusoid of amplitude 0.7 and a 120 Hz sinusoid of amplitude 1.

S = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t);
% convert row vector to column vector
S = S(:);
%---------------------------
% Windowing
%---------------------------
win = hamming(L); % Window
winS = S.*win; % Windowed signal

% Calc window scaling factor
coherentGain = sum(win)/L;
scaleFactorWin = 1/coherentGain;
%---------------------------
% FFT
%---------------------------
Y = fft(winS);
Y = scaleFactorWin*Y;
% Compute the two-sided spectrum P2. Then compute the single-sided spectrum P1 based on P2 and the even-valued signal length L.
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
Ampl = P1;
Ampl(2:end-1) = 2*P1(2:end-1); % Mult by 2 because we discarded the left side of the spectrum (skip first (DC) and last (Ny) values)
%---------------------------
% Power Spectrum : Magnitude^2
%---------------------------
% Calculate with already scaled windowed data
PS = 2 * P1.^2;     % Mult by 2 for discarded left side

% or with rms amplitude spectrum:

PSrms = 2 * (P1.*1/sqrt(2)).^2;

% PLot
f = Fs*(0:(L/2))/L;

plot(f,PS)
title('Single-Sided Power Spectrum of X(t)')
xlabel('f (Hz)')
ylabel('P²')

plot(f,PSrms)
title('Single-Sided Power Spectrum of X(t)')
xlabel('f (Hz)')
ylabel('P²')

This is my matlab script so far.
Im not sure how to calculate the power spectrum, either with PS or Psrms (see below)

Thanks for your help.
Mike

%-----------------------
% Signal time domain
%-----------------------

Fs = 2000; % Sampling frequency Hz
dT = 1/Fs; % Sampling period
L = 1024; % Length of signal [points]
t = (0:L-1)*dT; % Time vector

% Form a signal containing a 50 Hz sinusoid of amplitude 0.7 and a 120 Hz sinusoid of amplitude 1.

S = 0.7*sin(2*pi*50*t) + sin(2*pi*120*t);
% convert row vector to column vector
S = S(:);
%---------------------------
% Windowing
%---------------------------
win = hamming(L); % Window
winS = S.*win; % Windowed signal

% Calc window scaling factor
coherentGain = sum(win)/L;
scaleFactorWin = 1/coherentGain;
%---------------------------
% FFT
%---------------------------
Y = fft(winS);
Y = scaleFactorWin*Y;
% Compute the two-sided spectrum P2. Then compute the single-sided spectrum P1 based on P2 and the even-valued signal length L.
P2 = abs(Y/L);
P1 = P2(1:L/2+1);
Ampl = P1;
Ampl(2:end-1) = 2*P1(2:end-1); % Mult by 2 because we discarded the left side of the spectrum (skip first (DC) and last (Ny) values)
%---------------------------
% Power Spectrum : Magnitude^2
%---------------------------
% Calculate with already scaled windowed data
PS = 2 * P1.^2;     % Mult by 2 for discarded left side

% or with rms amplitude spectrum:

PSrms = 2 * (P1.*1/sqrt(2)).^2;

% PLot
f = Fs*(0:(L/2))/L;

plot(f,PS)
title('Single-Sided Power Spectrum of X(t)')
xlabel('f (Hz)')
ylabel('P²')

plot(f,PSrms)
title('Single-Sided Power Spectrum of X(t)')
xlabel('f (Hz)')
ylabel('P²')


texts: Should I increase the mean or the variance of Gaussian noise?

texts:
I am new in DSP.

I am currently working with a project and I want to try which algorithm works better given different image noises.

To test it out, I have been using the Gaussian noise in images.

I have been using this equation for imnoise:

f = imnoise(im,'Gaussian',0,0.001)

And then increasing by 0.001 and testing. As can be seen, I am using zero-mean noise. Is this the right approach to test the algorithm on different noises by
changing the variance, or should I just change the mean?


texts: Edge detection in spatial or in frequency domain?

texts: I want to ask why everywhere is used to do edge detection in spatial domain and not in frequency domain with FFT of image and apply HP filter. Thank you.


texts: estimate the impulse response after equalization

texts:
I want to estimate the impulse response of the channel at the receiver. Assuming
some arbitrary impulse response: h=[1  0.2 -0.4 0.0 0.6]. Once the equalizer is constructed, I get the equalizer weights in  w1.

Question1) How to get back the channel estimates: Say,

x is the input to the channel,

d is the output of the channel which is the input to the equalizer,

w is the equalizer,

y is the output of the equalizer or is the equalized signal.

Ideally, y and x must be same. Then would the channel estimates be the least square solution i.e., inv(w^Tw)*w*d' or  inv(w)?

Question2) Can I apply the LMS method for any source input --- Gaussian and non-Gaussian?

Below is the code for LMS equalizer

N=1000;
% number of samples

np = 0.01;
% noise power is 0.01
sp = 1;
% signal power is 1 which implies SNR = 20dB
h=[1  0.2 -0.4 0.0 0.6];  %unknown   impulse  response

x = sqrt(sp).*randn(1,N);
d = conv(x,h);
d = d(1:N) + sqrt(np).*randn(1,N);
w0(1) = 0;
% initial filter weights are 0
w1(1) = 0;
mu = 0.005;
% step size is fixed at 0.005
y(1) = w0(1)*x(1);
% iteration at “n=0”
e(1) = d(1) - y(1);
% separate because “x(0)” is not defined
w0(2) = w0(1) + 2*mu*e(1)*x(1);
w1(2) = w1(1);

for  n=2:N

    y(n) = w0(n)*x(n) + w1(n)*x(n-1);
    e(n) = d(n) - y(n);
    w0(n+1) = w0(n) + 2*mu*e(n)*x(n);
    w1(n+1) = w1(n) + 2*mu*e(n)*x(n-1);
end


texts: why is the DFS of a delta function equal to 1

texts:
I have a x[n] = $\delta$[n].
By formula is should be

$$

X[k]= \sum_{n=0}^{N-1} \delta[n]W_N^{kn}

X[k]= \sum_{n=0}^{N-1} e^{-j2*pi*kn/N}

$$

The formulae isn't showing for some reason. I took a screenshot of what I got here: http://imgur.com/6j0Ibgu

basically why is the summation of an exponential term going to 1. $W_N=e^{-j*2*pi*k*n/N}$ in this case. I tried to prove it via $\sum\alpha^k=\frac{1-\alpha^N}{1-\alpha}$ but it doesn't work for me.


texts: harmonic algorithm for extending compressed bandwidth of audio file

texts:
I just found this site and feel at home :)

I've been playin around with python and some compression ideas and was wondering if there's anything like what I am trying to figure out...

I've noticed that mp4's uploaded to instagram go from ~0 to 13100Hz and on Youtube from ~0 to ~14000Hz... which is OK, but what if there was an algorithm to extend or at least soften this highcut, by taking the amplitudes of the frequencies within the respective bandwidths and applying a harmonic function, (possibly with averaging/weights), to add in implied values to the currently non-existant values (i.e. the and additional bandwidth from 14,000 to 20,000Hz would be harmonically synthesized and still have realistic values.) It could be a plug-in or something.

Anyways, I'm going to check out the board some, and hope to hear some smart replies! thanks!


texts:
I have taken an image in height 2 meters of the ground as the original image and I want to simulate another image in 50 cm's of the ground, such as the simulated image.

The focal length of my camera is 3.7 mm (f = 3.7).

I have used from planar homography: H = k'*(eye(3) - ((t*n')/d))*inv(k) for each plane in original image in MATLAB where in the camera coordinate system:

x axis is perpendicular to buildings in the right side of the ground,
y axis is parallel to the ground,
and hence z axis is upward perpendicular to the ground.

The origin is camera location in original image and k is 3x3 camera calibration matrix:

[f 0 0;
  0 f 0;
  0 0 1]

t is the transfer matrix [0 0 -1.5], n is the normal vector of a plane, d is the plane distance to the origin, and eye(3) is camera rotation.

For example, the plane of the ground has n=[0 0 1] and d = 2, is it true?

Please help me with implementation of other planes (sky, buildings) in MATLAB.

Original image:

Simulated image:


texts: Minimum Sampling Rate of Bandpass signal

texts:
Consider a DSBSC signal s(t) = m(t).cos(2.pi.20k.t) where say m(t) is multi-tone signal with frequencies less than 5 kHz.

The spectrum gives a typical bandpass say triangle centred around 20k and -20k.
Now, when I calculate the minimum sampling rate;

Since the highest frequency is not an integer multiple of bandwidth of the bandpass signal; I end up getting the sampling rate as 25 KHz.

But when I draw the spectrum for the sampling rate of 10 kHz, the message signal centred at origin was doubled in amplitude and there was no aliasing as well.

Does that mean the minimum sampling rate is indeed 10 kHz?

I have used the typical bandpass sampling rate procedure many times, this was the first time it failed in giving me the minimum sampling rate!

edit: I doubt that the 10 kHz would be the minimum sampling rate to recover the message signal while the bandpass sampling rate of 25 kHz is the one which recovers the bandpass spectrum!

Also, if the above conclusion is correct, is it independent of the carrier frequency (Fc>>Fm)?


texts: Varying speed of shaft Spectral Kurtosis

texts:
I want to perform Spectral Kurtosis (Kurtogram) on 10-second vibration signal of 25.6 kHz sample rate, where the shaft is varying in speed.

Will a preprocessing step of resampling the signal be necessary (upsamling; so that there is an equal amount of samples within each shaft rotation), or will the SK handle the non-stationary behaviour?

I'm using Matlabs implementation:

  "The fast kurtogram algorithm uses bandpass filtering along with a
  simplified computation to approximate the spectral kurtosis for each
  window size and frequency ..

Update

(The kurtogram looks much better on the resampled signal, but I am not sure if this is the correct approach.)

The objective is to diagnose the health of a bearing by first applying an optimal band-pass filter based on SK. After that, I wish to find the enveloped signal. Then, I want to decompose the filtered signal using EEMD to look for fault characteristics. The bearing is on an operating wind turbine.

Much like this paper: Faulty bearing signal recovery from large noise using a hybrid method on spectral kurtosis and ensemble empirical mode decomposition


texts: How to calculate pitch fundamental frequency f( 0) ) in time domain?

texts:
I am new in DSP, trying to calculate fundamental frequency ( f(0) ) for each segmented frame of the audio file. The methods of F0 estimation can be divided into three categories:

based on temporal dynamics of the signal time-domain;
based on the frequency structure frequency-domain, and
hybrid methods.

Most of the examples are estimating fundamental frequency based on the frequency structure frequency-domain, I am looking for based on temporal dynamics of the signal time-domain.

This article provides some information but I am still not clear how to calculate it in the time domain?

https://gist.github.com/endolith/255291

This is the code, I have found, used so far :

def freq_from_autocorr(sig, fs):
    """
    Estimate frequency using autocorrelation
    """
    # Calculate autocorrelation and throw away the negative lags
    corr = correlate(sig, sig, mode='full')
    corr = corr[len(corr)//2:]

    # Find the first low point
    d = diff(corr)
    start = nonzero(d > 0)[0][0]

    # Find the next peak after the low point (other than 0 lag).  This bit is
    # not reliable for long signals, due to the desired peak occurring between
    # samples, and other peaks appearing higher.
    # Should use a weighting function to de-emphasize the peaks at longer lags.
    peak = argmax(corr[start:]) + start
    px, py = parabolic(corr, peak)

    return fs / px

How to estimate in time domain?

Thanks in advance!


texts:
Hi I have a FFT which is quite noisy. How to apply to my code Hamming window to make it less noisy. Look at my code:

# Loop for FFT data
for dataset in [fft1]:
    dataset = np.asarray(dataset)
    psd = np.abs(np.fft.fft(dataset))
    freq = np.fft.fftfreq(dataset.size, float(300)/dataset.size)
    plt.semilogy(freq[freq>0], psd[freq>0]/dataset.size**2, color='r')

for dataset2 in [fft2]:
    dataset2 = np.asarray(dataset2)
    psd2 = np.abs(np.fft.fft(dataset2))
    freq2 = np.fft.fftfreq(dataset2.size, float(300)/dataset2.size)
    plt.semilogy(freq2[freq2>0], psd2[freq2>0]/dataset2.size**2, color='b')

I put an answer in the questions.


texts:
I'm attending this course (Coursera: Audio Signal Processing for Music Applications) in which the professor derives a general equation for Discrete Fourier Transform (DFT) for a complex sinusoid. The following is the screenshot of the slide he used:

While the derivation is totally fine, I'm having trouble understanding the concluding statement

  if $k \neq k_0, denominator \neq 0$ & $numerator = 0$
  thus $ X_1[k]=N$ for $k=k_0$ & $X_1[k]=0$ for $k \neq k_0$

I did try value substitution and all, but all that can't seem to justify the statement. My understanding is that if $k=k_0$, both the numerator and the denominator would be zero and there would be no way the result is $N$ and I have no clue about the first part of the statement either. Or I'm missing something here (or forgotten some basic school math here). What's going on here?


texts: Find cepstrum of a function with dirac function

texts:
$$x[n] = d[n] + a\cdot d[n − N_p]\ ,\ \ |a| < 1$$

I need the complex and real cepstrum of this function as well as what would happen for $N=N_p/6$ for the DFT approximation $$x_p[n] = \frac{1}{N}  \sum_{n=0}^{N-1}\log(X[k]) e^{j2kn/N}$$.

It is a school assignment in which I am missing a lo of the functions required to proceed so  also any info you can link is greatly appreciated.


texts:
Image registration algorithms are usually based on point features such as SIFT (Scale-Invariant Feature Transform).

I saw some references to line features, but I was wondering if it would be possible to match image segments instead of points. For example, given source and transformed image:

I can do edge detection, blurring and Watershed Transform on each:

Regrettably, the segmentation turned out to be too different on each image to match individual segments.

I saw some papers on matching shapes and shape descriptors which are invariant to affine transforms, so this area seems to be promising...

Are there any segmentation methods more robust to affine (or even projective) deformations of the image?


texts: What does opencv `mean` function do when passing a mask of zeros?

texts:
INTRO

I am looking at this project to implement some changes.

The project is written in C++ using openCV library. Its goal is to detect cars in parking area to calculate the number of free spots. It does this in real-time (takes the input of a camera and processes it frame by frame).

The detection is performed on areas of the image that are rectangles whose vertices are loaded from a file. For each park (i.e. rectangle) to be monitored there is an instance of the class Parking defined in the files Parking.*** and all of them are collected in a vector (named parking_data in the following lines of code)

CODE

As far as I understood the processing starts on line 83 of the file main.cpp. The interesting lines are these:

cv::cvtColor(frame, frame_gray, cv::COLOR_BGR2GRAY);
cv::GaussianBlur(frame_gray, frame_blur, blur_kernel, 3, 3);
//...
for (Parking& park : parking_data)
{
  // Check if parking is occupied
  roi = frame_blur(park.getBoundingRect());
  cv::Laplacian(roi, laplacian, CV_64F);
  delta = cv::mean(cv::abs(laplacian), park.getMask());
  park.setStatus( delta[0] < atof(ConfigLoad::options["PARK_LAPLACIAN_TH"].c_str()) );

cvtColor converts the frame to grayscale then a gaussian blurring filter is applied. Then it iterates to the elements of the vector parking_data containing the instances of the class Parking which represent the area in the image corresponding to the parking spots.

getBoundingRect should return the vertices of the rectangle so to get a subimage of the blurred frame. Then a Laplacian filter is applied.

QUESTION

The line I do not understand is delta = cv::mean(cv::abs(laplacian), park.getMask());.

I understand it is calcualating a mean of the absolute value of the laplacian subimage. The mask is defined only in line 49 and used in the followings of the file Parking.cpp as

mask = cv::Mat::zeros(bounding_rect.size(), CV_8UC1);
for (cv::Point p : contours_points.at(0))
{
  polygon_points_in_bounding_rect.push_back(cv::Point(p.x - bounding_rect.x, p.y - bounding_rect.y));
}
vector<vector<cv::Point>> contours;
contours.push_back(polygon_points_in_bounding_rect);
cv::drawContours(mask, contours, -1, cv::Scalar(255), CV_FILLED);

so it is a matrix of 0s.
Then the value is compared to a given threshold.

OpenCV on the function says about the mean function.

When all the mask elements are 0's, the functions return Scalar::all(0)

So what is the mean function calculating when passing a mask of zeros?


texts: How good is the Leech lattice for AWGN coding or MSE quantization of Gaussian sources?

texts: Anyone know where to find graphs showing how well it works vs other schemes? thanks


texts: Reconstruct DWT for each cD1,cD2,cD3 and cA3 signals

texts:
This question must be basic for this forum, but I'm only start working with DWT recently, and I was working with CWT before.
I think I can decompose my signal to L level and then get the scalogram with the following scripts.
 my question is how can I reconstruct the signal for each details{cD1,cD2,cD3} and approximation CA3?

Here is what I have tried so far!

t=0:0.1:300;
f=[0.4 0.6 1];
y=sin(2*pi*f(2)*t)+sin(2*pi*f(2)*t)+sin(2*pi*f(3)*t);
len=length(y);level=10;
[c,l]=wavedec(y,level,'coif5');
cfd=zeros(level,len);
for k=1:level
    d=detcoef(c,l,k);
    d=d(ones(1,2^k),:);
    cfd(k,:)=wkeep(d(:)',len);
end
image(flipud(wcodemat(cfd,255,'row')));

P.S: Actually I want to estimate the frequency components in a signal, so I thought I must decompose my signal to few levels then take FFT from the details to estimate the components that involve is the main signal. Is that right? or a better idea?
for the CWT case I was estimating the frequency components directly.


texts: Confusion with regards to units in Sinc pulse

texts:
I am currently reading a section in the Handbook for MRI Pulse Sequence book where it says that the time-bandwidth product for a SINC pulse is given by:

$$
T \Delta f = Z
$$

where $T$ is the length of the pulse, $\Delta f$ is the bandwidth and $Z$ is the number of zero crossings of the SINC pulse. Now, my question might be very naive but it says that this is a dimensionless quantity but I do not understand why that will be:

For example, if the bandwidth is described in radians/sec, does it not have units of radians? In that case, how does it relate to the number of zero crossings, which is indeed a discrete number.


texts: Enhancement by adaptive wavelet shrinkage denoising

texts: I was wondering does anyone know any programming library or a way that can help me implement this algorithm.


texts:
I've been reading about Gabor filters, and there is a concept of scale. While I understand that Gabor features help in analyzing both the frequency and orientations of an image, what is the idea of scale? It appears that the original image is resized and then the filter is applied, for each scale.
Some elaboration on this would be much helpful.


texts: Visual regression test - Extract elements of GUI-Screenshot - MSER OpenCV

texts:
For a visual regression test I need to compare screenshots of webpages (different release-versions). I started with pixel by pixel compare. Actually I split the screenshot in different parts (maybe 10x10 rectangles) and compare the histograms of correspondending parts.

Now I try to find an approach to find regions of interest (with coordinates and/or bounding box) automatically (tried edge detection and so on) like this:

Example of image segmentation

I found different articles on the web dealing with blob-detection. One method for that is Maximally stable extremal regions (MSER) which is provided in OpenCV.

Feature Detection and Description

Now my questions:

Is MSER applicable for the given problem
Are there simpler approaches then MSER
Can someone give me some introductions to use MSER from OpenCV in Java (unfortunately I haven´t experience with OpenCV or JavaCV. Anyway I installed it on my machine.)


texts: Lost power while using FFT with apodizing mask

texts:
Previously I asked this question about FFT artifacts.

The solution that is most easily implemented is using an apodizing mask on the original image, tapering the borders to zero making the image left-right and top-bottom continuous.

The rest of the question is best illustrated with a series of pictures. First, we have the original image:

Now, we apply the cosine-tapered apodizing mask:

Note that we use this mask to do operations in Fourier Space on the second image; to get the power spectrum we will have to use a similar mask again but this is applied to both images so that is not the issue here. (i find it difficult to phrase this properly; please place a comment if you want me to elaborate).

So, now we use a FFT on the second image. In Fourier space we may do things with the image, but to trace our steps and see if the machinery works, we leave it be and just use the inverse-FFT back to real space.

Back in real space, we take the part of the image that was not tapered (so, the inner square) and obtain:

Now, for both the first and the third image I now want to find the power spectrum. In essence, this comes down to (1) using an apodizing mask again, (2) transforming to fourier space once more, and (3) taking the azimuthal average as a function of radius.

As a last subtelty, the power spectrum of the third image is multiplied by a factor $(N_1 / N_3)^2$, where $N_1$ is the number of pixels per side of the first image and $N_3$ the number of pixels per side of the third image. This factor is there to correct for the normalization $1/N$ in the definition of the FFT in numpy.

Doing all this, we obtain the following power spectra:

Clearly, the power spectrum of the third image (in red) has lost power with respect to the power spectrum of the original, first image (in blue). What is the cause of this and can it be compensated for? Or will I have to rely on more elaborate solutions than an apodizing mask (as given here) to not lose this power?


texts: I am building a morse code decoder in Python.  I am using SciPy and Numpy for FFT, filter, etc.  It is coming along ok, but the signal is noisy.  I am extracting the dominant frequency from the FFT of the input signal (which is audio from a wav file), and sampling at a rate great enough to detect whether a signal is present at these intervals.  Would I benefit from implementing a filter around this dominant frequency before processing the signal?  If so, what sort of filter, how wide, etc?  The dominant frequency is variable, but I determine it at each invocation programmatically (fft).


texts: How can we know the geometry of a blob detected by OpenCV

texts:
I am performing blob analysis using connectedComponentsWithStats, which returns elementary statistics about the blobs (AABB and area). It is also possible to compute geometric moments.

Anyway, I need to compute my own features and for this I need the detailed geometry, be it a chain of outline pixels, a run-lengh-coded description or similar.

I know that the contours can be obtained by other means (findContours), but I would like to know if the connected components engine can help me achieve my goal or if it is a black box forever.


texts: Which step response matches the system transfer function

texts:
A system has the following open loop bode plot: -

Which one of the plots below describe the closed loop step response for the entire system?

My attempt
My initial thought was to look at the static gain from the open loop bode plot, which is $-4\text{dB} = 0.63$. I interpret this as this will be the final value of the step response when it is settled. The only graph that sort of matches this description is graph 3.
However, it turns out that this is not the correct answer. And therefore my reasoning about the static gain can not be correct either. So my question is, how can I identify the closed loop step response from the open loop bode plot?


texts: How Many Samples Are Identical In Linear and Circular Convolution Result?

texts:
Suppose we have two sequences $x[n]$ of length $L$ , and $y[n]$ of length $P$. If we are doing a $L$ length circular convolution between the two sequences ($L>P$), how many of the resultant samples will be uncorrupted, i.e. same as linear convolution values ?

I think by intuition it should be $P-1$, but I'm not sure.


texts: Getting frequency content at different times from discrete wavelet transform coeffs

texts:
After being away from DSP for a long time, I am trying to familiarize myself with wavelet transform. Here is what I (think) have understood so far:

Wavelet transform provides you high time resolution at higher frequencies and high frequency resolution at lower frequencies.
DWT can be calculated by using QMF pair and subsampling. When used recursively, filter pair increases frequency resolution and subsampling decreases time resolution.
Result of level 1 DWT [cA, cD] = dwt(signal, Lo_D, Hi_D) essentially gives low frequency and high frequency splits of the signal called approximation and details.

If my understanding is correct, in this algorithm, the end result is essentially in time domain. If I want to know what frequencies are present, I have to take an fft of the coefficients. Is that correct? If so, is there a way I can get the frequency domain representation without extra step of taking fft?

Question 2: If scalogram is the answer to question 1, how is it generated using results of repeated use of dwt?
I want to avoid using ready MATLAB functions for better understanding except maybe dwt

Thanks


texts: 2D Image Convolution: Spatial Domain vs. Frequency Domain Convolution in the Computational Complexity Sense

texts:
I want to convolve 2 images I1 and I2. Let's assume they have the same size.
I know 2 methods:
1) Use direct definition of 2D convolution
2)  i. Computer FFT of I1 (call it I1_FFT) and FFT of I2 (call it I2_FFT)
   ii. Multiply I1_FFT and I2_FFT (call it I1I2_FFT)
  iii. Compute inverse FFT of result I1I2_FFT

I want to know the number of multiplications and additions for each case.  Note I don't want the Big O type of complexity.

It would be great if the answer includes a derivation.

In the end, I want to know for a specific NxM sized image, when I should use method 1 or 2.


texts: Solving Convex Optimization Problem Used for High Quality Denoising

texts:
The highest voted answer to this question suggests that to denoise a signal while preserving sharp transitions one should

  minimize the objective function:

  $$ |x-y|^2 + b|f(y)| $$

  where $x$ is the noisy signal, $y$ is the denoised signal, $b$ is the
  regularziation parameter, and $|f(y)|$ is some L1 norm penalty.
  Denoising is accomplished by finding the solution $y$ to this
  optimization problem, and $b$ depends on the noise level.

However, there is no indication of how one might accomplish that in practice as that is a problem in a very high dimensional space, especially if the signal is e.g. 10 million samples long. In practice how is this sort of problem solved computationally for large signals?


texts: Estimate/Compute parameters of image dataset

texts:
I have a dataset of images taken by Kinect (RGB and D)

How can I extract the following information only from the images:

Frame rate
Camera height

The dataset containes images of person (s) moving in the scene like here. So I do have the bounding box or/and the skeleton of the person.

How can I estimate them, is there any way to do it?


texts: What's the difference between a DSP and microcontroller?

texts:
What are the core attributes that make a DSP and microcontroller different from each other?

I have in mind a machine vision device that needs to run image processing algorithms but also control some networking and hardware I/O. How good would the DSP be at controlling hardware output pins, and how well can microcontrollers run image processing algorithms?


texts: Recognize irregular repetition of unknown feature (stamps at random locations)

texts:
I have an image which is generated by convolution between (1) a simple small patch and (2) larger image composed of localized point-like features (e.g. slat-and-pepper noise, sum of shifted $\delta$-functions).

In other words the image seems like the patch (1) stamped at random locations without any affine transforms.

Se example of images here 1
2
3

I would like to find out what is the stamp in order to eventually conduct de-convolution.

signal/image processing is not my profession, so I don't know if this fast is easy and solved long time a go, or if it is still difficult. But I guess it is easier than what is done in image processing today e.g. here or here since at least I don't need to consider any affine transforms.


texts:
Might be a simple question, but I keep second guessing my algorithm, so a confirmation would be good. I first do bandpass filtering, then PCA, then FFT of the top principle components.

I also do some sine wave fitting on the filtered data, and on the principle components, but it doesnt make sense to do it on the fft data? Right?

Thank you for your patience on this stupid question.


texts: non-consistent normalization problem in fft

texts:
I have a csv file of data sampled with Ts=1ns which looks like this:

This signal is a step response of some system which responds to a step of value 1.
I'm trying to get the impedance profile of the signal by using
$$
Z\left( f \right)=\frac{V\left( f \right)}{I\left( f \right)}=\mathcal{F}\left\{ \frac{1}{A}\frac{\partial }{\partial t}v\left( t \right) \right\}
$$
The signal is noisy and needs filtering, so i'm separating the droop part from the non-droop part and filtering them independently using savgol filter.

Than, i'm applying FFT on each of the signal parts above using the following code:
n = 1024
from scipy.fftpack import diff
frqs = np.fft.fftfreq(n, d=1.0 * 1.e-9)[range(n // 2)]
Y = np.fft.fft(diff(data_y, order=1), n=n, norm='ortho')  # fft computing and normalization
Y = abs(Y[range(len(frqs) // 2)])
return frqs, Y

Next, I take the high frequencies from the FFT applied on the droop part and the rest from the non-droop part and I stitch the results and plot them in log scale.
We get a result with expected shape as below:

Now, the problem arise when we change the droop part to be longer (and thus, the non-droop part to be shorter) as seen here:

The output of this looks like this:

The noisy shape is expected, but the amplitude is changing drastically and I can't figure out why.
can someone point what the problems are are help with a solution?
Thanks!


texts:
Consider the following DFT plot below, showing the output spectrum of an 9-bit ADC.

I want to calculate the THD (Total Harmonic Distortion). I am wondering how to do it correctly. I am aware of the following formula.

$$
THD=10\log\Bigg(\frac{P_{HD}}{P_S}\Bigg)=10\log (P_{HD}) - 10\log (P_S)
$$

So from my understanding $P_{HD}$ is the sum of all harmonics. But the signal power is shown in dBc, so I can only say how far the power of the harmonics are below the carrier, but I do not have the information about the absolute power levels. Is this correct? Then is it actually possible to calculate the THD, if yes, how?


texts: What would be the best solution to synchronize data from two inputs. (Camera and Acceleration sensors)

texts:
For a video stabilization project, a USB webcam (25 fps) and USB acceleration sensors are being used.

the sensors gives a speed of approximately 50 samples/second.

the derivative of the acceleration aka. 'Jerk' is done in order to localize the start and the end of the shake in the video.

the purpose is to compensate the camera movement in real-time, if it were you how would you :

Match the frame (from camera) with the data of the sensor.
Calculate the number of pixels to compensate from acceleration sensors data.

Thank you.


texts: Calculating $E_b/N_0$ from given SNR ratio

texts:
I have the following formula:
$$\frac{S}{N}=\frac{E_b \cdot R_b}{N_0 \cdot B }$$
Where $S/N$ is of course SNR in dB, $E_b$ is energy per information bit, $ R_b$ is information bit rate in [bit/s], $N_0$ is noise spectral density and $B$ is bandwitch in [Hz].
So let's say the $S/N$ ratio = $70\ \rm dB$

$R_b = 250\ \rm kb/s$
$B = 1\ \rm MHz$

So:
$$70\ \mathrm{dB}=\frac{E_b \cdot 250 000\ \mathrm{b/s}}{N_0 \cdot 1000000\ \mathrm{MHz}}$$
$$70\ \mathrm{dB}=\frac{E_b \cdot 1\ \mathrm{b/s}}{N_0 \cdot 4\ \mathrm{MHz}}$$
Now, of course, I can't simply multiply because this 4 is not on a logarithmic scale. But can I just logarithmize 4 like that? Won't there be a problem with units? How to do it to get $E_b/N_0$ ratio in decibels, because I keep doing something wrong.


texts: Is it due to vocal length normalization during MFCC? Please be kind enough to explain.


texts:
I have some Points and Lines from original image and the corresponding points and lines from top view image.

-> I tried to find homography using points

Used OpenCV function findHomography(objectPointsPlanar, imagePoints) :

It takes 4 points from original image and 4 points from top view image(required image) and it returns homography matrix of 3x3.

Used mathematical formula to find homography matrix (Didn't use any library) :

Please look at the attached images for the formula to find homography matrix using 4 points by your own.

I have tested the formula by writing it in c++. It gives the homography matrix same as given by opencv function findHomography().

PS : I used the determined homography matrix to warp the image, but I am not happy with the output. It doesn't come to be accurate.
So I am thinking to use lines to find homography matrix.

-> Homography determination using Lines :

I tried to look for any library which can find homography matrix using lines, but I didn't get.
I also looked at research papers to get the mathematical formula to find the homography matrix, but was not able to understand the formula.
I am sharing the link of a research paper  : A new normalized method on line-based homography estimation

So could someone please help me to find the homography matrix using lines.

Thanks in advance,
Regards.


texts: Why transforming the data to a high-dimensional feature space in which classes are linearly separable leads to overfitting?

texts:
I read in my book (statistical pattern classification by Webb and Wiley) in the section about SVMs and linearly non-separable data:

  In many real-world practical problems there will be no linear boundary separating the
  classes and the problem of searching for an optimal separating hyperplane is meaningless.
  Even if we were to use sophisticated feature vectors, $\Phi(x)$, to transform the data to a
  high-dimensional feature space in which classes are linearly separable, this would lead
  to an over-fitting of the data and hence poor generalization ability.

Why transforming the data to a high-dimensional feature space  in which classes are linearly separable leads to overfitting and poor generalization ability?


texts: EKF smoothing for prediction at t=0 when no there is no measurement

texts:
I have a simple first-order reaction batch system for which I have some discrete measurements ($0<t_{k}\le t_{endbatchsample}$).

I have an initial guess for $x_0$ and $P_0$ and from here I perform the time update step (in a continuous-discrete EKF fashion) until $t_1$, where the measurement update kicks in. This process is repeated until $t_{endbatchsample}$ In this particular application it is useful for me two things:

Prediction from $t^+_k$ (posterior) until the $t_{endbatchsample}$ with covariance (in this case I just have one state, variance) propagation through the continuous Riccati equation. This is implemented.
Correction of  $x_0$ and $P_0$ after the batch. For this, I have thought of the (Extended)RTS smoother, but I am not sure how to perform it since there was no measurement at $t_0$, hence no EKF estimate. Is there a conventional way to predict (extrapolate) the new $x_0$ and $P_0$ (at $t_0$) after a series of measurements are available?

Sorry for any confusion, this isn't really my field.

Best regards


texts: The current state of the art in audio signal interpolation

texts:
Three questions:

What are all the metrics one can use to measure audio interpolation quality, objectively? (but also in terms of psychoacoustics if possible)
By those metrics, what is the current state of the art in audio interpolation?
Suppose I were to render two files from a sequence of notes from virtual instruments in two resolutions and then compare an upsampling of one file with the high-frequency rendered version, what software could one use to compare these objectively? - ideally using the before-mentioned metrics

So far, from I've been able to gather, these resamplers provide some of the best quality

http://www.mega-nerd.com/SRC/
http://sox.sourceforge.net/SoX/
http://www.izotope.com/tech/src/

One of the problems that these resamplers seem to have is pre- and post ringing.

I should note that of key interest is signal reconstruction (insofar as that term makes sense), so upsampling more than downsampling.

Edit: Interpolation time efficiency is irrelevant in this context.

Best regards,
The curious :-)


texts:
I have got a question concerning the definition of the frequency vector for an fft operation.

Generally, I work with a frequency vector, f, with power of 2 elements (2048, 4096, 8192, ...).

Given a certain simulation analysis time, time (e.g. 600s), I should define f as follows:

% Frequency defition
t = 0:dt:(time-dt);
df = 1/(time);
fn = Nfft/time;

$$ f = -f_{n}/2:df:f_{n}/2-1; $$

where $ f_{n} $ represent the Nyquist cut-off frequency.

Actually, for:

computational reasons
symmetry of the power spectra along f axis
not throwing away real or imag part of the fft

I aim to define only half of the frequency range, as for example

$$ f = 0:df:f_{n}/2-1; $$

After calling the fft of my input signal, I would get the desired time series as

ouput = [real(fft) imag(fft)];

But, this way, I count the 0 frequency term twice and the -fn/2 is completely discarded.

How would it be possible to emcompasses the whole standard frequency range starting from only half of it?


texts: Role of Power Metric and LPF in AGC(Automatic Gain Controller)

texts:
I am trying to understand & implement simple toy AGC system . Few questions while reading about it,

1) How POWER metric came into picture, why not other metric like absolute value to control amplitude of input signal?

2) What's the role of LPF(Low Pass Filter) in AGC? LPF is an attenuator or blocker to high frequency components. How such filtering effect can control input amplitudes?

Any reading references will be helpful.


texts: what is the correct answer for shifting and scaling this signal?

texts:
The discrete signal is provided in the following table:

$$\begin{array}{|c|ccccccc|}\hline
 n& -3 & -2 & -1 & 0 & 1 & 2 & 3\\ \hline
x[n] & 3 &  2 &  1 & 0 & 1 & 2 & 3 \\ \hline
\end{array}$$

The book is asking to determine $x[3n-1]$. This involves shifting by one and then scaling by $3$. My answer to this problem is:

$$\begin{array}{|c|cc|}\hline
 3n-1 &    0 & 1\\ \hline
x[3n-1]&  1 & 2 \\ \hline
\end{array}$$

However the answer in the book is:

$$\begin{array}{|c|ccc|}\hline
 3n-1 &  -1 & 0 & 1\\ \hline
x[3n-1]& 0 & 1 & 2 \\ \hline
\end{array}$$

Which answer is correct?


texts: Understanding where the constant $2/N$ comes from in Fourier transformation

texts:
I'm implementing Fourier transformation in my analysis and I wanted dig a bit deeper on the reasons why the absolute value of Fourier transformation is usually multiplied by the constant $2/N$ to get the peak amplitude value of a sinewave with certain frequency.
In the book Understanding Digital Signal Processing by Lyons, the author states the following relationship between the peak amplitude $A$ of a sinewave and the output magnitude $M_r$ of the discrete Fourier transformation (DFT) for that particular sinewave is:
$$M_r=AN/2,\;\;\;\;\;\; \tag{1}$$
where the $r$ stands for real input values to DFT and $N$ is the number of input values to DFT. From this relationship, I trivially get the amplitude I want to as $A=2M_r/N$, which I see is many times done in many fft-examples in Matlab found throughout the web.
Now my big question was, why is the relationship in $(1)$ true? I started to read more from the book Fourier Analysis and Its Applications by Folland and I found the following in his book (in section about DFT):
$$\widehat{f}\left(\frac{2\pi m}{\Omega}\right)\approx \frac{\Omega}{N}\widehat{a}_m,\;\;\;m=0,1,...,N-1\;\;\;\;\;\;\tag{2}$$
where $\widehat{f}$ is the amplitude function, $\widehat{a}_m$ is the $m^{th}$ output of DFT, $N$ is again the number of inputs and $\Omega$ is the length of the time interval $[0, \Omega]:$
$$\widehat{f}\left(\frac{2\pi m}{\Omega}\right)=\int_0^\Omega e^{-2\pi i m t/\Omega}\;f(t)\;dt,$$
where $f$ is the wave function. Now when I look at $(1)$ and $(2)$, there seems to be a connection between them:
$${\color{red}{\widehat{a}_m}} \approx{\color{blue}{\frac{N}{\Omega}}}{\color{green}{\widehat{f}\left(\frac{2\pi m}{\Omega}\right)}},\;\;\;\;\;\;{\color{red}{M_r}}={\color{blue}{\frac{N}{2}}}{\color{green}{A}}.$$
These two results are almost satisfying but I wondered why it seems to be the case that $\Omega=2$?
My questions: Where does this $2$ come from? Why in Lyons's book there is $2$ instead of $\Omega$?
I thought could it be somehow related to the symmetry of the DFT output? One time unit to left and right: $[-1,1]$ so the length of the interval would be $\Omega=2$? A bit vague this last part but could I be onto something here?
UPDATE:
The definition for DFT in book Understanding Digital Signal Processing is given as:
$$X(m)=\sum_{n=0}^{N-1} x(n) e^{-2\pi i nm/N},$$
where $x(n)$ is some continuous time-domain signal. In the book Fourier Analysis and Its Applications the corresponding definition is:
$$\widehat{a}_m = \sum_{n=0}^{N-1}a_n e^{-2\pi i mn/N}\;\;\;(0\leq m<N),$$
where $a_n = f\left(\frac{n\Omega}{N}\right)$.


texts: Does a RHP zero imposes some limitation in the time domain?

texts:
I am studying control systems and I have seen that a non-minimum phase zero imposes a limitation on the bandwidth, which is that the bandwidth has to be lower than the frequency which is half the frequency of the RHP pole.
I have also seen that it imposes in the step response an undershoot, but as far as I know there are no other influences of RHP zero in the time domain.
Does a RHP zero impose any other limitation in the time domain?


texts: Advice on generating stable FMCW waveform

texts:
I am trying to derive parameters for a triangular FMCW waveform such that the phase of the signal has consistency from one period to the next. Perhaps this is arbitrary and feel free to tell me so, but I have a bistatic FMCW system and as such I will be implementing a delay and freq offset calibration procedure on start up to sync the Tx and Rx sides.
So I wish to set my triangle FMCW period, FMCW bandwidth, and sample rate such that the resulting time domain output is clean and steady, however for some reason I'm not doing well with deriving it via some simple trig.
$$
\cos(2\pi f_0 t + \pi Rt^2) \\
\text{where } R=\frac{B}{T_p}
$$
so at $t=0$ the argument is also $0$, so I would think I need to solve for a value of $T_p$ such that when $t=T_p$, the argument is equal to some large multiple of $2\pi$ while keeping in mind the sample rate should yield a whole number of samples per period.
This approach has not been working however, and I'm at a loss of what the issue is.
For additional info
$T_P$ is half of the full period of the triangular modulation waveform
for current parameters I have set
$$
sample rate = 8M\\
B = 1.6M \\
Tp = 1.25ms \\
$$
EDIT:
Images to illustrate point

Note


texts: Specifying the frequencies when designing a filter

texts:
I understand the general method of designing high pass and low pass filters with a specific number of poles and zeros. However, if given a specific sampling frequency, let's say 2.4 kHz, and asked to make

the peak of the high pass filter to be 1.2 kHz
and
the peak in the passband region of the bandpass filter 1.2 kHz.

How would I specify these frequencies for the peak and band pass region of the filters?


texts: Difference between frame rate and sampling rate?

texts:
Below are the two different methods of generating time to create an audio sine wave. While framerate is being used as steps in one method, the inverse of sampling rate is being used as steps in the other method. Where is the difference?

Method 1

# Create a sound
framerate = 44100
time = np.linspace(0,5,framerate*5)         #Multiplying framerate

signal = np.sin(2*np.pi*220*time) + np.sin(2*np.pi*224*time)

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

Method 2

#Simulation parameters
srate = 500 #Sampling Rate in Hz
time = np.arange(0.0, 2.0, 1.0/srate)       #Dividing Sampling rate

#Sine Wave Parameters
freq = 3         #Frequency in Hz
ampl = 2         #Amplitude in a.u.
phas = np.pi/3   #Phase in Radians

#Generate the Sine Wave
sine_wave = ampl*np.sin(2*np.pi*freq*time + phas )


texts: Confusion related to PN sequence : Terminology and application in statistical signal processing

texts:
A Pseudo-random Noise (PN) sequence is a sequence of binary numbers, e.g. ±1, which appears to be random; but is in fact perfectly deterministic. The sequence appears to be random in the sense that the binary values and groups or runs of the same binary value occur in the sequence in the same proportion they would if the sequence were being generated based on a fair "coin tossing" experiment. In the experiment, each head could result in one binary value and a tail the other value.

Consider a model, $$ x(t) = h^T u(t) +w(t)  \tag{1}$$
where $u$ is the input and $w$ is Additive White Gaussian noise.

In Book : Fundamentals of Statistical Signal Processing by Steven Kay, Chapter 4 : Eq(4.21) finds the variance of the estimated coefficients when the input is PRN as

$$var(\hat{h_i} )= \sigma_w^2/Nr_{uu}[0] \tag{2}$$ under the assumption of the input to be a PN .
where $\sigma_w^2$ is the variance of the measurement noise; $\sigma_u^2$ is the variance of the input
The MVU estimator using a PRN sequence is mentioned as $\hat{\theta }= (H^TH)^{-1} H^T x$

It is stated that for the deterministic sequence $u$ the $H^TH$ becomes a symmetric Toeplitz autocorrelation matrix. For $H^TH$ to be diagonal we should use a PRN sequence as input which yileds $H^TH = Nr_{uu}[o]I$.
These expressions can be found in document Pg 7 Section 4.1.3

My problem is that I am confused if PRN means binary or $\pm1$ sequence or only $iid$ real valued White Gaussian noise. If PRN is different from PN or binary, then what will be the estimates and the CRLB?


texts: multiplications per second for FIR filter

texts:
In my filter the fs=300MHz. and no. of coefficients is 31. I have following queries.
(1) what will be the number of multiplications per second?
(2) can multiplication be done in one clock cycle?
(3) why normally multiplication per second is calculated and not addition/substraction per second in the filter implementation?


texts: Understanding the Eye Patterns: analysis

texts:
I got the following  Eye Patterns:

how to understand them? what kind of information can I get?


texts: Reconstruction of Audio Signal from its Absolute Spectrogram

texts:
I have the absolute Spectrogram of an audio signals.
I lost the phase data of the Spectogram because of various processing applied on the original spectrogram of the signal.

I'm trying to reconstruct the audio signal in a meaningful (Audibly) manner from teh absolute value only of the Spectrogram.
The obvious inverse won't work (The DFT inverse of the absolute, Since the Phase is significant).

The Spectrogram is a result of fusion of few audio signals as I'm trying to create a smooth transition between audio signals.

Anyone has experience with the problem?  Anyone has experience with this procedure?
Could anyone refer me to a code, article, etc...

Thanks.


texts:
Based on another question in  this subforum:
SIFT - why s+3 scales per octave?
"For s=3 this means you will have s+3=6 blurred images (the Gaussian images shown in the paper in Figure 1 on the left). Having 6 Gaussian images will result in 5 DoG images (shown in Figure 1 on the right). This will allow you to do the extrema detection on s=3 scales (using the method shown in Figure 2)."
The above mentioned is pretty clear, now I was left wondering what the used sigma for the 6 images is? Let's say we go for 3 scales per ocatave like suggested, meaning we have to detect extrema points in 5 DoGs, making it necessary to smooth 6 images in the octave. According to the given formula $k = 2^{1/s}$, with $s=3$ we would need to more than double our initial sigma for 6 images due to $\sigma_n = k^n \sigma_0$, where $n$ is the number of images in one octave and $\sigma_n$ the used sigma for image $n$ of the octave. Since actually the first formula for $k$ is used to compute the number of intervals until $\sigma_0$ doubles, it feels odd. Is this intended or should the formula actually be $k = 2^{1/(s+3)}$.

[1] Distinctive Image Features from Scale-Invariant Keypoints
D. G. Lowe
Int. Journal of Computer Vision 60(2) (2004), pp. 91--110


texts: Comparing the spatial frequency between groups of images

texts:
I have 2 groups of images, and I want to compare these groups based on the spatial frequencies of the images.
In the end, I'd like to be able to say "group A is/is not significantly different from group B based on spatial frequency" or something to that effect.

Can anyone think of a way to do this? This question seems to touch on the topic, but as far as I can see it doesn't provide the solution I'm looking for (particularly not for groups of images).
MATLAB solutions would be the most valuable to me.

Thanks for any help!


texts: phase offset extraction using matched filter

texts: I'm looking for an application of detecting the drift of a received signal carrier frequency from ADC clock of the receiver (in case of no drift, ADC clock should be optimally an integer multiple of signal's carrier frequency), which can be observed as a varying offset component of the instantaneous phase. one proposed solution that I found was to use a matched filter to extract the variable offset. Yet I didn't find any supporting description to understand the mechanism. did anyone come across such an implementation?


texts:
What algorithms are available for optimizing estimated bit clock recovery of an audio FSK modulated signal over an extremely low S/N ratio channel? (e.g. barely audible due to high ambient noise and low max available signal power)

Assume the audio signal is very long in duration and the payload is suitably pseudo-randomly pre-scrambled.  Even if the data cannot be recovered (or perhaps without massive super-redundant FEC), can the clock be recovered accurately?


texts: Relative beat detection for two audio recordings

texts:
I'm trying to synchronize (overlap) two audio recordings of the same song: a complete HQ version, and an incomplete, noisy version (phone recording for example). The noisy recording may have had its tempo increased/decreased.

Currently I have something working based on this 2001 paper: "The Beat Spectrum: a New Approach To Rhythm Analysis"

In short:

Get a self-similarity vector for each song, which will contain peaks at every multiple of the song bar duration (every 4 or 8 beats for latin music for example)
For each r (relative tempo ratio) from 0.7 to 1.3, scale the second vector and then do a dot product on the two vectors, using v1[i]·v2[i] = d1[i]*d2[i] (where d1[i] is the sign of v1's derivative; same for d2)

When the beats are clear in the recordings, I get a clear maximum at the correct ratio (usually 1.0). The problem is that when the beats are lost (due to EQ, bad mics, noise) the self-similarity chart stops showing the periodic bars and I get maximums in the wrong places in the dot product.

My question is, what else could I use instead of the dot-product-of-the-derivative-sign, and what other robust ways are there to determine the tempo change ratio in the presence of noise, besides doing BPM(song1) / BPM(song2).

In case someone can take this further: computing self-similarity by overlapping the sample twice (at [0, ofs, ofs+ofs] instead of just [0, ofs]) gives much more accurate maximums for each beat duration. But I haven't tested it yet on all the samples I have.