texts: CIC Filter implementation for Normalized Correlation

texts:
This is a follow-on to the original post ...

I am trying to implement real-time correlation of two signals. I noticed the CIC Filter discussion here:

How to find correlation/cross-correlation of two signals in real time?

... however I need to use a normalized correlation for detection.  What would the complete block diagram or equation be for a normalized CIC Filter?

Additional:
I don't have to use the CIC - it just appeared to be computationally cheap and simple to implement.  Would a moving window sum of the signal product, divided by the product of stdevs (to normalize) also work?

I.e., if I took the most recent N samples of each signal, and computed the normalized cross correlation, would that be the typical approach?

Follow-on:  I have experimented with various ways of measuring signal similarity, and am getting some odd results.  My scenario is two signals, X and Y, aligned in time.  The X series counts from 0 to 20 over 21 samples (N).  The Y series counts from 20 to 0.

I calculate the Pearson product-moment correlation coefficient over the entire series and get an R value of '-1', as expected.  But since I'm interested in getting (near) real-time indication of similarity, I implemented the CIC filter approach, as discussed below, and normalized the value by dividing by the square root of the product of two other CIC filters (that each square the X and Y inputs, respectively).

So I am computing the normalized cross correlation of X and Y, using a 'window' of 3, by:

R(normal) = Rxy / (Rxx * Ryy)^0.5, using three CIC structures suggested by Boschen in:

How to find correlation/cross-correlation of two signals in real time?

Note I am not scaling this by N, as suggested by the answer below (I'm still not clear on this).

Now the problem ... the normalized value I get is '+1', not '-1'.  Also, the result using the three CIC filter structures seems to be sensitive to when the data series is near, or crosses, zero.  My guess is I'm still not implementing this correctly ... or is the CIC approach just sensitive this way?

As a further experiment, I implemented a 'windowed' Pearson correlation calculation, computing normalized value 'R' for the most recent 3 samples, and obtained correlation values near '-1', as expected, and it does this consistently regardless of whether the data series is near or crosses zero.

The CIC approach would seem to be preferred as a faster computation, but is giving inconsistent results.


texts: Where to get a proper dataset collection for bikes and cars

texts:
I am doing a project on image classification and I don't know where to get the data-set collection for cars and bikes.Any suggestions to sites for this purpose is welcome.Also please suggest the general size of images to act as a database.

Thanks and regards in advance.


texts: difficulties implementing windowed sinc interpolation (C++)

texts:
I'm currently trying to implement a windowed sinc interpolation. I've already written some code for that, but it only seems to work for cases where phaseInc <= 1.0 (phaseInc = outSampleRate / inSampleRate). By "not working" I mean that something is causing a lot of distortion or aliasing. I hope that the comments in the code do explain what I'm trying to do:

// SINC_FILT_THRESH = 0.8, SINC_WINDOW_SIZE = 16

float sincStep = phaseInc > SINC_FILT_THRESH ? SINC_FILT_THRESH / phaseInc : 1.00f;
/* ^^ this should in theory detect the transition between phaseInc <= 1.0 and phaseInc > 1.0
 * but with a little headroom for the transition band. "sincStep" then shall be used later
 * to make the sinc function wider to narrow down the rectangle in frequency domain
 * The code below seems to work as long as sincStep = 1.0 */

int i = 0;
do {
    float sampleSum = 0.0f;
    float kernelSum = 0.0f;
    for (int wi = -SINC_WINDOW_SIZE + 1; wi <= SINC_WINDOW_SIZE; wi++) {
        float sincIndex = float(wi) * sincStep - phase;
        float windowIndex = (float(wi) - phase) / SINC_WINDOW_SIZE;

        float s = sinc(sincIndex);
        /* ^^ this is a normalized sinc function with zeroes at whole numbers */

        float w = windowFunc(windowIndex);
        /* ^^ windowFunc = 0.5 + (0.5 * cos(pi * t))
         * I know there is stuff like blackman window, but I'd rather like to start with
         * something simple */

        float kernel = s * w;
        sampleSum += kernel * fetchBuffer[i + wi + SINC_WINDOW_SIZE - 1];
        kernelSum += kernel;
    }
    phase += phaseInc;
    int istep = static_cast<int>(phase);
    phase -= static_cast<float>(istep);
    i += istep;
    /* ^^ This simply advanced phase and i by the appropriate amount.
     * I used it with a simple linear interpolator and didn't have problems with that */

    *outData++ = sampleSum / kernelSum;
    /* ^^ Store the result in the output buffer with unity gain */
} while (--numSamples > 0);

So, as I've mentioned it only seems to work right in cases where sincStep = 1.0 and I don't know why. My intention is to widen the sinc function to lower the cutoff. Doesn't seem to work quite as I want.

What I do know is that there are no memory access violation problems (according to fsanitize) and my fetchBuffer always has enough data (so no out of bounds thing).

I'm sorry if I'm missing something quote obvious but I really don't want to just copy someone elses code. I'd like to actually understand it. Most examples however just seem to focus on the equations and not on the actual thing with all the possibly faulty boiler plate code.

Thanks in advance!

EDIT:
I want to use that interpolation for resampling different waves at different pitch for music. I thought I could add this for clarification for my problem:


texts: What transform can I use to "zoom" in a specific portion of a signal's spectrum?

texts: If I recall correctly, there is a variation of the DFT that can be used to analyze a specific band of the spectrum of a signal. How is it called?


texts: Creating a prediction from data in MATLAB

texts: How can I create a closed loop system that can guess the input value that is going to return 0 from unknown transfer function. The output of this transfer function is  either a positive constant, a negative constant, or 0. I am using MATLAB.


texts:
I know that FSK occupies a larger bandwidth compared other keying schemes,

$$
BW= f_1 -f_2 +2R_b
$$

where $f_1$ and $f_2$ are carrier frequencies and $R_b$ is the bit rate.

For achieving a particular bit rate say 8 Kbps, how can I find the minimum bandwidth if I am free to choose the carrier frequencies?


texts:
I'm working on an echo implementation on the board adsp 21262. I've added a buffer to have more samples to work with. I end up then with 4 buffers of size 2048 each.

The problem is that I need to access at least 7000 samples back. (An echo can only be heard with a delay of 1/15th of a second, we have 2*48K~100k samples per second. So we need 100k/15~7k samples).

The problem is that 2048 samples per buffer is the max we can have because afterwards we're out of memory.
I'm trying to access the ROM but I have no idea as to how to do that. I've spent a lot of time looking for this information finds nothing.

Thank you !


texts: Steps for recognizing hand written signs

texts:
What are the steps to built a system which recognizes hand written signs in an image and returns the most similar sign in Python? I know that depends on the image and many variables, but isn't there a sequence of steps, a main path to follow for recognizing images with at most 2 or 3 signs per image?

I thought the steps should be:

Take an image and convert to binary
Apply find_contour to binary image
Extract scale- and rotation-invariant properties
Create a feature vector with scale/rotation-invariant properties (if there is much more than 1 sign, like for example a star and a square not overlapped, is that a problem?)
Repeat steps 1 to 4 for every image
Put them into a k-NN or SVM

These next steps, I do not know how to implement:

I choose an image
System spits out the result (most similar image)


texts: Why are analytic signals so important in Time-frequency analysis?

texts: I am little confused about why we need analytic signals so bad in time-frequency analysis. What might happen if I use non-analytic signals to do time-frequency analysis?


texts: Resolving two nanosecond pulses with filter

texts:
Problem:

I have to design a filter that should be able to resolve two ~20ns pulses as well as get rid of some ringing. I was able to find a filter kernel that if convolved with the noise signal is able to do the required filtering. Below are the graphs of the pulses before and after applying the convolution.

The graph below is the filter kernel that I used to do the filter in both the time and frequency domain.

Attempted solution:

According to my understanding -- which is very superficial -- the discrete Fourier transform of the filter kernel is the frequency response of the filter. Based on that I went on and tried to design a filter with a similar frequency response which according to the graph is a low-pass filter. However, when I design and apply a low-pass iir filter to the signal it is not able to resolve the pulses. What's wrong with my understanding? How could I go about designing a filter that has a similar filter response as the one I require? Is it possible to do it with an analog rather than a digital filter?

PS: I was not able to upload more pictures but if possible I will upload them to the comments as well as the code and data use to produce the graphs.


texts: MATLAB's $\tt cpsd$ and $\tt pwelch$ - different results for cross spectral power density

texts:
I am wondering why the following code does not yield the same results for Sxy1 and Sxy2, where Sxy1 is calculated using the MATLAB's pwelch function, and Sxy2 is calculated using cpsd function.

fs = 1000;
t = 0:1/fs:5-1/fs;
win=hann(500);
overlap=300;
nfft=500;
fs=1000;

x = cos(2*pi*100*t)+randn(size(t));
y=sin(2*pi*160*t)+cos(2*pi*150*t)+randn(size(t));

[Sxx,f] = pwelch(x,win,overlap,nfft,fs);
[Syy,f] = pwelch(y,win,overlap,nfft,fs);

Sxy1=sqrt(Sxx.*Syy); %******

[Sxy2,f] = cpsd(x,y,win,overlap,nfft,fs);

Sxy2=abs(Sxy2); %******


texts: Root raised cosine pulse shaping peak value

texts:
I am using a FIR filter to generate Root Raised Cosine pulse for an IQ modulator.

How can I predict the peak value at the output, depending on the roll-off factor of the pulse ?

Knowing that :

FIR taps are root raised cosine pulse, with a peak value of 1.0
Input of the filter is I or Q channel, valued +1 or -1, upsampled by adding 0s between each symbol

By experimenting, I can see that depending on the input sequence and the roll-off factor, the peak output varies. But I can't predict it theoretically.

This plot illustrates my question (while this particular sequence may not generates the maximum peak that I want to predict) :


texts: Digital Equations problem - The fundamental ones that everyone should know

texts:
I think I am looking at conflicting equations from a few sources, or maybe I just dont understand it. But please can someone help me with these fundamental equations that everyone should know

$  \Omega $ is Discrete frequency
$  N $ is Discrete Period in number of samples
$  \omega $ is Continuous frequency
$  T $ is Continuous time period
$  F_{s} $ is Sampling frequency
$  T_{s} = n $ is Sampling Time Instant
$  k $ is Discrete harmonic index which goes from 0 to N/2 for a real signal

First equation forms I have seen
$$
\Omega = 2\pi f/F_{s} =\omega / F_{s} = \omega T_{s} = \omega n = 2\pi n / T
$$
Second equation I have seen
$$
\Omega N = 2\pi k
$$
$$
\Omega  = 2\pi k / N
$$
Third equation I have seen
$$
\Omega_{k} = k F_{s} / N
$$

First and second look the same with if $k=T_{s} $ and if $N=T$ and first and third seem contrary.

Can someone explain please?
Also if you can point in in a direction where all this is neatly derived and defined?

Thanks


texts: Why Normalized Cross Correlation Score of different templates is high?

texts:
I am using OpenCv's built in template matching function to search for an object in  image. I am using Normalized Cross Correlation Method.The function is returning a value which I think indication of similarity so the larger value the more similar template.The problem is NCC value when object is matched is 0.93 whereas NCC value when different regions are found is 0.91 .

Is the value returned by NCC percentage of similarity ? If so how can different templates produce high score ?


texts: How do I implement an Infinite Impulse Response (IIR) filter using Finite State Machine(FSM)? The problem is FSM has finite memory, but an IIR filter has infinite outputs. The closest match is this paper  . Unfortunately I do not have access to the paper.


texts: Correct method for processing FFT in realtime

texts:
I am attempting to build a 'live' FFT for continuous analysis of an audio signal.  I am relatively new to this, but have attempted the full thing, and then updating the graph at 20Hz refresh rate.  My FFT is working with regards to identifying the correct frequencies, however is very 'jumpy' with regards to the amplitude. Before I analyse my code to see where I have gone wrong, I would like to just confirm the process I have used is correct.

Read audio samples into circular buffer.
Read half of the window size from last half of previously processed array.  (50% because I am trying to achieve 50% overlap).
Read half of the window size from the buffer to give the complete window size.
Perform Window Function on the window to process.
Zero pad the end of this to achieve a n^2 bit FFT Size.
Perform FFT on the full array
Save Array to List of arrays to average (add to end, delete first one to always retain x no. of averages).
Average the arrays in List of arrays (treat real and imaginary values as two independent values when averaging).
Obtain the magnitude for each bin in the averaged array by obtaining the square root of the imaginary^2 + real^2. This is then divided by the number of bins, and the entire value multiplied by 2 as we have halved the FFT spectrum to only look at the positive values. This is then multiplied by a window correction factor.
This thread should then sleep for however long it takes to capture half the window length needed (the next set of audio samples) subtract the computational time to perform this FFT step.

That is currently the method I am working to.  I know my FFT algorithm is good as I have tested it on known input and output values.  I am sure it is a coding mistake, but as this is my first project I wanted to check that the overall method I am following is sound.  If anyone can see any potential issues, please let me know.  Thanks.


texts: Calculating the power spectral density on a different bandwidth

texts:
Let us say I have the power spectral density for a given bandwidth F1 at a transmit power T. I would like to calculate the power spectral density for bandwidth F2 for 1) the same transmit power T
2) for different transmit power T2

How can I calculate this

For ex: Let's say I am looking at the spectrum of an OFDM signal in the baseband with a BW of 20Mhz transmitting at power of 23dbm and i know the power spectrum of this signal. Now I would like to derive the power spectrum of the OFDM baseband signal with a BW of 10MHz transmitted at the same power


texts: Why removing noise increases my audio file size?

texts:
I took some mp3 files of people talk with noise in the background and then removed the background noise using online voice seperator tools (such as https://vocalremover.org/). I was surprised to find out that the size of my files increased (about 5 times bigger).
Why did it happen? Is there another audio format or background noise removing method that would reduce my file size when removing the background noise?


texts:
I know non-sinusoidal waveforms can be decomposed into a sum of higher integer harmonics, but are there any mechanisms that cause interharmonics?
I collected some noisy data from a source that was forced at a known frequency. I applied a Blackman-windowed FFT and found not only many integer harmonics, but also 1/5 integer harmonics (ie, 0.2, 0.4, 0.6, 0.8, 1., 1.2, 1.4, etc). Are there any signal processing causes of fractional harmonics, or should I be looking for a physical mechanism?


texts: Developing DSP software, am I 10 years too late?

texts:
I'm interested in developing software related to digital signal processing. The most obvious ideas are audio and video editing. There's an overwhelming amount of audio and video editing software out there. I feel like developing such software would be in vain.

What are the most lucrative frontiers in DSP software today?
Do you think audio/video software still have room for improvement?


texts: Apply timeshift to undersampled signal after bandpass after sampling?

texts:
I have sampled my continuous signal $x(t)$ at $x[n]=x(n/f_{\rm samp})$.
The signal is undersampled, and I know I will not be able to recover my original signal anymore. After sampling a bandpass is applied to $x[n]$. The signal is now bandlimited, but is affected by aliasing. I now have

$$
x[n]=\textrm{Bandpass}\left[x\left(\frac{n}{f_{\rm samp}}\right)\right]
$$

However, can I somehow reconstruct

$$
z[n]=\textrm{Bandpass}\left[x\left(\frac{n}{f_{\rm samp}}+\Delta\tau\right)\right]\quad\text{from my}\quad x[n]\quad?
$$

$z[n]$ would also be an undersampled version of $x(t)$ but sampled with a timeshift.
I have no intentions of recovering $x(t)$, but I need my signal shifted to synchronize it. Even though it is undersampled.

Any ideas or is that a lost case?


texts: What is the Assumptions behind obtaining Better Frequency Resolution based on Transformation of the Random Variables?

texts:
If the sequence is short (N is small), we cannot resolves frequencies being closed together (even by some tricks like zero padding).

One way to solve it is to create a pseudo-sequence from the data in hand by linear transformation of the RVs. Next, overlapping (by about 25%) these two series, we can process the new series with longer (pseudo) data. This method is based on the linear transformation of RVs.

This theorem is proposed in 8.2 by 1. I can realize the simulation using the code shown in below. However, I don't really understand the assumptions behind this operation.

tes.m

n = 1:64;
x =  sin(0.3*pi*n) + sin(0.32*pi*n) + 0.5*randn(1,64);
[ax,ay,w,y] = linear_modified_periodogram(x);

subplot(2,1,1);
plot(w,ax,'k');
xlabel('w');
ylabel('PSD');
subplot(2,1,2);
plot(w,ay,'k');
xlabel('w');
ylabel('Modified PSD');

linear_modified_periodogram.m

function [ax, ay, w, y] = linear_modified_periodogram(x)
    [sx, ax, px] = general_win_periodogram(x,2,512);
    y1 = [x zeros(1,48)] + 0.2*rand(1, 112);
    y2 = [zeros(1,48) x] + 0.2*rand(1,112);
    y = (0.2*y1+(y2)*0.2);
    [sy, ay, py] = general_win_periodogram(y,2,512);
    w = 0:2*pi/512:2*pi-(2*pi/512);
end

general_win_periodogram.m

function [s, as, ps] = general_win_periodogram(x, win,L)
    if(win == 2) w = rectwin(length(x));
    elseif (win ==3) w = hamming(length(x));
    elseif (win ==4) w = bartlett(length(x));
    elseif (win ==5) w = tukeywin(length(x));
    elseif (win ==6) w = blackman(length(x));
    elseif (win ==7) w = triang(length(x));
    elseif (win ==8) w = blackmanharris(length(x));
    end
    xw = x.*w';
    s = zeros(L,length(x));
    for m = 1:L
        n = 1:length(x);
        s(m) = sum(xw.*exp(-j*(m-1)*(2*pi/L)*n));
    end
    as = ((abs(s)).^2/length(x))/norm(w);
    ps = atan(imag(s)./real(s));
end

Understanding Digital Signal Processing with MATLAB and Solutions


texts:
After the estimation of channel impulse response (using the training sequence in burst), to cancel the effect of channel on the received signal, I think, we should convolve the received signal with the inverse of channel impulse response.

But,in GSM, it convolves the received signal with the conjugate of channel impulse response and calls this process as "matched filtering". Could anyone help me with this?


texts: Would someone be able to give me some info or a link etc ... regarding haarlets aka Haar wavelet-like features. I'm reading several papers for my master dissertation and several of these papers briefly mention haarlets. I can't seem to really get what they are and how they are used. I apologize if this is a basic question. But i would really appreciate some clear, structured info on the subject.


texts:
I am fairly new to the signal processing world and that being said I have little to no experience. The problem that I am having is that I am not quite sure how to use dsp.RLSFilter. So far I have only used highpass filter and it was straight forward - I had to just decide on the cutting frequency, type of highpass and sampling frequency, whereas for the RLS filter I have a ton of parameters to choose from. For example, how to decide on the method to calculate the filter coefficients? Furthermore, in the documentation it is stated:

  Call step to filter each channel of the input according to the properties of dsp.RLSFilter. The behavior of step is specific to each object in the toolbox.

and in the example step is not used:

rls1 = dsp.RLSFilter(11, 'ForgettingFactor', 0.98);
filt = dsp.FIRFilter('Numerator',fir1(10, .25)); % Unknown System
x = randn(1000,1);                       % input signal
d = filt(x) + 0.01*randn(1000,1); % desired signal
[y,e] = rls1(x, d);
w = rls1.Coefficients;
subplot(2,1,1), plot(1:1000, [d,y,e]);
title('System Identification of an FIR filter');
legend('Desired', 'Output', 'Error');
xlabel('time index'); ylabel('signal value');
subplot(2,1,2); stem([filt.Numerator; w].');
legend('Actual','Estimated');
xlabel('coefficient #'); ylabel('coefficient value');

So what is the difference between using step and using the method in the example? When I try using the same the same way as in the example I get the following error message: Array formation and parentheses-style indexing with objects of class 'dsp.RLSFilter' is not allowed.  Use objects of class 'dsp.RLSFilter' only as scalars or use a cell array. I tried using num2cell on x and d, however, I had 0 success.


texts: How can the I and Q components of a quatrature signal be different from each other (apart from the 90º phase shift)

texts:
I've spent many hours reading and doing research about this. I don't attend any university classes or any type of classes that may approach such a topic. I am trying to learn this by my own without any help what so ever. So please understand if I say something inaccurate or if I'm asking an obvious question. I have been studying quadrature because I need it to develop in software defined radio platforms. The problem is that I have found a lot of contradictory information.

Let's start from the beginning. Given a real signal to a quadrature circuit, two components are generated, I, the real component, and Q, the imaginary component, which is the same as I, only it has a 90º phase shift. So, we have I and Q, which are the same signal but 90º shifted from each other. I understand this, it's simple.

Now let's see an example. Have a look at this paper about QAM on ni. I read the paper and I understood it. The problem that is if you take a look at the baseband graph, the I and Q components are completely different. They shouldn't, right? They should be equal to each other apart from the 90º phase shift. At least that is what I understood after reading many papers and articles on the subject.

Now let's take a look to this article about IQ on tek. Having a look at the last two images, I can visualize what's happening. The I and Q are being transmitted and then received. I can see that the I and Q components transmitted are the same being received. What I don't understand is that Q isn't by any means equal to I with a 90º phase shift. That is what's supposed to happen right? Maybe not? I can't seem to understand it. I know that the I is obtained by multiplying the real signal by 2cos(w(c)t) and that Q is obtained by multiplying the real signal by -2sin(w(c)t). So, what I saw before clearly makes no sense.

Again, I'm sorry for my ignorance in the subject but I am doing my best trying to understand this. I'm just confused by this contradictions. Thanks in advance for anyone who tries to explain this to me, really.


texts: How to normalize audio data for feature extraction

texts:
I recorded the same piece of music from different instruments and I want to compare their timbre by means of feature extraction. As you can see from the image, the recording is composed of different sections with very different amplitude values. Furthermore, the recordings present a slightly different RMS energy.

I am going to divide the recording into short frames and extract the features from each frame. I want the instruments to be comparable. For this reason I performed a RMS normalization on the whole recording before windowing and feature extraction, in such a way that each recording has a RMS value equal to 1. Am I right, or is there a better approach?


texts: How can we prove the correctness of the integration property of the Laplace transform?

texts:
I was going through an Electrical Engineering textbook for understanding the Laplace transform and came across the following proof for one of the properties of the Unilateral Laplace transform.
Integration property of the unilateral Laplace transform:

In the proof, it is stated that:
$$ e^{-st} \rightarrow 0 \mbox{ as } t \rightarrow \infty \ \ \ (i)$$
and therefore the term:
$$ -\frac{e^{-st}}{s}\int_{0-}^{t}{f(\tau) \ d\tau} = 0 \mbox{, as } t \rightarrow \infty \ \ \ (ii)$$
But my doubt is: Isn't there a case where the real part of $s$, gets cancelled by the function obtained after the integration of $f(\tau)$?
If that happens then we can't guarantee that (ii) would hold true right - i.e. the LHS in (ii) would not be zero right?
So, isn't the proof that has been provided wrong?


texts: Why do we care about the multiplicity of poles and zeros in rational Z-transforms?

texts:
In the DSP class I'm taking, a lot of the questions ask me to list the multiplicity of a pole or zero in a rational Z-transform. For example, the multiplicity of the zero at $z=1$ in :$\frac{(z-1)^2}{(z-5)}$ is $2$.

What knowledge do we gain by knowing the multiplicity of some pole or zero?

Of course, if we know the multiplicity of all the poles and zeros, we can describe the z-transform of a rational function, but this is not a very satisfying answer.


texts:
I'm using a spectral subtractor to eliminate some continuous tones equi spaced in frequency through overlap-add method, some example code in Matlab:
%Windowing parameters
NFFT   =N;
WINLEN =N;
WINSTEP=N/2;

overlap=max(WINLEN-WINSTEP,1);
win=@hamming;
ACw=(WINLEN/sum(win(WINLEN))); %Amplitude Correction of window effect

frames = frame_sig(sig,WINLEN, WINSTEP, win); %frame the signal
cspec = fft(frames,NFFT,2); % complex spectrum
pspec=abs(cspec*ACw/N).^2;
phase = angle(cspec);%save the phase

%Using precalculated subtractor nFilt
load('nFilt.mat')
gSub=1.5;subtraction gain, some audio has more level than others

%Subtraction
noise_est =  gSub*nFilt;% noise estimator
clean_spec = pspec-noise_est;%repmat(noise_est,size(pspec,1),1); % subtract noise_est from pspec
clean_spec(clean_spec < 0) = 0; % negative power spectrum is not allowed

%reconstruction
reconstructed_frames = ifft(  N/ACw*sqrt(clean_spec).*exp(1i*phase) ,NFFT,2);
reconstructed_frames = real(reconstructed_frames(:,1:WINLEN));
enhanced_signal = deframe_sig(reconstructed_frames,length(sig),WINLEN,WINSTEP,win);

The code does the job correctly except when a loud click sound appears, when this happen that frame doesn't work properly as you can see in this image (top original, bottom subtracted):

I guess this has to be related with the sharp phase change that this "click" sounds imply but, any clue on how to solve it?
All the best


texts: Is applying a gaussian kernel a fancy way of saying convolving with a gaussian mask?

texts: When someone says apply gaussian kernel, is that the same thing as saying, convolve it with a gaussian mask?


texts: How to correct phase of Diagonal Volterra kernels obtained with exponential sweep sine?

texts:
I'm trying to model a non-linear system using non-linear convolution with Novak's (2010) synchronized exponential sine weep (SESS) that models them with a Generalized Hammerstein (Volterra diagonal).

I get pretty clear kernels in terms of SNR but clearly phase-distorted. It's a common problem in real estimation.

I can theoretically correct them with a time compacting inverse Kirkeby filter (Farina, 2007).

[EDITED AFTER FIRST ANSWER]
I use Novak's himself code from http://ant-novak.com/Swept_Sine_Simulation.php

Changes: Doubled both number of kernels and sampling frequency and raised a lot SESS duration for better kernels separation.

fs = 192000;         % sampling frequency
f1 = 10;            % start frequency
f2 = 12000;         % final frequency
T_ = 50;            % approximative time duration
N = 8;              % order of nonlinearities

L = round(f1/log(f2/f1)*T_)/f1;
T = L*log(f2/f1); % this definitive sweep duration in theory makes Novak's SESS more robust to phase issues (synchronized)
t = (0:ceil(fs*T)-1)./fs;
x = sin(2*pi*f1*L*(exp(t/L)-1)); % sweep

Novak uses an elegant frequency domain definition alternative of the regular inverse sweep with amplitude correction. I've tested both with similar results, although Novak's frequency method gives better SNR.

The other change is that I model a "real" system, well, only a guitar amplifier VST plugin. So I use the plugin's output to the SESS for the deconvolution stage.

I'm not 100% sure, but the kernels look extremely "right-sided" like in Farina's 2007 paper.
The frequency magnitude of the kernels is coherent with a high gain guitar amplifier VST (even harmonics stronger). The kernels phase is a mess, triangle-saw like.
When I test the model with a dry guitar sound as input the output barely improves a regular (first kernel) equalization matching.

Both in the SESS (plugin input/output) and simulation (dry guitar signal) stages I apply a steep bandpass in the range 100hz-10Khz. The plugin is configured to do x2 oversampling (384Khz) to avoid aliasing artifacts.


texts:
Operate on an image by performing Median Filtering in a 3x3 window.
Operate on the resulting image by performing, again, Median Filtering in a 3x3 window.
Can the resulting image be obtained from a single Median filtering?

my initial thought is that it can be done with the right mask. maybe a median next to a median. but i'm not sure.


texts: Time response of complex frequency response function

texts:
In the frequency domain I have the following function:

$$
H(f)=A(f)\cdot e^{i\beta(f)}
$$

Where both $A(f)$ and $\beta(f)$ are known functions. The thing I would like to do is to model the response of this filter in the time domain to an incoming discrete signal ($x[n]$) according to:

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

Where:

$$
h[n]=\frac{1}{f_{s}}\int_{-\frac{f_{s}}{2}}^{\frac{f_{s}}{2}}H(f)\cdot e^{2\pi i f n/f_{s}}df
$$

When I use a real valued filter, e.g. $h[n]=0.5$ for all values of n, the resulting signal is also real valued. However, when I add an imaginary part to the filter, e.g. $h[n]=0.5+0.5i$ the resulting real part of the signal is the same as for $h[n]=0.5$ and the output signal contains an imaginary part, even though this should only be a phase shift.

How do I convert this to a real signal with a phase shift?

Thanks in advance


texts: How to use Expectation maximization to estimate poisson noise matlab

texts:
I delved into the state of the art of algorithms for Poisson noise estimation in order to estimate the variance, I found that the Expectation Maximization algorithm is very used and it is very effective and easy to use,

But I have not found the matlab code, I found only one that is was used for classification,

can you give me the code or the algorithm for estimating and not for classification please.

I would be grateful if you answer me


texts: Pre-emphasis: What is connection between parameter a and 6db/oct?

texts:
I know that the standard pre-emphasis function is $y[n] = x[n]-\alpha\cdot x[n-1]$ and $\alpha$ is often set as $0.97$.

However, when I was learning speech processing using Praat, the pre-emphasis factor in it is $6\ \rm dB/oct$ which means the signal emphasise $6\ \rm dB$ per octave increase. That's what confused me.

I've read lots of questions and answers, most of them are explaining how to calculate the value of $\alpha$ or the $\mathcal Z$-transformation. But What I'm confused is the connection between $\alpha$ and $6\ \rm dB/oct$.

I found this equation in the notes of Praat:

  The pre-emphasis factor $\alpha$ is computed as:
  $$
\alpha = \exp(-2\pi F \Delta t)
$$
  where the $F$ above which the spectral slope will increase by $6\ \rm dB/oct$, $\Delta t$ is the sampling period of the sound

According to above notes, it seems that the $\alpha$ only influences the "beginning frequency" of pre-emphasis but not "how much it changes".

Hope I explained myself clearly enough.


texts: Definition of difference between compression and filtering

texts:
I'd like to explain the motivation.

It sounds like a neuroscience topic but I am after the strict Signal Processing POV.

(I don't come from a Signal Processing background).

Interested to learn whether it's a rookie question in Signal Processing or  touches on a deeper paradox.

Motivation for Question
A scientist, many years ago, told me that the most advanced data compression we have is in the optic nerve. That visual data is compressed there and then bought to our brain for more processing.

(I'm guessing visual data is defined as a 3D field that maps light, based on a variety of factors e.g. intensity, wavelength, Cartesian coordinates etc).

If my dog and my cat and I  all look at my cup of coffee- presumably, we will all compress the visual data distinctly. For sure, we all use visual data very differently, which must mean we extract different meaning from it, ie. have different filters for what's relevant in the visual field (e.g. color, movement, sensitivity to small or large shapes).

Question

I am interested to learn whether, in Signal Processing, compression and filtration are separate concepts or they are strictly equivalent?

Because the more I thought about it, the less I can understand the difference between compression and filtering.

It seems like compression involves filtering automatically. (i.e. compression is unconscious filtering, whereas filtering is conscious).

Thanks!


texts:
I am doing an experiment, where I am trying to get the same results from correlating two matrices in the fourier domain and then scaling the data by some constant.

Method_1: Do FFT, Correlate the data, Do IFFT, Scale data

Method_2: Scale data, Do FFT, Correlate the data, Do IFFT

In other words, I am expecting that if I scale the data before I do correlation I should get the same results if I do correlation first, and then scale the data afterwards(provided that the scaling factor is the same). What's different in my case is that I am doing the correlation in the fourier domain instead of the time domain.

Below illustrates an example of scaling data before or after fft, which will yield the same results

data_1=[1,2,3,6;4,5,9,10;7,8,9,4;10,14,13,12];
%scaling constant
h=magic(4);

%METHOD_1

    %go to fourier domain
    method1_data1_fourier=fft2(data_1,8,8);             %zero padding

    %do ifft
    result1=ifft2(method1_data1_fourier,8,8);

    %scale
    result1=result1(1:4,1:4)./h;

%METHOD_2

    method2_data1 = data_1;

    %scale
    method2_data1 = method2_data1./h;

    %go to fourier domain
    method2_data1_fourier=fft2(method2_data1,8,8);              %zero padding

    %do ifft
    result2=ifft2(method2_data1_fourier,8,8);
    %result2=result2(1:4,1:4);

Below is a matlab example of what I am trying to do in C++ with 32 bit floating point precsion and fft length of 47x47 (padded with zeros to 94x94 before transform)

%correlating data_1 against data_2
data_1=[1,2,3,6;4,5,9,10;7,8,9,4;10,14,13,12];
data_2=[7,2;9,8;3,1];
%scaling constant
h=magic(4);

%METHOD_1

    %go to fourier domain
    method1_data1_fourier=fft2(data_1,8,8);             %zero padding
    method1_data2_fourier=fft2(data_2,8,8);

    %do correlation
    result1 = method1_data1_fourier.*conj(method1_data2_fourier);

    %do ifft
    result1=ifft2(result1,8,8);

    %scale
    result1(1:4,1:4)=result1(1:4,1:4)./h;

%METHOD_2

    method2_data1 = data_1;
    method2_data2 = data_2;

    %scale
    method2_data1 = method2_data1./h;

    %go to fourier domain
    method2_data1_fourier=fft2(method2_data1,8,8);              %zero padding
    method2_data2_fourier=fft2(method2_data2,8,8);

    %do correlation
    result2 = method2_data1_fourier.*conj(method2_data2_fourier);

    %do ifft
    result2=ifft2(result2,8,8);

Unfortunately I am not getting the results I expected to get, as in result1 and result2 do not match at all. Can anyone tell me what is wrong with the above method I am using? I believe it has something to do with the result be being shifted after the correlation in the fourier domain because indice (3,3) of result1 and result2 are the same, and the first matlab example works fine.


texts: Comparing a clean, noisy and enhanced signal

texts:
I was doing some experiments on analyzing few audio samples and I'm stuck with this. Suppose we have a clean audio signal ($y_{clean}$), a noisy version of it ($y_{noisy}$) and an enhanced version ($y_{enhanced}$), by applying some speech enhancement algorithm to $y_{noisy}$.
What would be the right way to quantify the relative increase/decrease in the quality of $y_{enhanced}$ (compared to $y_{clean}$) in say, particular sections of the audio ? I'm not interested in SNR (or is it wrong to ignore it?)
I was thinking something like computing the spectrogram for both and compute the distance between them, say using $L_{\inf}$ or $L_2$ norm. I don't know if it's even correct way to do this. Any help would be much appreciated !
(PS: I'm from CS background, so I don't have much knowledge about signal processing. Apologies if this question seems too basic)
Thanks.


texts: Slow Down Music Playing While Maintaining Frequency

texts: Playing a piece of music audio at a slower speed would lower its pitch (frequency). Is there a tool and theory to slow down the song playing while keep the frequency the same? I suppose one can do windowed Fourier transform or wavelet transform. It seems one have to either pre-select the window size or dynamically choose the wavelet basis. Is there any specific and detailed theory and application for doing that?


texts:
I was looking at the answer Whitened Matched Filter. It mentions that in case of ISI channel,  we can use $g^∗(−t)$ as matched filter. For ISI channel, $r(t)=c(t)∗s(t)+n(t)$ where $r(t)$ is the received signal, $c(t)$ is the channel response, $s(t)$ transmit signal and $n(t)$ is the noise signal. The receiver structure for ISI channel consists of matched filter $g^*(-t) = c^*(-t) * s^*(-t)$ followed by a sampler (A/D).

My question is, that if this receiver structure practically implemented or is it used for analysis only? Since in a real receiver, the matched filtering is implemented at the front-end does not include the channel response but only signal response ($g^*(-t) = s^*(-t)$). This is because we don't know the channel response before the sampler and we are not adapting analog filters to match to the channel response.
How does a practical receiver structure looks like with adaptive filters after the sampler on the digital side? What am I missing here?

Thanks
strong text


texts: What is differences between structured image and textured image

texts: I have been working on image processing. But there is one thing that has made me confused. What is difference between structured image and textured image. In my idea, both are the same. But recently, I read a paper that claimed that they are different. Somebody let me know what diffrerences they have?


texts: How to get a Standard Deviation weighted average in image?

texts:
Now I'm trying to implement a SDLWGW algorithm by the Matlab.
currently I'm referencing from
Automatic White Balancing Using Luminance Component and Standard Deviation of RGB Components
the problem is I can't understand  how to get a $SDLWA_{R}$ and the bar of $L_{red}(k)$?
Also what if the $L_{weight}(i,j)$ is $0.9$ then how to calculate the bar of $L_{red}(k)$?
I can't understand it.
I want to know that what is the different between $SD_{red}(k)$ and $SD_{red}(l)$?
$$SDLWA_{R}=\sum_{k=1}^{n}\frac{SD_{red}(k)}{\sum_{l=1}^{n}SD_{red}(l)}\times red_{i,j}(k) $$
What is the $SD_{red}(l)$? Especially, what is the $l$?
$$ \sum_{l=1}^{n}SD_{red}(l) $$
update 1
I have got a question one more,
The paper was written like this
L_weight(i,j) is a positive single-peak function (e.g
triangular function or Gaussian function) value for the
luminance value at i-th row, j-th column of the k-th block.
but I can't understand what is the mean? How can I make L_weight(i,j) ?
Would you please give me any hint please ?
update2.
I'm trying to understand a below answer.
Q1. What if I've got a this block then I think that the k is to
be 12 and p is 5, and q is 5. am I correct?
Q2. How to get a $\overline{L\_red(k)}$ value as below the case?
Q2-1. Especially, I want to know that what does exactly indicate the $L\_weight(1,1)red_{1,1}(k)$ and $L\_weight(1,1)$, as below the case?
I think we know the value of the $red_{1,1}(k)$ But I'm confused that how to know the value of $L\_weight(1,1)$ ?
[![enter image description here][2]][2]
update 3
Does anyone know the relationship between Gaussian weight and one Block ?


texts: Finding the correct $x[n]$ from a given DFT $x[k]$

texts:
I'm struggling to solve this question from an exam. The question gives you a set of IDFTs from a given DFT, only one of those is the correct.

The first thing I did is notice that $2\sin(\pi k)$ will always be zero for every $k$ so I ignore it.

Then I proceed to do the IDFT with the rest as $x[n] = x[k]e^{2j(2 \pi /N)n}$

And there is where i get lost and confused cause i know that the cos can be expressed in the Euler form and also:

l -> lambda;   cos(wn) = l[n-m] + l[n+m] <- Modulation.

So, when i see the given answers... i think i am missing something.

Thanks for your time.


texts:
I implemented a filter similar to the one below, within a FPGA. Two substractors and an accumulator.

The recursive equation is given by:

My inputs are exponentially decaying pulses of the form  sampled every n

I want to express y[n] as a function of the pulse amplitude A (which I need to retrieve), and the constants , L, and K.

Any hints would be greatly appreciated

Edit:

I followed Andy Walls' advice and this is what I get

$$H(z) = \frac{(1-Z^{-K})(1-Z^{-L})}{1-Z^{-1}}$$
$$X(z) = \frac{A}{1-aZ^{-1}}   $$
$$with: a = e^{\frac{-nT}{\tau}}$$
$$Y(z) = H(z).X(z) = \frac{A}{1-aZ^{-1}} .  \frac{(1-Z^{-K})(1-Z^{-L})}{1-Z^{-1}} $$
Using partial fraction expansion I try to find B and C so that
$$ Y(z) = A.(\frac{B}{1-aZ^{-1}} + \frac{C}{1-Z^{-1}})$$
I get $$ B = \frac{A.(1-a^{-K})(1-a^{-L})}{(1-a^{-1})}; C=0$$
So $$Y(z) = \frac{A.(1-a^{-K})(1-a^{-L})}{(1-a^{-1})(1-aZ^{-1})}$$

I don't know if I am right at this point, but this is where I get stuck. This filter is trapezoidal, but because the input is an exponential decay, I would expect at least 2 exponential terms corresponding to the trapeze rise and fall times.

Am I on the right path?


texts: I am doing an NMR experiment by applying an rf signal to a sample which is mounted near a magnet and in the end I get an FID (free induction decay) and do Fourier Transform to this FID . So we have a signal which is a convolution of 3 things! An excitation RF pulse , a signal from the magnet and the noise which is absolutely low. I need to deconv this signal into its components to get the pure signal from the magnet and plot it! So get red the rf pulse and the noise! I wonder whether this is possible or not! Any help would be appreciated!


texts: Noise problem when transforming and quantizing difference between images

texts:
I have a program which does forward and inverse DCT and quantization and dequantization using JPEG matrix. If I apply it to a single image it works fine, but if I try to apply it to a difference between current and previous frame it gets extremely noisy. Album with examples ("high frequency components" line in descriptions means any component starting from anti-diagonal and beyond it).

What is the correct way to compress the difference between two images?

uncompressed

Compressed, high frequency components removed

Difference between this frame and previous one was coded. High frequency components kept

Difference between this frame and previous one was coded. High frequency components removed


texts: How can I prove stability of a biquad filter with non-zero initial conditions

texts:
Ok, so the situation is that I have a DFII biquad with some filter coefficients:

  \begin{align}
w[n] &= x[n] - a_1*w[n-1] - a_2*w[n-2]\\
y[n] &= b_0*w[n] + b_1*w[n-1] + b_2*w[n-2]
\end{align}

While the filter is running, I change the coefficients and I believe my filter is going unstable.  I want to build in a check that it is going to be stable.  Currently, I simply check that the poles lie inside the unit circle by checking my $a_n$ coefficients.  But it doesn't work.  I believe that the problem is that I have non-zero values in my states (i.e. $w[-1]$ and $w[-2]$) and my stability check doesn't even consider this.  I cannot reset my states while the filter is running without producing an audible click.

The questions: How do I check for stability of a filter with non-zero initial conditions?

I have a couple thoughts for solutions:

Simply take the $\mathcal Z$-tranform of the impulse response (of the feedback portion), which looks something like this:
$$
h[n] =
\begin{cases}
    0 & n < 0\\
    1 - a_1*K_1 - a_2*K_2    & n = 0\\
    0 - a_1*h[0] - a_2*K_1   & n = 1\\
    - a_1*h[n-1] - a_2*h[n-2] & n \ge 2
\end{cases}
$$
Solve the difference equation à la differential equations:
$$
h[n] + a_1*h[n-1] + a_2*h[n-2] = 0, n \ge 2
$$
Or prove that the sequence $h[n]$ (from above) converges
$$ \sum_{n=0}^{\infty}\lvert h[n]\rvert  < \infty $$

Am I on the right track? I can't seem to figure out the solution for any of these methods.


texts: Reconstruction of bandpass filtered signal from decimated version of itself

texts:
I know how to up-sample a discrete, real signal (by an integer factor n) that is band-limited to frequencies between $f_1 = 0 \mathrm{Hz}$ and and $f_2$: Just insert $n-1$ zeros between every original sample (thereby multiplying the sampling rate by $n$) and send the result through a brick wall (sinc) filter to get the interpolated values.

In my DSP book I've read about a generalization of the sampling theorem:
If I have a discrete, real signal that is band-limited to frequencies between $f_1 = f_0  k$ and and $f_2 = f_0  (k + 1)$, with $k$ even, I can just subsample it with a frequency of $2 f_0$ (i.e. decimate) it, and I will retain the complete information (basically the frequencies get shifted down to the baseband between $0 \mathrm{Hz}$ and $f_0$).

My question to you: does the following method let me correctly reconstruct the original signal (in the frequency range between $f_1$ and $f_2$ as given in the previous paragraph) from the decimated signal?

Just insert $k-1$ zeros between every original sample (thereby multiplying the sampling rate by $k$) and send the result through a brick wall filter (bandpass filter with passband $f_1$ … $f_2$) to get the interpolated values.

If this method is not correct, what needs to be changed to get the correct results?

I assume the problem could also be solved by converting to an analytic signal, upsampling in the baseband and then shifting the frequencies up by $f_1$. I want to avoid this because of the additional cpu load this would cause. Please let me know whether you think the method above would work as well, or some other ways to save some processing costs.

Looking forward to your tips,

Flo


texts: Find which groups of ants are coexisted in the location from their signal profiles using NMF

texts:
I try to solve biological problem to tell which groups of ants are co-existed at this location by their sounds in specific time frame. My chosen framework is Non-negative Matrix Factorisation (NMF) where $V \approx WH$ (dimension is $(m*n) \approx (m*k)(k*n)$ where $m$ is frequency, $k$ is type of ant, $n$ is time in seconds) using Kullback-Leibler Divergence.

My prior knowledge are:

1) I know 3 specific frequency with intensities of each type of ant sounds. For example, Ant 1: frequency 1100, 1155, 1210 with probable intensities 0.45, 0.35, 0.20. So I supervise $W$ with intensity value at frequency $m$ and fix it (not update). Each column of $W$ thus are frequency profile of each type of ant.

2) I know which kind of ants can stay together. So I group these ants in $W$, and apply group sparsity according to group position in $W$ on $H$.

3) I can predict which duration of this time frame I will hear a specific ant sound (80% accuracy). So I initialize $H$ with peak at the predicted time with all equal intensity.

The problems that I don't have techniques to solve yet are:

1) the specific frequency and intensity can be slightly shifted or disappeared.

2) I have about thousand kinds of ants for about 10 groups (at this stage). It is overdictionary because I have to include every ants that are possible into the group. I don't know can NMF handle this?

3) Initialed peak in H with the same intensity is not good because so many ants don't exist at this time frame. But I don't have prior information of which ants appear in each time frame.

4) Peak in $H$ should be Gaussian, some tailing is ok. If it is not Gaussian it may be not ant sound.

5) While each ant can be in more than one group is rare, the overlapped frequency is much more possible. Particularly, not only exactly overlapped but partially like Ant1 1100, 1155, 1210 and Ant B 1155 1210 1345 and Ant C 900 1000 1100 are often found.

Any suggestion for overdictionary of supervised NMF with slightly shift should be allowed from supervising $W$ and $H$ by prior knowledge? The algorithm should guarantee the convergence. Is this possible?


texts: Proof of fourier transformation of multiplication of two signals

texts:
I've been trying to find a proof of the following, but still I m unable to proof it, can someone help me?
$$
ℱ[x(t)g(t)] = \frac{1}{2\pi}
 [X(\omega)*G(\omega)]
$$


texts:
I have two distance maps that describe distances from blue segment and the red segment. I need to compute a single weight map that describes transition from blue segment to the red segment. The weight map should have values in the 0-1 range where 0 is one segment and 1 is the other one.

There are the distance maps with the respective segments highlighted:

Having $d_{1}$ and $d_{2}$ the distances in the first and second image, I tried the following formula for computing the transition weight map:

$$d=\frac{d_{1}-d_{2}+d_{2max}}{d_{1max}+d_{2max}}$$

where $d_{1max}$, $d_{2max}$ are maximum distances found in the first and second image.
The numerator have range $<0, d_{1max}+d_{2max}>$ so the result is in $<0,1>$ range.
This resulted in the following map:

However, there are sharp transitions between each segment and the map.

I have also tried multiplicative formula:

$$d=\bar{d_{1}}\cdot (1 - \bar{d_{2}})$$

where $\bar{d_{1}}, \bar{d_{2}}$ are distances normalized to $<0,1>$ range. This resulted in little better transition, but still suffers from the same problem:

Do you have any other ideas how to combine the distance maps to obtain transition map?


texts: Is there a known correlation between human voice (360Hz to 5kHz) and the frequencies above/below that?

texts:
Hello SignalProcessing!

Given the important frequencies of a human voice (360Hz to 5kHz) to communicate and the frequencies above and below that make the voice sound way better.

Imagine we are recording a phone call... the provider normally transmitts only the frequencies needed for moderate communication. I am interested in reconstructing the information for the lost frequencies above and below by training deep neuronal nets.
To improve the training progress it would be helpful to figure out if there is sort of a correlation.

Any reference about that topic would be great!
Thank you.


texts: why a limit on passband in spectral masks?

texts:
I was wondering why there is a usually flat limit on passband of every spectral masks.

For sidebands it makes sense but for passband not so much.


texts: Filtering a square signal with a median filtering to preserve the edges

texts:
If needed, you can find my first post for this problem here.
I am trying to clean the following signal:

As proposed in the comment, I tried to use this post proposing 2 methods: median filtering and total variation denoising. I'm a python user, and I didn't find a good implementation of the second, thus I only kept the median filtering, implemented with scipy.
from scipy.ndimage import median_filter
window_size = 10
data_filtered = median_filter(data, size=window_size, mode='nearest')

The result is, at first look, not too bad:

However, the ramp is not completely preserved. The first one is slightly shifted to the right; while the second remains at the correct position. If I zoom on the low amplitude signal, I get:

My goal is to measure the width of the square. This is directly impacting my measure. What is causing this phenomenon, why is it not symmetric and can it be fixed?
N.B: This is impacting my measure because I would like to measure the width at a us scale. Obviously, as those first measurements were done with a sampling frequency fs = 1e6, this resolution is not possible. Would this problem be solved if I increase the sampling frequency? (Max is 1e9). Which sampling frequency would you recommend?
Thank you for the guidance,
Mathieu


texts: Time-invariance, causality and stability of $h(t)$ of four given systems

texts:
Question:

  The impulse response functions of four linear systems $S_1,\ S_2,\ S_3,\ S_4$ are given respectively by

  \begin{align}
h_1(t)&=1\\
h_2(t)&=u(t)\\
h_3(t)&=\frac{u(t)}{(t+1)}\\
h_4(t)&=\exp(-3t)*u(t)
\end{align}

  where $u(t)$ is the unit step function. Which of these systems is time-invariant, causal and stable?

  A) S1
  B) S2
  C) S3
  D) S4

My Approach:

we know if $h(t)$ is not equal to $0$ for $t<0$ then the system is non-causal; so option $A$ is non-causal. We also know if the impulse response function of the system is finite then the system is stable; so all the systems are stable. Now what is the condition to check Time-invariancy??


texts: Normalized LMS with a posteriori Error and Woodburry's Matrix Inversion

texts:
I was going through this paper and the author mentioned that we can prove the following using the Matrix Inversion Lemma (AKA Woodburry's Matrix Inversion Identity): Using matrix inversion lemma we can show that  $$w(k+1) = w(k) + \mu(k)e_p(k)x(k)$$  and $$e_p(k) = d(k) - x^T(k)w(k+1)$$  where $e_p(k)$ is the a posteriori error, is equivalent to the NLMS coefficient update equations $eq.(1)$ and $eq.(2)$ which employs a priori error in its update equation.  The way I did it, I did not use the Matrix Inversion Lemma and still proved the equivalence I want to know how to do it with the Matrix Inversion Lemma (which apparently is giving me a hard time).
I would be glad if anyone could help me out here. Even a hint would do.


texts: MATLAB and Radon Transform-cannot obtain a sinogram

texts:
I'm quite new to MATLAB and I'm trying to obtain an implementation of the computed tomography without using any built-in functions such as radon() and iradon(). I guess I understood what Radon Transform does: it simply adds up projections taken at some angles (in matrix columns in MATLAB) and produces a sinogram as a result.  But I am at a loss at this point. My attempt never produces anything that is similar to a sinogram. What am I doing wrong? I cannot get ceil(128*sqrt2) as the maximum value for the detector size as well (of one unit each on the projection axis).

Side note: Well, I admit that the code is lousy (my image, too), but I tried my best in the last week to decipher MATLAB. Thank you very much in advance for any kind of help and suggestions.

A=zeros(128,128);
for k=1:128
    for n=1:128
        if k>70&&k<100&&n>70&&n<100
            A(k,n)=1;
            else if k>20&&k<35&&n>30&&n<55
            A(k,n)=0.7;
                else if k*k+n*n>100*40&&k*k+n*n<100*50
                    A(k,n)=0.8;
                    else if (k-50)^2+(n-50)^2<=400
                        A(k,n)=0.5;
                        else if (k-100)^2+(n-20)^2>=100&&(k-100)^2+(n-20)^2<=400
                                A(k,n)=0.8;
                            end
                    end
                end
            end
        end
    end
end
imshow(A)
M=size(A);
[ilength, iwidth]=size(A);
idiag=sqrt(128^2+128^2);
lengthpad=ceil(idiag-ilength)+2;
widthpad=ceil(idiag-iwidth)+2;
pad=zeros(ilength+lengthpad,iwidth+widthpad);
pad(ceil(lengthpad/2):(ceil(lengthpad/2)+ilength-1),ceil(widthpad/2):(ceil(widthpad/2)+iwidth-1)) = A;
figure, imshow(pad)
N=100;

% angles through which we will perform the projection
for j=0:N-1
            thetaj=j*(pi/N);
            costhetaj=cos(thetaj);
            sinthetaj=sin(thetaj);

    %finding the midpoints of the pixels on the image
    for iy=1:size(pad)
        for i=1:size(pad)

            y_img=iy-pad/2-0.5;
            x_img=i-pad/2-0.5;
            img=sum(costhetaj*x_img*pad(i,iy)+sinthetaj*y_img*pad(i,iy),1);
            last_img(i,iy)=sum(img);

        end
        end
end
    figure, imshow(last_img)
    figure, plot(last_img)

Edit: I added zero padding to the code, but I still don't know why I get half of the image as output.


texts:
I have a vibration signal sampled at 25.6 kHz, and I want to perform Envelope Analysis on it.

"Vibration-Based Condition Monitoring of Wind Turbines" by T. Barszcz has the following figure:

I understand why the frequency is changing in the FFT plot during the HP and LP filters, but why is the frequency range changing from (b) -> (c) [0 -> 8 kHz] when performing the abs (rectifier) step?
Cheers!

Edit:
Here is an FFT for "Gearbox High Speed shaft" from my data.


texts:
I am trying to follow this tutorial to build a flow graph with GnuRadio 3.8. However, I am unsure what settings should I use for Chunks to Symbols and get some errors when I execute the flow graph:

Traceback (most recent call last):
  File "/home/user/top_block.py", line 186, in <module>
    main()
  File "/home/user/top_block.py", line 164, in main
    tb = top_block_cls()
  File "/home/user/top_block.py", line 129, in __init__
    self.digital_chunks_to_symbols_xx_0 = digital.chunks_to_symbols_ic(1+1j, 1-1j, -1+1j, -1-1j, 1)
TypeError: make() takes from 1 to 2 positional arguments but 5 were given

>>> Done (return code 1)

My flow graph is:

My settings for Chunks to Symbols are as follows:

Input type: int

Output type: complex

Symbol table: 1+1j, 1-1j, -1+1j, -1-1j

Dimension: 1

Num ports: 1

Am I doing something wrong here? Any help is appreciated.


texts:
I would like to understand the "justification" for the bilinear transform. The basic idea as I understand it is that by integration rule of Laplace transform we have for continuous $y(t)$:

$$\mathcal{L}\{y\}(s)=\mathcal{L}\{\int_0^ty'(x)dx\}(s)=\frac{1}{s}\mathcal L\{y'\}(s)$$
$$\Leftrightarrow \frac{\mathcal L\{y\}(s)}{\mathcal L\{y'\}(s)}=\frac{1}{s}$$

Using the trapezoid rule gives the approximation ($h$ is the sampling period):
$$y(t)=y(t-h)+\int_{t-h}^ty'(x)dx\approx y(t-h)+\frac{h}{2}(y'(t)+y'(t-h))$$

And taking the $z$ transform using the approximation gives
$$\mathcal Z\{y\}(z)\approx \mathcal Z\{ y[n-1]+\frac{h}{2}(y'[n]+y'[n-1])\}(z)$$
$$\approx z^{-1}\mathcal{Z\{y\}(z)}+\frac{h}{2}(\mathcal Z\{y'\}(z)+z^{-1}\mathcal Z\{y'\}(z))$$
$$\Leftrightarrow\frac{\mathcal Z\{y\}(z)}{\mathcal Z\{y'\}(z)}\approx \frac{h}{2}(\frac{1+z^{-1}}{1-z^{-1}})$$

Putting the above together we finally end up with the estimate
$$\frac{\mathcal L\{y\}}{\mathcal L\{y'\}}(\frac{2}{h}(\frac{1-z^{-1}}{1+z^{-1}}))\approx \frac{\mathcal Z\{y\}(z)}{\mathcal Z\{y'\}(z)}$$

Question 1. How do we get from the above estimate to the desired estimate

$$\mathcal L\{y\}(\frac{2}{h}(\frac{1-z^{-1}}{1+z^{-1}}))\approx \mathcal Z\{y\}(z)$$

Question 2. It is possible to use trapezoid rule on any linear ODE directly (iterating the rule for higher order derivatives as needed). However, interestingly this gives a different answer than the bilinear transform. The resulting coefficients $a_i$ are the same, however the $b_i$s are different.

So it appears that the bilinear transform distorts the start of the impulse. For example, for 1st order Butterworth lowpass filter we have from bilinear approximation that $b_0=b_1$. As such the impulse response ramps up before ramping down, while the continuous response and trapezoidal estimate only ramps down (see example below). Why are the $b_i$s given by the bilinear transform the "best" (e.g. most widely used) estimates (note: we only consider $b_i$s up to the order of the filter)?

Example:
Consider the 1-pole Butterworth filter given by the transfer function:

$$H(s)=\frac{1}{1+s/\omega}$$

The differential  eq. is given by:
$$y=-y'/\omega$$

By trapezoid method and the above eq.:
$$y(t+h)\approx y(t)+\frac{h}{2}(y'(t+h)+y'(t))=y(t)-\frac{h\omega}{2}(y(t+h)+y(t))$$
Collecting the terms we get the difference equation
$$y[n+1]=\frac{(1-\frac{h\omega}{2})}{(1+\frac{h\omega}{2})}y[n]$$

Which is the same as with bilinear transform, but lacking the term $b_1$.


texts: How to find the minimum no of bits required to express the filter coefficients in fixed point format from Histogram plot?

texts: I have designed a filter in MATLAB and now I have to convert the filter coefficient in fixed point format from floating point format. For that my professor asked me to find histogram first and then find how many bits I required for fixed point representation. But I am not sure how to find that from histogram plot.


texts: Step-by-Step calculation of stationary wavelet transform gives shifted results?

texts:
Good morning,

I am doing an SWT by hand to understand it better, and have a couple of questions. Maybe someone here has experience with this?

I am using the book "Conceptual Wavelets" by Fugal, where they also perform an SWT by hand. This is the page I'm referring to:

Now the details coefficients we get in the book are
-80 0 0 0 80 0 0 0 0.

When I did the convolution by hand I got this same result. Matlab agrees also.

However, when I do a level 1 SWT in Python or Matlab, the details coefficients look like this:
0 0 0 56.6 0 0 0 -56.6

I understand that the numbers are different - the book uses [1 -1] as wavelet filter, and most algorithms use [1/(sqrt(2)) -1/(sqrt(2))]. (I don't know why, it's not explained in the book.)

But the Matlab result is shifted to the left! Somehow, Matlab gets rid of the delay in the result during the computation of the details coefficients. They also reduce the convolution back to the right length. In the book, all of these transformations are done at the end, to get the original data back.

So my questions are:

Is the delay always 1, for every 1D SWT?
How does Matlab calculate its result, and how can I replicate this in my calculation? Is it always a shift to the left? Which element do I drop to reduce the size of the result? A simple convolution seems to be not enough.
(And if someone happens to know a source that mentions why haar filter values are so often 1/(sqrt(2)), I'd be happy to learn about that, too.)

Thanks in advance


texts: Total Harmonic Distortion THD for the rectifier

texts:
I have a circuit with rectifier and this cause the non-linearity and consequently harmonic. I show it on the oscilloscope and then collect the signal from the oscilloscope to the computer by using the Labview. The procedure is to find the Total Harmonic Distortion. Could anyone please guide me what I have to do in order to find the THD?
I know that the f(input),f(sampling), N(number of samples)and M(number of cycles) are very important.
Thanks in advance!


texts: What does it mean to rotate an image by 0.1 degrees?

texts:
In the old Tangent Prop paper (you don't need to know this paper), the tangent vector of the rotation transformation (think of this as a vector-valued function) at some input image is approximated by (1) rotating this input image by a very small amount, (2) subtracting this rotated image from the original image and (3) dividing the result of subtraction by the degree rotated.

However, I'm not understanding what it means to rotate an image (i.e., a matrix) by a very small amount (less than one degree).

Also, how can this be implemented in Python using libraries like Python / Scipy / Skimage?

In addition, how am I supposed to know that the image has actually been rotated, since the rotation is so small?

Note: I tried to do rotation with Skimage, but it seems like Skimage allows rotations of at least 1 degree.

Since I'm new to these concepts, I wouldn't mind if someone points me to some resources or references since it might be hard to provide a succinct answer.


texts: Deriving Frequency Response for 2-pole Zero-Delay Feedback State Variable Filter

texts:
I have an existing zero-delay feedback (ZDF) 2-pole state variable filter implementation (along the lines of the theory presented in VA Filter Design by V. Zavalishin), and I wish to determine the frequency response so that I can plot it. My usual approach with other filters (e.g. biquads) has been to determine the difference equation by inspection, apply the Z-transform, and then substitute $z=e^{j\omega}$ and plot phase and magnitude against $\omega$ from $0$ to $\pi$.

However I am having difficulty determining the difference equation of this filter.

Here is a diagrammatic representation of the implementation I have:

Here's the implementation of an iteration, corresponding to a single input sample, in pseudo-code. Fs is the sample rate, F is the filter frequency of interest, R is the resonance, x is the input sample, and lp, bp, hp represent the low, band and high-pass outputs respectively:

# z1 and z2 contain state from previous iteration

k = tan(pi / Fs * F)
b = R + k
a = 1 / (1 + b * k)

hp = (x - (z2 + z1 * b)) * a
x1 = hp * k  # temporary
bp = x1 + z1
x2 = bp * k  # temporary
lp = x2 + z2

z1 = x1 + bp
z2 = x2 + lp

For analysis, I assign $v[n]$ to be the output of z1, and $w[n]$ to be the output of z2. From this I am able to determine by inspection the following equations:

$$hp[n] = x[n] - (b.v[n] + w[n])$$
$$v[n] = k.hp[n-1] + bp[n-1]$$
$$w[n] = k.bp[n-1] + lp[n-1]$$
$$bp[n] = k.hp[n] + v[n]$$
$$lp[n] = k.bp[n] + w[n]$$

I'm not sure how to handle such a system of equations algebraically. For example, if I wish to analyse the low-pass response, then I will need to determine a difference equation for $lp[n]$ in terms of $x[n]$. However I have not had success manipulating these equations by rearrangement and substitution even with the help of the Z transform. Is it possible to do this or are there too many unknown variables?

If an algebraic solution is impossible, is there an alternative method to determine the frequency response of this filter?


texts: Suppose we denote the frequency of the note $\textbf{C}_{0}$ by $f$ then we know that $f,2f,3f,...$ are all resonating frequencies. However the frequencies $2f,4f,...,2^{m}f$ belong to the same pitch class $\textbf{C}$ having the "same" sound(higher pitch) whereas the frequency $3f$ is the note $\textbf{G}_{1}$ and belongs to the pitch class $\textbf{G}$. What is the reason(physical) for this ?


texts: How to prove this theorem about the Z transform and final value theorem?

texts:
Claim:

  If $\lim_{k\rightarrow\infty} x[k]$ exists and is finite then $X(z)$, the Z-transform
  of $x[k]$, has no poles in the region $|z|>1$ and at most  1 pole at $z = 1$.

Attempt:
\begin{align*}
   X(z) &= \sum_{k\ge 0} x[k]z^{-k}\\\
   H(z)/G(z) &= \sum_{k\ge 0} x[k]z^{-k}\\\
\end{align*}

First prove that no pole can be in $|z|>1$
\begin{align*}
   H(z)/G(z) &= \sum_{k\ge 0} x[k]z^{-k}\\\
   H(z)/[(z-a)G'(z)] &=  \sum_{k\ge 0} x[k]z^{-k},  a>1\\\
   H(z)/[(z-a)G'(z)] &= x[0] + x[1]z^{-1} +x[2]z^{-2}+x[3]z^{-3}...\\\
   H(z)/G'(z)     &= (z-a)(x[0] + x[1]z^{-1} +x[2]z^{-2}+x[3]z^{-3}...)\\\
\end{align*}

Can someone help carry this further?


texts: Can Principal Component Analysis (PCA) Solve the Cocktail Party Problem?

texts:
I'm looking into the cocktail party problem and trying to figure out whether something like Principal Component Analysis is enough to separate out all the various voices at the cocktail party into its constituent sound sources.
If its not enough, why? What other techniques should be used in conjunction with it, so that I wind up with distinct signals for each cocktail party patron's voice? Spatial filtering, such as various beamforming methods, has been suggested. But in researching PCA it seems that should (possibly) be enough to "split" the total signal of the cocktail party into individual signals for each human voice attending the party. Beamforming and similiar methods appear to be filters for subsequently focusing on just one of those voices and filtering the rest out.
Can anybody with PCA experience weigh in here with whether or not PCA can solve this, or if its requires additional processing?


texts: Describing a certain quality of an image

texts:
Generally I guess I’m asking for a source of image processing terminology.

Specifically here I’m trying to find the word for the quality that would identify the main difference between the Mondriaan on the left and the Monet on the right. Graphical/strong edges vs no solid areas and no well-defined edges.

I know a good amount about image processing but this vocab part is a glaring omission for me.


texts: How to detect the singular value with signal processing technology?


texts: Estimating an Image Filter from Examples

texts:
I am sure this question must be fairly easy and could have been answered on this site but perhaps I am not finding the right keywords to search for it.

I would like to know what are the common techniques to estimate an image filter from a set of input (unfiltered) and post-filtered images. From a machine learning point of view, this seems easy to do under some assumptions and constraints. Assuming the filters are linear, we might have to estimate the parameters of a convolution kernel using an optimization algorithm such as gradient descent. A similar approach may be done for non-linear kernels assuming a suitable non-linear model. However, I have still not found lots of literature discussing this.

As such, my questions are the following: How I can perform this operation in practice, how it is typically referenced in the image processing literature, and where can I find more information about it?


texts: How can one know the polynome equation based on the output of a system?

texts:
I came across this situation in my textbook:

However I have no clue about how you can (starting from the stepresponse on the left), get the polynome equation on the right.

Could someone please explain?

EDIT: what about this case? http://postimg.org/image/726noo2k3/b8bf257


texts: Fourier transform of dirac comb with function: The scaling factor

texts:
Multiplication in the time domain corresponds to convolution in the frequency domain:

$$
f(t) \cdot x(t) \iff F(j \omega) * X( j \omega)  \tag*{No scaling factor}
$$
I know the fourier transform of the dirac comb is:
$$
\mathcal{F} \big \{ \text{III}_{T_{s}} (t) \big \} =  \omega_{s} \cdot \text{III}_{T_{s} } (j \omega)
$$
But according to oppenheim and others:

$$
\mathcal{F} \big \{ f(t) \cdot \text{III}_{T_{s}} (t) \} =  \underbrace{\dfrac{1}{T_{s}}}_{\text{scaling factor}}  \cdot \displaystyle \sum_{k \to - \infty}^{ k \to \infty} X(j( \omega - k \omega_{s} ))  \text{ (Scaling factor)}
$$
So my question is, for the above, why does the scaling factor not remain $$ \omega_{s} $$ why does it become $$ \dfrac{1}{T_{s}} $$
Which means that
$$
x(t) \cdot y(t) \iff \dfrac{1}{2 \pi}  X(j \omega) * Y(j \omega) \tag{ ??? }
$$
Why does the scaling factor not remain the sampling angular velocity?


texts: Confused in Smoothing a signal in the Frequency domain

texts:
First of all i am new in DSP, so i have no solid education in the field.

I have convert my time domain acceleration-seismic data to frequency using DFT and i am trying to smooth the DFT data, in order to have a more appropriate view of them. I didn't know how exactly to approach my problem until i found a python code in the Scipy-Cookbook that smooths signals with convolution of them with a specific window (e.g. hanning, hamming).

You can find the code here: https://scipy-cookbook.readthedocs.io/items/SignalSmooth.html

I tried this code and the results seem quite promising as you can see.

So my questions are:
1) Is this a theoretically accepted method to smooth data in DFT domain or only applies in time domain ?  (references are welcome!)
2) Are there different methods to smooth DFT data?

Thank you!


texts: Sampling a signal with varying frequency

texts:
Question
I'm trying to figure out the sampling rate for my ADC to sample essentially signal essentially of the form:
$$y(t) = \sin(\max(t, \omega_{max})\times t) + n$$
where $n$ is noise.
Context
This signal is from a motor, where the rotational speed linearly ramps up to some $\omega_{max}$. To figure out the required sampling rate, the only method I know of is to use the Nyquist rate of 2 times the maximum frequency content. The above signal can be re-written as
$y(t) = rectangularPulse(t \epsilon [0, \omega_{max}]) \times \sin(t^2) + heaviside(t-\omega_{max}) \times \sin(\omega_{max}t)$
This isn't bandlimited since a rectangular pulse is a $\text{sinc}$ in the frequency domain. What I want to do is low pass filter the signal first as described here. But I'm at a loss as to figure out what I should set my cut off frequency to, as the fourier transform of the signal doesn't die off at any frequency (unless I'm mistaken).
I feel like this should be a classical problem, but I can't find anything online. If you know of how to solve this problem, or you could point me to where to look, it would be greatly appreciated!


texts:
I have been trying to calculate the inverse of a sweep (by inverse I mean something that when convolved with my original signal will yield a unit impulse) and although I managed to find a time-domain solution, in the process of trying to do it in the frequency domain I encountered an issue which puzzles me a bit.

So, in order to calculate the inverse of the signal I tried the same thing in two slightly different ways. The first was to synthesize the inverse filter by setting its magnitude and phase responses manually and then taking it back to the time domain. I calculated the magnitude as the reciprocal of the original signals magnitude and the phase as the negative of the original signal's phase.

The next try was to calculate the complex spectrum of the original signal, calculate the reciprocal of that and then take this back to the time-domain.

The plots below show the results. I have also plotted the differences between the two approaches (the time domain difference doesn't provide much info, it is just there for visualization purpose).

So, now my question is, why do I get such a different time-domain signals? I know that the correct is the second on but, assuming that most of the residuals originate from round-off errors of the implementation (I used MATLAB's build-in functions for that), I can't think of a good reason.

Any insights would be most welcome here.

UPDATE:
For completeness I also provide the MATLAB script I used to calculate the signals.

%% Create signal and calculate the inverses
time = 0:1/1000:10; % Create a time vector
sweep = chirp(time, 1, 10, 500, 'logarithmic'); % Create the sweep

% Artificial inverse filter
artMag = 1./(abs(fft(sweep))); % Magnitude of filter
artPhase = -angle(fft(sweep)); % Phase of filter
artInv = pol2cart(artPhase, artMag); % Convert to Cartesian to IFFT it
artInv = ifft(artInv); % Time-domain signal
artInv = ifftshift(artInv); % Shift the signal to centre it

% Original inverse filter
invSig = 1./fft(sweep); % Calculate the inverse filter
invSigMag = abs(invSig); % Calculate its magnitude
invSigPhase = angle(invSig); % Calculate its phase
invTime = ifft(invSig); % Time-domain signal

%% Plot results
figure(1)
subplot(3, 3, 1)
plot(artMag); grid on
title('Artificial inverse magnitude', 'FontSize', 16)
grid on

subplot(3, 3, 2)
plot(artPhase); grid on
title('Artificial negated phase', 'FontSize', 16)

subplot(3, 3, 3)
plot(artInv); grid on
title('Artificial inverse IFFT', 'FontSize', 16)

subplot(3, 3, 4)
plot(invSigMag); grid on
title('Reciprocal signal magnitude', 'FontSize', 16)

subplot(3, 3, 5)
plot(invSigPhase); grid on
title('Reciprocal signal phase', 'FontSize', 16)

subplot(3, 3, 6)
plot(invTime); grid on
title('Reciprocal signal IFFT', 'FontSize', 16)

subplot(3, 3, 7); grid on
plot(artMag - invSigMag)
title('Magnitude residuals', 'FontSize', 16)

subplot(3, 3, 8); grid on
plot(artPhase - invSigPhase)
title('Phase residuals', 'FontSize', 16)

subplot(3, 3, 9); grid on
plot(artInv - invTime)
title('Residuals in time domain', 'FontSize', 16)


texts:
I am trying to create an adaptive noise canceller using the RLS algorithm. The dsp toolbox from matlab offers the RLS adaptive filter already implemented, so this saved me some trouble.

My goal is to filter out the hearth beat signal from the muscle signal, however, so far I have had zero to no success. After reading a decent amount of articles on the internet, I found out that the adaptive filters are suitable for this purpose as I have a reference signal of the electrocardiogram, which I want to remove. I chose the RLS algorithm, as it converges faster than the LMS algorithm, so it should be the better option. The filter order is set to 60 and the ForgettingFactor is 0.999, which according to multiple sources, should result the removal of the noise signal, however, in my case the hearth beat signal is still significantly noticeable.

From the attached images, you can see that the signal is not filtered. Should I process the signal before using the adaptive filter or is the method wrong (or badly implemented) and I should try something else?


texts: Is it possible to observe more area in azimuth than physically traveled in stripmap SAR?

texts:
I have been examining stripmap SAR algorithms, mainly Range Migration Algorithm. When I check the textbook examples so far I always saw the data trajectory is larger than scene that is observed. An example from Melvin, Scheer et. al page 309.

As you can see, the platform trajectory roughly travels -20 to 20 and scene is half of it. Similarly, if you check this RMA implementation, this RMA Matlab implementation from the book "Spotlight SAR" by Carrara et. al., you can see the along-track axis of the reconstructed image is actually some subset of trajectory points. Just like the case above.

My question is, what happens in the opposite scenario? For instance in image above, we travel a small distance [-5 to 5] and try to produce, maybe not so precise, the image of the scene [-10 10]. If our antenna does not have high directivity, we should also able to get reflections from targets that are not directly in front of the antenna, but how we should map the frequencies in RMA?

I know that in Stolt mapping of RMA, we shift and strech each azimuth column, but in that case I should retrieve more azimuth columns than I have.. I don't know, if there is a such imaging type in SAR literature, please refer me to correct terminology, I am new to that area so I might be missing correct nomenclature.


texts: Image registration / fusion optimal similarity metrics

texts:
I have a theoretical question about optimal similarity metrics for comparing data sets.

In reading this linked paper, pp. 488-489 [1], I read the following 2 interesting statements. On page 488:

  In the ... scenario in which $A$ and $B$ differ only by Gaussian noise, then it can be shown that $\rm SSD$ [Sum of Squares Difference] is the optimum measure

And on page 489:

  If the intensities in images $A$ and $B$ are linearly related, then the
  correlation coefficient $\rm CC$ can be shown to be the ideal similarity
  measure

Are there are other conditions under which other similarity metrics [e.g. Mutual Information] can be proven to be optimal, or under which the above 2 metrics ($\rm SSD$ and $\rm CC$) are also optimal? The linked PDF focuses on image registration. I am thinking mostly in terms of image registration/image fusion similarity metrics, but would also be interested in more general similarity between any 2 or more datasets (don't have to be images).

[1] P. A. Viola, Alignment by maximization of mutual information. Ph.D. thesis, Massachusetts Institute of Technology, 1995.


texts: Functional analysis for signal processing engineers

texts:
It seems that the most advanced pure mathematics course most EE engineers take is Fourier analysis, and after that it's basically 'applied' courses. There's probably a good reason for this, but I'm not sure what it is. Functional analysis seems like it pops up all the time in signal processing.

Would exploring functional analysis be worthwhile for an aspiring communications engineer?


texts: How to convert the mean and variance of a processed received signal into a SNR or BER?

texts:
So I'm trying to eventually get some sort of signal quality metric, and to do that, I'm trying to take the mean and variance of my bit determination signal and use that to come up with a Signal to Noise Ratio, but I think I'm missing something. What I'd like to do is take my calculated SNR value, determine the Probability of error, $P_e$,  from this, and then transmit a bunch of known data, and at the receiver see how many bits actually are received in error to see if it matches with my calculated $P_e$.

What's going on is when we get our signal, we multiply it with the previous bit to get a correlation, and from that we sum it with the correlation from the previous bit to get our detection signal. This detection signal is what I'm taking the mean and variance of. If the mid-bit point of this detection signal is above a threshold (0 for now) then we say the bit is a 1, if it's below the threshold, we say it's a 0. I'm taking the mean and variance of all the 1's and the mean and variance of all the 0's, and am trying to come up with an average SNR type quantity. My thought was I could just say $\mu = Mean$, $\sigma ^2 - variance$ so say my $SNR = \frac{\mu}{\sigma^2} $ since the mean of this decision signal should be directly related to my signal power, and the variance of these points (which should relate to the spread of my decision points for 1's and 0's) will be directly related to my noise.

When transmitting Binary symbols, I know that if I increase my data rate, then the power per bit will decrease, and I see this when calculating $\mu$ and $\sigma^2$. For instance the $\mu$ is cut in half if I double the bit rate, but $\sigma^2$ stays the same, which I'm interpreting as the power/bit is being halved, but the noise is staying the same so the mean should be half and the variance of the noise should be roughly constant.

Our received signal is received in Volts, and then put through and ADC so I realize there is probably something to take into consideration there... I'm just not quite sure what yet. Am I thinking about all of this correctly?


texts:
We're doing a project in which we're sending an OFDM frame, modulated to some appropriate carrier, over a channel. Before modulating, we are instructed to upsample and low-pass it. Is there a good explanation this should be done? We have full control over sample rates at both Tx and Rx.

I guess it has something to do with an interpolated signal resembling a smooth signal more, which makes it somehow more appropriate for multiplying by a modulating function. But this is very imprecise and hand-wavy.

Can anyone explain why one would want to upsample, given that there are no differences in sample rate between transmitter/receiver?


texts:
I am new to DSP, and I am trying to recreate some analysis of a signal,

They use a chebyshev I filter to go from a 100Hz signal to a 50Hz signal. Whereas I am currently using a chebyshev II filter, however I have very little understanding of how they work, and the differences between them, and even if they do the same thing.

So my question is, firstly, how exactly does the chebyshev I filter decimate my signal?. And secondly what is the difference between that process and the chebyshev II filter and are they interchangeable? (my background is in maths if it affects the explanation)


texts: How to convert a sound signal function to wave file?

texts:
I am new to understanding sound signals. A sound signal, as far as I know, is a real function of time $x:[0,\infty)\to \Bbb R$. For example, $x(t)=\sin(t)$ is a sinusoidal signal.

What is $x(t)$. For example, what are the numbers $x(1)$ or $x(2)$? What do these numbers represent?
How can I produce a sound file e.g. sound.wav from a signal function, e.g. $x(t)=\sin(t)$. Can you give me a simple code (in c#, c++,..) for converting $x(t)$ to a sound file? (I need the code just to see and understand the details).


texts:
I wonder what is bit reversal for mixed radix FFT. Is there any algorithm that compute the bit reversal for various mixed radix FFT?
Here should be mention why I want that. For example if I have signal with N=48 points. There is some method (which I still understand very little) to make DFT by prepare matrix 2x2x2x2x3, and make 5 dimensional DFT calculations. And because I can’t understand that idea but I understand radix-2 and radix-4 I wonder if I can make the matrix like 3x16 and on each 16 point row make the radix-2 or radix for. But not sure how to prepare bit reversal for input data (or for output if I make decimation in frequency).


texts: statsmodels STL not finding all frequencies

texts:
I am trying to do a signal decomposition with statsmodels.tsa.seasonal's STL to find a seasonality. I have 3000 datapoints spanning over several years. so I know there are monthly, yearly and weekly seasonalities in my data.
When using STL it only finds the peaks in the below fft plot, i.e. periods of ~7,3.5,2.5
Is there anyone who can help me understand why I am not seiing several more seasonalities?

EDIT:
If I use two superpositioned sine series as input data, the STL only finds one of them as seasonalities:


texts: To compute fft with FFTW algorithms I think that it needs to factorizes the array size in prime factors, but it is an hard problem. which Algorithm uses FFTW to do it?


texts: Is it better to use sampling rates which are a power of 2? 32-64-128-256-512Hz? instead of using for example 30- 60-100-250-500Hz?

texts: I am recently working with a sensor that I its sampling frequency can only be configured to either 32hz- 64hz- 128hz-256hz... and I was wondering whether using such sampling frequencies have any benefits. For example is it better to use 128 Hz than use 100 Hz? or 64Hz rather than 60 Hz?


texts: When I modulate a signal $x(t)$ with $\cos(2 \pi f t)$ and the modulated signal passes through a HPF, what output do i get in the frequency domain?


texts: What should be the output of the pdf for multivariate normal distrubition?

texts:
Let's think we have a multivariate normal distribution with
$$\mathbf x_1=\begin{bmatrix}x_{11}\\x_{12}\end{bmatrix},\quad \mathbf x_2=\begin{bmatrix}x_{21}\\x_{22}\end{bmatrix},\quad\mathbf m=\begin{bmatrix}m_1\\m_2\end{bmatrix}, \quad\text{and}\quad\mathbf S.$$
Where $\mathbf m$ is the mean, and $\mathbf S$ is a $2\times 2$ covariance matrix.

What should the result of the PDF be? Should it be a $2\times 1$ array or a single value?
Can  any of its values be bigger than $1$?


texts:
I have a basic question about sweep times and number of points.

My frequency span of interest is: $f_{span} = 30\ \rm MHz$
My resolution bandwidth is: $RBW = 10\ \rm kHz$
The minimum number of points is: $n = 2\frac{f_{span}}{RBW} =  6000\ \rm points$
My minimum dwell time per point is: $T_s = 15\ \rm msec$

My scan time for this measurement would be:
$$t_1 = n T_s = 6000 \cdot 15\ \rm msec = 90\ \rm sec$$

If I exceed the number of points in the required span by 2, can I decrease my required dwell time by 2?

$$t_2 = 90\ \rm sec = 12,000 \cdot 7.5\ \rm msec$$

Can these 2 measurements on a Spectrum Analyzer be considered "equivalent"?

What are the trade-offs?

To give some background, I am using a software tool that fixes the number of points in a span to default value. Since I cannot set number of points to exactly what is required, I would like to set the sweep times to the smallest allowed time. Is my reasoning correct?


texts: Comparing MFCC Features ,What do they represent?

texts:
I know that MFCC features are the spectral envelope of the input signal but
I can't understand what do they mean and what do they represent  . and if I have two people saying the same word how can I compare the resulting features.

For example I have Mel Frequency Cepstral Coefficient  for the word "please":
first person:

second person:

I see these signals seem different but thay are representing the same word which is 'please' .

I've done these steps to extract the features:

framing with overlapping percentage equals to 50%.
Windowing : hamming window was used
perform  fft(signal);
take the squared magnitude of the ff(signal);
apply mel filter bank :number of filters = 26;
take the logarithm of  squared sum for each filter
perform dct .


texts:
I'm trying to understand the derivation of LS algorithm for channel estimation in OFDM.
Could anyone explain why do we need complex conj transpose in the below equation?


texts: rate 1/128 channel code for 5 bits of information

texts: I need to construct channel code for 5 bits of information for a wireless channel. Possible code rate could be as low as 1/128 (very reliable). At decoder I assume to use channel LLRs as decoder input to take into account channel conditions. At a first glance I'd use a repetition code (or spreading sequence) thus there will be 2^5 coded word in total and use ML decoder (correlation with each code word). These seems to be simple in realization. Are there any other types of channel codes for such requirements.


texts: Detect ending of frame in H.264 video stream

texts:
I am using JM reference software for encoding and decoding of video. What I want to do is to add Redundant bits using RS codes to the encoded video so that if some bits get lost during transmission we recover it on decoder side.

The problem I am getting is to detect end of frame when reading H.264 file.
So my questions are as follows

How to detect end of frame when reading h.264 file?
Also, how to detect end of GOP in h.264 stream?


texts: How to take the FFT of an undersampled signal

texts:
I am sampling a 5mhz blocks at an IF of 40 MHZ. I am sampling at 11MHz per bandpass sampling theory, but when I run the FFT I get all the aliases.

So I am trying to figure out how to constrain the fft interpretation to the appropriate bands.

IE i may be sampling the band from 80-85MHz by mixing it down to 37.5-42.5 MHZ and sampling at 11 MHZ. How do I take the FFT and plot it showing the peaks appropriately?


texts: What is limit on performance of edge-preserving filters for low signal-to-noise regimes?

texts:
Let's say I have a signal, $x(t)$, defined as such,

$
x(t) =
\begin{cases}
 0 &\text{if} \,\,\, t < -\alpha/2 \\
 \frac 1\alpha t+\frac12 &\text{if} \,\,\, -\alpha/2 \leq t \leq \alpha/2 \\
 1 &\text{if} \,\,\, t > \alpha/2
\end{cases}
$

Thus, $x(t)$ is some signal that exhibits a jump from $0$ to $1$ over some timescale $\alpha$ around $t = 0$. For $\alpha = 0$, $x(t)$ has a discontinuous jump or an edge.

Now, let's add Gaussian Noise with $\sigma = 1$ to $x(t)$. Thus, $x(t)$ has a signal-to-noise of 1.

Now, say, I want to use some sort of edge-preserving smoothing filter on the noisy $x(t)$ to determine $\alpha$. It is obvious that any smoothing filter will wash away this edge a little and make it difficult to distinguish between a signal with $\alpha = 0$ and a signal with $\alpha = \alpha_{\mathrm{min}}$ for some minimum alpha. But there should be some minimum $\alpha$ that is distinguishable from a sharp edge.

For example, a Gaussian smoothing filter of some scale will make it impossible to distinguish between $\alpha$ less than that scale, but should make it possible to detect $\alpha$ larger than the scale of the filter.

What is the best possible performance a smoothing filter can have with regards to edge-preservation in this case? That is, what is the best possible smoothing filter that I can use for this low signal-to-noise regime that will be able to distinguish between $\alpha = 0$ and $\alpha = \alpha_{\mathrm{min}}$ for the smallest possible $\alpha_{\mathrm{min}}$.

Now, obviously, in practice we have discrete signals and $x(t)$ will be sampled at some sampling rate. So, hopefully, the answer to this question will address this. There's no real use to answering this question as if it was referring to a continuous signal.


texts: In case of Complex DFT spectrum, why the x axis range from N/2 to N point mean a negative frequency?

texts:
When we transform  a complex signal into frequency signal by using a complex DFT,

The range from N/2 to N point on the X axis of the spectrum mean a negative frequency...

But i cant understand why it is..

In addition, we should not use the range over the N/2 point because of the aliasing issue.  But i heard that the negative frequency region is used in case of complex dft spectrum of a complex signal.

Very confused..

Could you clearly explain it for me?

Thank you.


texts: Noise Reduction: How can I filter out saw or other tool noises?

texts:
I recently came across old recording from the 1930's. Unfortunately not only is the quality low, it's a bunch of carpenters talking as they saw things and make all sorts of other noises.

I've tried:

Various filters
Audacity's native Noise Reducer

and yet no luck, I find that either what I'm doing is not working, or when I get the speech to sound louder than the tools the amplitude is too low to make out anything.

Is there any well defined technique for this?


texts: Python audio analysis: which spectrogram should I use and why?

texts:
I am doing my final project at university: pitch estimation from song recording using convolutional neural network (CNN). I want to retrieve pitches existed in a song recording. For CNN input, I am using a spectrogram.
I am using MIR-QBSH dataset with pitch vectors as data label. Before processing the audio to CNN (each audio has 8 sec duration in .wav files of 8 KHz, 8 bit, mono), I need to pre-process the audio into a spectrogram representation.
I have found 3 ways to generate a spectrogram, the code are listed below.
Audio example I am using in this code is available here.
Imports:
import librosa
import numpy as np
import matplotlib.pyplot as plt
import librosa.display
from numpy.fft import *
import math
import wave
import struct
from scipy.io import wavfile

Spectrogram A
x, sr = librosa.load('audio/00020_2003_person1.wav', sr=None)

window_size = 1024
hop_length = 512
n_mels = 128
time_steps = 384

window = np.hanning(window_size)
stft= librosa.core.spectrum.stft(x, n_fft = window_size, hop_length = hop_length, window=window)
out = 2 * np.abs(stft) / np.sum(window)

plt.figure(figsize=(12, 4))
ax = plt.axes()
plt.set_cmap('hot')
librosa.display.specshow(librosa.amplitude_to_db(out, ref=np.max), y_axis='log', x_axis='time',sr=sr)
plt.savefig('spectrogramA.png', bbox_inches='tight', transparent=True, pad_inches=0.0 )

Spectrogram B
x, sr = librosa.load('audio/00020_2003_person1.wav', sr=None)
X = librosa.stft(x)
Xdb = librosa.amplitude_to_db(abs(X))
# plt.figure(figsize=(14, 5))
librosa.display.specshow(Xdb, sr=sr, x_axis='time', y_axis='hz')

Spectrogram C
# Read the wav file (mono)
samplingFrequency, signalData = wavfile.read('audio/00020_2003_person1.wav')
print(samplingFrequency)
print(signalData)

# Plot the signal read from wav file
plt.subplot(111)
plt.specgram(signalData,Fs=samplingFrequency)
plt.xlabel('Time')
plt.ylabel('Frequency')

Spectrogram results are displayed below:

My question is, from the 3 spectrograms I have listed above, which spectrogram is best to use for input to CNN and why should I use that spectrogram type? I am currently having difficulty to find their differences, as well as their pros and cons.


texts: Getting blobs from similar colored contours in grayscale image

texts:
I want to select the staff lines as individual blobs in either of these images:

Is this possible?


texts: How can I decompose a signal into square waves?

texts:
I am dealing with signals that are a superposition of different square waves with different amplitudes and phases. Normally, one would decompose a signal into sine waves with help of the Fourier transform, but in this particular case a decomposition into square waves would be much more effective. A Fourier transform would produce a very complicated spectrum, while a square wave decomposition should give just a few clear lines.

I know that such a decomposition is possible. In fact, I could use any periodic function as a basis for the decomposition and this is mentioned in many texts on the subject. But I could never find a formula or an explicit example for a decomposition into a non-sinusoidal basis.

My approach to decompose a signal consisting of the $N$ samples $x_k$, was to use a DFT-like formula
$$ u_k = \sum_{n=0}^{N-1} x_n \, \mathcal{R}_k(n)$$
where $\mathcal{R}_k$ is a real-valued square wave with a frequency $k$ times the base frequency. But this is certainly not complete, since I do not obtain any phase information for the constituent square waves, and I couldn't invert the procedure.

How can I decompose my signals into square waves with well-defined amplitude and phase?


texts: I'm looking for a way to shield or protect produced signal or wave. In mechanics suppose you run the wave in tunnel. Is there a way to do it in transmitting the signals or waves created by mobile phone or wifi or in the internet?


texts: How range of the signal increases when modulation is done?

texts: I am very confused because we know by increasing frequency we decrease the wavelength hence the distance covered will be very less but when we modulate a signal we are basically increasing the signal frequency which increases the range of the signal  , so how is that possible??


texts:
I have a curve like this one:
I need a function to approximate this curve, but I need that the function be a low order function (less than 5). What is a good way to obtain what I expect? Thanks!