texts: Choosing correct filter parameters for IMU sensor datas

texts:
I have a IMU sensor breakoutboard using ITG3701, LSM303D LinkProduct. This involves accelerometer, gyroscope and magnetometer. I have connected the sensors to a Arduino and sending the datas to another computer over Xbee. My sensor is placed on a wheel along its radius. The controller sends data at every 100 ms delay.

The frequency parameters for sensors , i have chosen to be

// Specify sensor full scale
uint8_t OSR = ADC_8192;  // set pressure amd temperature oversample rate
uint8_t Gscale = GFS_4000DPS; // gyro full scale
uint8_t Godr = GODR_190Hz;  // gyro data sample rate
uint8_t Gbw = GBW_low;  // gyro data bandwidth
uint8_t Ascale = AFS_16G; // accel full scale
uint8_t Aodr = AODR_200Hz;  // accel data sample rate
uint8_t Abw = ABW_50Hz;  // accel data bandwidth
uint8_t Mscale = MFS_2G; // mag full scale
uint8_t Modr = MODR_25Hz; // mag data sample rate
uint8_t Mres = MRES_HighResolution; // magnetometer operation mode

It can be changed to several other allowable limits given here

I want to process the sensor datas through two filters, first HPF (to remove DC components noise) and then LPF (to remove high frequency noise)

I want to do the following

Accelerometer_Raw ->  /HPF/  -> Accelerometer_HPF_datas -> /LPF/ -> Accelerometer_LPF_datas

Similarly

Gyroscope_Raw ->  /HPF/  -> Gyroscope_HPF_datas -> /LPF/ -> Gyroscope_LPF_datas

On the receiver side, i did the following on the sensor datas

   FSa=200; % accelerometer sampling
        FSg=190; %gyroscope sampling

d = fdesign.highpass('N,F3dB',2,3/(FSa/2));
            H_a1 = design(d,'butter');
            H_a1.PersistentMemory=true;
            H_a2 = design(d,'butter');
            H_a2.PersistentMemory=true;
            H_a3 = design(d,'butter');
            H_a3.PersistentMemory=true;
            d = fdesign.lowpass('N,F3dB',2,30/(FSa/2));
            L_a1 = design(d,'butter');
            L_a1.PersistentMemory=true;
            L_a2 = design(d,'butter');
            L_a2.PersistentMemory=true;
            L_a3 = design(d,'butter');
            L_a3.PersistentMemory=true;
            d = fdesign.highpass('N,F3dB',2,3/(FSg/2));
            H_g1 = design(d,'butter');
            H_g1.PersistentMemory=true;
            H_g2 = design(d,'butter');
            H_g2.PersistentMemory=true;
            H_g3 = design(d,'butter');
            H_g3.PersistentMemory=true;
            d = fdesign.lowpass('N,F3dB',2,30/(FSg/2));
            L_g1 = design(d,'butter');
            L_g1.PersistentMemory=true;
            L_g2 = design(d,'butter');
            L_g2.PersistentMemory=true;
            L_g3 = design(d,'butter');
            L_g3.PersistentMemory=true;

Gyroscope_Z axis in degrees per second

Accelerometer_X axis in g

I chosed the cutoff frequncy to be HPF to be 3Hz and for LPF to be 30 Hz. I found there is a very big difference with the gyroscope datas after filtering when the wheel rotates. The output datas are significantly reduced for gyroscope.There are also difference at the peaks of accelerometer datas measured. Hence, when i input this to Kalman, I am getting  lower velocty and the lower angle rotated. I have also verified that the output is wrong, as rotating a complete rotation gives me less than 2Pi radians after using filtering.

My questions

If my method is correct for filtering and choosing the sampling rate, cutoff frequency, order of filter etc.
If not, what filter parameters should I chose instead to get better results?
If HPF for gyroscope is necessary?

Additional:

After suggested the lower frequency limit should be much lower: I did the following High pass filtering at 0.05 Hz, The accelerometer datas are good, however i am somehow filtering out the content from Gyroscope datas

The lowpass filter at 3 HZ gave me much better result.

PS: I found the these 3 and 0.05 HZ from the FFT analysis


texts: Scaling the input vs scaling the impulse response for an LTI system

texts:
Two different cases:

We pass $x(t)$ to an LTI system with impulse response $h(2t)$ and get the output $y(t)$.

We pass $x(2t)$ to an LTI system with impulse response $h(t)$ and get the output $z(t)$.

I'm looking for possible restrictions (e.g. on the input / impulse response) such that the $y(t)$ and $z(t)$ are equal, if possible.
So far I have this:
$$
y(t) = x(t) * h(2t) \implies Y(j\omega) = X(j\omega) \frac{1}{2} H(j\omega/2)
$$
$$
z(t) = x(2t) * h(t) \implies Z(j\omega) = \frac{1}{2} X(j\omega/2)H(j\omega)
$$
If we take, for example, $X(j\omega)$ from the first equation and plug into the other one, we get
$$
Z(j\omega) = \frac{Y(j\omega/2)H(j\omega)}{H(j\omega/4)}
$$
I'm stuck at this point to comment on the relation between $y(t)$ and $z(t)$ and assert some restrictions to make them equal possibly.


texts: I have pressure data with respect to crank angle and I want to convert it into Sound pressure level with respect to frequency. Please suggest any known procedure .Thank you


texts: I have to find the Fourier coefficients of $$ \frac{1}{1+ t^{2}} $$ I tried with $$ \frac{1}{T}\int_{0}^{T} \frac{1}{1+it}e^{-i 2 \pi f_0 T } $$ but I should do at least two integrals by parts , so I tried to calculate first the even part and the odd part of x(t). I obtained that the even part is $$ x_p = \frac{1}{1 + t^2} $$ and the odd part is $$ x_d = \frac{-2it}{1 + t^2}  $$ now I should find the Fourier coefficients of $x_p$ and of $x_d$ ( $ X_p(k)$ and  $ X_d(k)$) and the Fourier coefficients should be the sum of these two. The only problem is that I that I still have integrals by parts.


texts: Is it necessary to downsample data for filtering of low frequencies?

texts:
My data is sampled at $f_s = 32kHz$, and for a particular analysis I'm only interested in signals with certain frequencies from $f_{highpass} = \text{0.1Hz}$  to $f_{lowpass} = \text{300 Hz}$.

When I build and apply a 4th order butterworth filter to my data I get an error and there's not much left of the signal. Here's some Matlab code:

order = 4;
[b,a] = butter(order, [0.1 300]/(32000/2), 'bandpass');
y = filtfilt(b,a,data);

The error message I get is: Warning: Matrix is close to singular or badly scaled. Results may be inaccurate. RCOND =  1.101227e-16.

I've also tried to apply a low pass filter first, then followed by a highpass filter:

FOI = [0.1 300]
SampleRate = 32000;
[b_l,a_l] = butter(order, FOI(2)/(SampleRate/2), 'low' );
[b_h,a_h] = butter(order, FOI(1)/(SampleRate/2), 'high' );
tempData = filtfilt(b_l, a_l, data);
output = filtfilt(b_h, a_h, tempData);

But I keep getting the same result.

My questions:

why does the filtering fail?
why does the filtering with order 3 work but not with order 4?
do I need to downsample the data (e.g. to 1kHz) before I can do the filtering?

EDIT: other people have hit the same problem but I'm looking for a definite answer.

EDIT2: if I change the highpass value to 1 and do the highpass filtering before the lowpass filtering, everything works. But why? Why is the cutoff of 0.1 causing problems?


texts: How to conduct spectrum analysis on a wav audio file using SciLab?

texts:
I am working on a wav audio file using SciLab. The wav file is to undergo the FFT (Fast Fourier Transform) by using the fft function to obtain the spectrum analysis of the wav file.
The wav file is sampled at 8kHz.
Here, I attach my code and the result that I obtained :
File = "E:\HDD DOCS\UPM SEM 5\DSP\Assignments\A5\message1.wav";

wavread(File, "size")
wavread(File, "info")

[y, Fs] = wavread(File); Fs

sample_rate= 8000
t=0:1/sample_rate:1 ;   //establishing period
N=length(File);

clf

z=fft(y);
f_vector=sample_rate*(0:(N/2))/N;
n=size(f_vector,'*');

clf();

subplot(2,1,1)
plot(y)
xlabel("Time Domain");
subplot(2,1,2);
plot(f_vector,abs(z(1:n)));
xlabel("Frequeancy(Hz)");

It seems that the result of the fft is not right and the plot for wav file has some errors.
Therefore, my question is :

How to change the $x$-scale of the original wav file plot into time ?

I'm assuming that my equation for fft plotting is not correct, how should I alter the equation?

Thank you in advance.
Regards.


texts:
my question is quite straight-forward to be grasped and handled.

I have an 1D array, F_11 representing a target power spectrum. It is sketched in the following figure

Then, I start generating a time series u which embodies the properties of F_11, as follows

nn = complex(normrnd(0,1,[1,1000]),normrnd(0,1,[1,1000]));
u = abs(fftn(fftshift(F_11.*nn)));

Finally, just to check that u does actually reproduce a spectrum as F_11, I try going backwards as follows

ps =  abs(ifftn(ifftshift(u)));

But, when plotting ps, it does not exactly match F_11

In the high frequency area, ps drifts upwards, which is quite weird.

Would you mind to sketch the flaw in my procedure?


texts: Stability of a logarithmic function and system

texts:
A system such that

$$y[n] = \log_{10}(|x[n]|) $$

was given as a stable one. I am not able to understand it though. When $x[n] = 0$ (which is bounded), the output $y[n]$ tends to get unbounded. So how is this system stable?


texts:
I'm new to DSP so the question might seem a little bit stupid: When using a low pass filter (say Butterworth filter) with Python, I came across an algorithm with which I can reconstruct my signal after filtering (from scipy.lfilter). I'm just wondering how the signal is reconstructed. I thought they used Fourier transform but it's not mentioned.
Thank you!


texts: What are the main advantages and disadvantages of the lifting scheme?

texts: Can anybody tell me the advantages of the lifting scheme and why it is better then a normal filter bank?


texts: Implement the 3D Structure Tensor in C/C++?

texts:
I am trying to figure out the details on how to implement the 3D structure tensor in C/C++ in an easy but efficient way and need some advice!

For a discrete function $ I(x_i,y_j,z_k)$ the 3D structure tensor is given by:
$$ S=\begin{pmatrix}
W \ast I_x^2 & W \ast (I_xI_y) & W \ast (I_xI_z)\\
W \ast (I_xI_y) & W \ast I_y^2 & W \ast I_yI_z \\
W \ast (I_xI_z) & W \ast I_yI_z & W \ast I_z^2 \\
\end{pmatrix}$$
where W is a smoothing kernel, $\ \ast $ denotes convolution and subscript denotes partial derivative with respect to.

The calculation of the structure tensor have two main steps:

calculate the partial derivatives of the function in a way that is
robust to noise over a window.
smooth products of the partial
derivatives over another larger window.

I start by looking at step 2.
I want to use a Gaussian as smoothing kernel. The normal distribution in 3 dimensions is separable:
$$ g(x,y,z) = g(x)g(y)g(z) $$
where
$$ g(x) = \frac{1}{\sqrt{2\pi}\sigma}exp(-\frac{1}{2}(x/\sigma)^2) $$
etc.

The Fourier transform of the normal distribution in 3 dimensions is given by:
$$ G(kx,ky,kz) = G(kx)G(ky)G(kz) $$
where
$$ G(kx) = \frac{\sigma}{\sqrt{2\pi}}exp(-\frac{1}{2}(kx\sigma)^2) $$
etc.

How do I implement the Gaussian smoothing?

The simplest way to implement the Gaussian smoothing would be to loop over a 3D Gaussian kernel for each point in I.

Since the Gaussian is separable however it should be more efficent to perform a convolution with a 1D Gaussian in the x direction followed by the y direction and z direction.

Another even more efficient approach (?) would be to do the convolution in the Fourier space (where it becomes a multiplication):
$$ g*a=\mathcal{F}^{-1}[GA] $$

Next I look at step 1.

The derivative of the Gaussian is given by:
$$ G_x(x) = -x/\sigma^2G(x)G(y)G(z)$$
etc., where subscript denotes partial derivative with respect to.

The Gaussian derivative can be used to estimate the partial derivatives of I in a way that is robust to noise:
$$ I_x=-x/\sigma^2G(x)*G(y)*G(z)*I $$

Here I am faced with the same decision:
implement it without using separability,
implement it using separability or implement it in Fourier space?

Any advice or comments?


texts: Modeling audio effects, reverse engineering

texts:
Okay so I'm relatively new to signal processing and don't know the specific terminology to help narrow my search, but I'm looking for advice in the form of books, papers, references, etc that can help me solve the following problem:

Lets say a dry audio signal is passed through a digital effect or VST in a workstation.
An effect is thus applied to the input sample generating the wet output.

Now is there a general method to analyze and determine the exact function of that effect so that we could reverse engineer it in the form of some physical I/O system we build?
(ex: find a really cool audio effect that works well on guitar, sadly it exists only in the DAW environment and isnt exactly practical for live performance, reverse engineer effect it into a pedal, can now plug in any instrument.)

bonus: would using Audio Neural Style Transfer be the best option? (in which case any machine learning references that might help I'd be grateful to read!)


texts: Processing and plotting quadrature data in the frequency domain

texts:
Ok, so I recently started playing around with an SDR in an effort to challenge and educate myself. So, let's say the sampling rate (Fs) is 3MHz and of those we grab 3M samples. Now we are left with 3M samples consisting each of an IQ pair. Complex data. My understanding is that I would perform an FFT using the entire sample, consisting of IQ as the input for the FFT which would yield another complex answer. From this, calculate the magnitude and plot starting from 0? I'm not entirely positive which samples to throw out, which ones to keep, and in which order to plot them.

Do I perhaps separate the real and imaginary prior to the FFT for the first 1.5M samples, perform an FFT for 1.5M of I and FFT for 1.5M of Q and plot the real to the right of 0 and the imaginary to the left of 0? Throwing out the rest of samples from 1.5-3M?

Any help is very much appreciated. Thanks all.


texts: when exposure time is longer, the Poisson noise become Gaussian

texts:
Can we have an explanation of why exposure time is longer, the Poisson noise become Gaussian?
Or any reference? Please provide a reference and I like to dig into more.
Thanks in advance!!


texts: How can I get the polyphase components of a signal

texts:
My question, What is the polyhase components of the signal? What does mean a polyphase component for a signal?

Is it true that example :

N = 100;
h = firls(N, [0 .2 .25 1], [1  1 0  0],[.001 .0001]);

% PolyPhase Componets
%First polypase component

for i=1:length(h),
h0(i)=h(i);
end

%% Second polyphase component

for i=1:length(h),
h1(i)=h(i)*exp((-j*2*pi*((i-1)))/4);
end

%% Third polyphase component:

for i=1:length(h),
h2(i)=h(i)*exp((-j*2*pi*(2*(i-1)))/4);
end

 and so on ..

Thank you


texts:
Let's suppose that I send the following sequence of bits:

$x = [1, 0]$

I use a 7 length Barker code ($b = [−1,−1,−1,1,1,−1,1]$) in a BPSK modulation scheme, resulting in the following sequence:

$x_B =\\
[−1,−1,−1,\:\:1,\:\:1,−1,\:\:1\\
\:\:\:1,\:\:\:1,\:\:1,-1,-1,\:\:\:1,-1]$

We send it, through a AWGN channel, and after decoding there are some errors in the sequence:
$y_B =\\
[\:\:1,−1,\:\:1,\:\:1,\:\:1,−1,\:\:1\\
\:\:\:1,-1,\:\:1,-1,-1,\:\:\:1,-1]$

How would I decode this sequence?


texts: Why can adaptive IIR filters result in unstable solutions?

texts:
For adaptive filtering, both finite and infinite impulse response (FIR/IIR) filters can be utilized. As an advantage of FIR filters in this context, guaranteed stability is often mentioned, while IIR filters do not share this property (see this related question and its answer).

I understand IIR filters are unstable when their impulse response diverges, meaning a transfer function pole lies outside the unit circle (has a radius $r>1$).

My intuition is that since in the example of system identification, only stable systems are practical to identify, it would make sense that estimates of such a systems' transfer function were stable as well.

My question is:

Why is this not the case? How can the adaption of IIR filters result in unstable solutions?

Please note that I am not asking about measures to ensure stability in IIR adaptive filtering.


texts: How to detect rectangular images inserted into another image?

texts:
It's just for a simple binary classification, so it's not needed to find the location of the inserted image or the content. I'm only interested in checking whether something like that is present or not.
Positive examples:

Negative examples:

The inserted image rectangles are always upright, i.e., they are not rotated by non-90-degree steps. They cover at least 1 % of the image area, i.e., their edges are at least 10 % of the image dimensions.
The images have undergone lossy compression to an unknown degree.

I've tried looking for straight horizontal and vertical edges using Scharr filters, Canny edge detection, and Hough lines, but the results are not very reliable. Any creative ideas would be appreciated. :)

My current attempt looks as follows:
import numpy as np
from cv2 import cv2

def get_value(image: np.ndarray) -> float:
    assert len(image.shape) and image.shape[2] == 3, "please provide a color image"

    scharr_scale = 1 / 16

    # Detect vertical and horizontal edges on all three channels (BGR) separately.
    scharr_x_abs = np.abs(cv2.Scharr(image, cv2.CV_32F, dx=1, dy=0, scale=scharr_scale))
    scharr_y_abs = np.abs(cv2.Scharr(image, cv2.CV_32F, dx=0, dy=1, scale=scharr_scale))

    # Reduce edge images to one channel.
    scharr_x_abs_gray = cv2.cvtColor(scharr_x_abs, cv2.COLOR_BGR2GRAY)
    scharr_y_abs_gray = cv2.cvtColor(scharr_y_abs, cv2.COLOR_BGR2GRAY)

    # Erode the edge images perpendicular to the filter direction to remove noise.
    scharr_x_abs_eroded = cv2.erode(scharr_x_abs_gray, np.ones((5, 1)))
    scharr_y_abs_eroded = cv2.erode(scharr_y_abs_gray, np.ones((5, 1)))

    # Get the average edge value per column/row part.
    scharr_x_row_abs = cv2.resize(scharr_x_abs_eroded, (scharr_x_abs_eroded.shape[1], 8), interpolation=cv2.INTER_AREA)
    scharr_y_col_abs = cv2.resize(scharr_y_abs_eroded, (8, scharr_x_abs_eroded.shape[0]), interpolation=cv2.INTER_AREA)

    # Build second derivative, but perpendicular to the original one,
    # to detect abrupt starts and ends of edges.
    scharr_x_row_abs_scharred_again = np.abs(cv2.Scharr(scharr_x_row_abs, cv2.CV_32F, dx=1, dy=0, scale=scharr_scale))
    scharr_y_col_abs_scharred_again = np.abs(cv2.Scharr(scharr_y_col_abs, cv2.CV_32F, dx=0, dy=1, scale=scharr_scale))

    # The maximum abrupt start/end of edges are the final result.
    scharr_x_row_max = np.max(scharr_x_row_abs_scharred_again)
    scharr_y_col_max = np.max(scharr_y_col_abs_scharred_again)
    result = float(np.maximum(scharr_x_row_max, scharr_y_col_max))

    return result

Output for the provides examples images:

100.99776458740234: A high value for this positive example. That's good.
37.71356201171875: A low value for this positive example. That's not good.
56.305137634277344: A not-too-high value for this positive example. That's ok.
27.330982208251953: A not-too-high value for this positive example. That's ok.

So there is no threshold that would divide the two classes correctly.


texts: Are there algorithms to invert a magnitude spectrum?

texts:
I have a speech signal and have calculated the magnitude spectrum for the signal. Is it possible to get a good estimate of the original signal from the magnitude spectrum, seeing that phase information has been discarded after the FFT?

I have looked at the Wiener-Khintchine theorem, but I am unable to find any code examples of it in use, and I have not managed to get it to work.


texts: Why Does the Kalman Filter Remove Only Gaussian Noise?

texts: What and where in the derivation of the Kalman filter is the assumption of Gaussian noise? Why and how does this assumption help?


texts: How can I obtain a Zero Phase Moving Average Filter? I read for example in Matlab that filtfilt give you zero phase doing forward and backward filter, I dont understand well how that work, I think taking the same number of values in the past and future can give me that, but it isnt causal right?


texts: Do we have to insert the DC null in OFDM when using RF direct sampling ADC architecture?

texts:
Generally, DC Nulls are added in OFDM because zero-IF receivers don't pass the DC component. But what if I use an RF direct sampling architecture? Will there be DC offset? If not, then is it okay not to insert a DC null?
Thank you in advance!


texts: If a Kalman filter can only receive information on $(x, y)$ position, is there any reason to have acceleration as part of the model?

texts: I'm new to working with Kalman filters. I have a Kalman filter that's working well for modelling the 2D position of an object. The model at the moment holds information on position ($x$ and $y$) and velocity ($x$ and $y$). I only have position information available for measurement. I was wondering if adding acceleration to the model would improves things in any way, where the only measurement information available is still just the $x$ and $y$ position. Would adding acceleration to the model make it so sudden jerks in movement are handled better or would there be no difference at all?


texts:
I am looking to some C algorithm(for MCU), where I have these signals :

First graph you can see that there is an event in the middle and it goes back to normal, second and third there is some event at a certain point.
Any of them can happen per session(= when I turn the device there will be 1 behavior)

I have to find  real time when the event begin and ends, and it has to be FAST.

All maximums in all 3 graphs are more or less the same (with 2% diff)
I don't know the values before hand but I can calibrate the "no event" state.
I need 1 algorithm to all of them.

I can do all sorts of basic stuff with averages and LPF, but I guess there is more clever and generic way to do so (?)

Something that can cover all cases and do what my eyes do - I can just tell where is the event.


texts: I think we don't have to do carrier phase recovery. Do i need to do carrier frequency recovery. My symbol rate is high around 2 MHz  so i am thinking carrier freq offset is less compared the symbol rate. Or does my gardner timing recovery can take care of this carrier freq offset recovery. I have another question should i do differential decoding then timing recovery or timing recovery followed by differential decoding . Because the pi/4 dqps will have valid symbols at zero level also ie symbol with zero I data and Q non zero data. Similary zero Q data and I non zero data. So i am thinking if do timing recovery first i may get proper output as garder looks for zeros crossing between symbols. Is my explanation right


texts: The essential bandwidth of a rectangular pulse

texts:
The problem is to estimate the essential bandwidth of a rectangular pulse

\begin{equation}
g(t) = \Pi(t/T),
\end{equation}

Where the essential bandwidth must contain at least $90\%$ of the pulse energy.

Here's what I'm trying to solve :

$$0.9 =\frac1\pi\int\limits_0^{2BT\pi} \mathrm{sinc}^2 \left(\frac x2\right) \,dx$$

I'm not sure how to get MATLAB to solve this numerically for the upper bound value.

first the energy for the pulse is $E_g$

\begin{equation}
E_g = \int^{+\infty}_{-\infty} g^2(t)dt\, = \int^{+T/2}_{-T/2} dt\, = T
\end{equation}

and

\begin{equation}
\ \Pi(t/T) \Longleftrightarrow T\mathrm{sinc}(\pi fT)
\end{equation}

and for the Energy Spectral Density of the pulse

\begin{equation}
\ \Psi(f) = \mid G(f) \mid ^2 = T^2\mathrm{sinc}^2(\pi fT)
\end{equation}

The energy $E_B$ within the band from $0$ to $B$ hertz is

\begin{equation}
\ E_B = \int^{+B}_{-B} T^2 \mathrm{sinc}^2(\pi ft)df\,
\end{equation}

setting

\begin{equation}
\ 2 \pi fT = x
\end{equation}

\begin{equation}
\ df = dx/((2 \pi fT)
\end{equation}
plugging into the integral
\begin{equation}
\  E_B = (T/ \pi)  \int_{0}^{2 \pi BT} \mathrm{sinc}^2(x/2)dx\,
\end{equation}
\begin{equation}
\ E_g = T
\end{equation}

\begin{equation}
\ 0.9 = (E_b/E_g) =  1/ \pi \int^{2 \pi BT}_{0} \mathrm{sinc}^2(x/2)dx\,
\end{equation}

I wanted to work through this example fully, but I'm not sure how to get this numerically solved, solving for the bandwidth B from the upper bound of the integral using MATLAB.


texts: Finding the output of a system where the input is a sum of complex exponentials

texts:
So, I have to find $H\{ x(t)\})$ (which is an LTI system), where $$x(t) = \sum_{k=0}^{\infty} a_ke^{ \ jw_kt}$$ and where the impulse response of the system is given by: $$h(t) = \frac{\delta(t+\tau)-\delta(t)}{\tau}$$

Here's is my attempt. I know that $$H\{x(t)\} = \int_{-\infty}^{+\infty}x(t-u)h(u)\ du \\ =\frac{1}{\tau}\left(\int_{-\infty}^{+\infty}x(t-u)\delta(u+\tau)\ du\ \  - \ \ \int_{-\infty}^{+\infty}x(t-u)\delta(u)\ du\right)$$
but I don't know where to go from there. If you could help, that would be really nice, thank you.


texts: Precise Centre frequency of an All-pole digital filter

texts:
I'm not an engineer and have essentially taught myself all I know - this present problem is giving me some problems.

I have an all-pole filter (a gammatone filter - impulse response is essentially a damped cosine) implemented as a cascade of identical 2nd-order pole pair sections. The poles are located in the z-domain at $\exp(-\beta + i\theta)$: $\beta$ is the bandwidth in radians/sample and $\theta$ is the 'ringing' frequency in radians/sample.

Given a value for $\theta$, this pole pair does not produce a filter with peak gain at the frequency $\theta$ - I'm trying to find out exactly what the offset is (i.e. at what frequency a filter with a given $\theta$ actually has peak magnitude), so that I can produce filters that peak at precisely a given arbitrary frequency.

I've found an equation for calculating what this offset is given poles located in the s-domain at $-\beta + i\theta$:

$$\mathrm{f_{c,actual}} = \theta  \sqrt{1 - \dfrac{\beta^2}{\theta^2}}$$

This formula provides results consistent with observations up to a certain frequency, though I believe it needs to be modified to produce results that work for arbitrary discrete spectrums.. I've tried a number of things and just can't quite figure it out.

Any help is greatly appreciated.


texts: Symmetry of real and imaginary parts in FFT


texts: Can you turn a square wave into a sine wave using a low pass filter?

texts: And if it could, would it make the sound of the square wave thinner than before because of losing its harmonics?


texts: How Are Unmeasured Properties (Velocity and Covariance of Velocity) Handled with a Kalman Filter?

texts:
I'm trying to understand how I can update a Kalman filter with a state variable for position and velocity when I only measure position. I have a covariance matrix of the position measurements. But what does the overall covariance matrix of the measurements look like where velocity is unmeasured?
It seems that zeros would send the gain with respect to velocity to favor the measurement, which would drive the velocity estimates toward whatever value I input for the velocity measurements.
I have a similar problem where I am trying to fuse alternating measurements where I either have X and Y but no Z or X and Z but no Y. I would expect a single filter to fuse the two sets of measurements if I do this right, but I don't have a full grasp on how to handle to covariance terms involving the unmeasured property.


texts: How to shift frequency bins in audio pitch shifter

texts:
I would like to ask you for help to understand shifting frequency bins in audio pitch shifter.
I have two issues that I can’t undersrand:

FIRST QUESTION
At first how to shift frequency bins when pitch factor is not rounded?
For example let’s imagine some not typical but easier to explain case:I have complex numbers array f with freq bins, where each member of array corresponds to one Hertz. So:
f[0] correspond to DC
f[1] correspond to 1Hz,
f[2] corresponds to 2Hz,
f[3] corresponds to 3Hz
and so on.
And I have only one value at freq bin with index 1:
f[1].real = 1;
f[1].imag = 0;

And all the rest members in f array are zero (both real and imag). So it means that my signal in time dmain is just single 1Hz sinusoid with no phase shift and with amplitude equal 1, and can be expressed by: $1 \cdot sin(1\cdot x)$
(Of course apmlitude depends on fft size, but the point is frequency so for simplify let's assume amlitude is 1 in that case).
And now let's assume I want to make pitch shifting with factor 1.75 to achieve one 1.75 Hz sinusoid expressed by: $1 \cdot sin(1.75 \cdot x)$
But I don’t have in my f array the member which correspond to 1.75 Hz. So I thought I should spread it on two bins with appropriate ratio. Let me illustrate it:

So after such pitch shift there are:
f[1].real = 0.25;
f[1].imag = 0;
f[2].real = 0.75;
f[2].imag = 0;

But in that case it looks like I have two sinusoids, first 1Hz and amplitude 0.25, and second 2Hz and amplitude 0.75. And it could be expressed by:
$0.25 \cdot sin(1 \cdot x) + 0.75 \cdot sin(2 \cdot x)$
So I've plot it on Desmos, to compare with desired $1 \cdot sin (1.75 \cdot x)
And definitely they are not even similar.
So it looks like I do something wrong? But then how should it be solved?

SECOND QUESTION
Second question corresponds to the first one. Let's say I want to pitch shift with factor exact 2, so I should achieve the signal with one octave higher than original. While it works for one sinusoid, I get very ugly results for normal audio signals where all freq bins in array are filled with some value. So after shift all bins by factor 2 then I have every two bins empty, let mi illustrate it:

And it seems (and it sounds) like it’s wrong way to shifting freq bins. But then how to do that properly?

For any help great thanks in advance.


texts: What does the 'B' parameter represent in Wikipedia's page on window functions?

texts: I was trying to learn some more window functions and was referencing the Wikipedia article on the subject. I've noticed on a lot of the plots, there appears to be a capitalized and italicized letter B. What is this parameter? I can't find any reference to it on the page. link to wikipedia article And here's a photo that features two plots, each having an associated 'B' parameter with them. imgur What is this?


texts: Several questions about power of periodic signal

texts:
I am currently learning Fourier series (periodic signals) and there are few things I am not sure about.

(1)

$$
P = \frac{1}{T}\int_{-T/2}^{T/2}|s(t)|^2dt = \frac{1}{T}\int_{0}^{T}|s(t)|^2dt
$$

Let's say, that period is $T = 2$ and my signal $s(t) = t$:

$$
P = \frac{1}{T}\int_{-T/2}^{T/2}|s(t)|^2dt = \frac{1}{2}\int_{-1}^{1}t^2dt = \frac{1}{3}
$$

but

$$
P = \frac{1}{2}\int_{0}^{T}|s(t)|^2dt = \frac{1}{2}\int_{0}^{2}t^2dt = \frac{4}{3}
$$

What am I missing here? Are these equations wrong? I've checked several sources and all of them state it this way. What kind of power do I calculate this way? Is it total power of periodic signal?

(2) I also have a problem, I don't understand how to calculate.
I am given:

\begin{align}
x(t) &= t\\T &= \frac{1}{3}
\end{align}
what is the power on interval
$$
\left(\frac{1}{3};\frac{2}{3}\right)
$$
How do I calculate this?

(3) Am I right if I say, that I calculate power of first $N$ harmonics with Parseval's theorem?

$$
P_n = a_0^2 + \frac{1}{2}(a_n^2 + b_n^2)
$$

I am really trying to understand it, but the more I read, the more confused I am. So I really need to point me in the right direction.


texts: List of all available ARM+DSP processors for surveillance

texts:
I have come across TI's C6000/Integra. I want to know is there any other similar or better ARM+DSPs available. Specially for surveillance/vehicle number plate recognition etc. Having camera input at least for VGA.

Link: TI's-Integra


texts: sensing and representation matrices from exhaustive data

texts: If I know y (samples) and x (signal), is there a practical matlab algorithm to calculate A in y = Ax?


texts: Fourier transform possible on non-rectangular part of an image

texts:
Dear Signal Processing readers,

I want to introduce 'noise' into parts of images.

Until now, I worked with rectangluar images for a similar purpose and did the following, using (inverse) Fourier transform (2D) with MATLAB:

fft = fft2(image);    % run fourier transform on image
amplitude = abs(fft); % extract power spectrum>
newPhase = ( (rand(size(amplitude)) * 2) - 1) * pi;          %generate random phase
newImage = real(ifft2(amplitude. *  exp(1i * newPhase)));    %use both to generate noised image

This randomizes the phase all over the image, preserving the power spectrum and should guarantee the same spatial frequency spectra in input and output.

What I want to do now, is applying this effect to only a part of an image, which is not rectangular (its the "inner" parts of a face and will fit in an ellipse).
So I want to "cut out" this ellipse from an image, turn it to "phase noise" (with power spectrum/spatial frequencies constant) as above, and but it back into the "hole" I cut out of the source image before.

My "prime solution" of course would be that FFT/iFFT could operate on ellipses as well, so I could keep using my method shown above. After some amount of google-ing I am quite pessimistic about that option - is there any hope?

I found the following traces to possibly promising approaches, but as I am not very familiar with the math behind FFT, I'd appreciate any hint/explanation (why) which of these ideas go into a good direction:

Somehow "stretch" the ellipse to a rectangle shape, apply my method as is, strech the resulting rectangle back to the input ellipse shape. Here I fear that the rectangle/ellipse-transformation might destroy the matched power spectrum. Do you know good way of transforming without loosing this information? (This led me to this idea)

Two more hints came from this post:

Using a 3D-fourier transform, in which (as far as I got it) the 3rd dimension gives information about which part of the image shall be used. Any idea / experience how to implement this?
Using "Elliptic Discrete Fourier Transform" - but I do not understand what (and how) it does and how it can be implemented.
In this forum, someone suggested to a related question some sort of window tapering (which, as far as I got it, is made to "filter" the influence that only comes into the image because the ellipse of interest is on a black background) that I did not understand as well, plus i am not sure if it could help me because it did not sound like that question would require an ellipse as output.

Thanks in advance for any suggestion!

Cheers. Jan


texts: Change pitch of a wav file without changing the Sample rate?

texts:
I'm creating a Python library I call PAL, (Python Audio Library) and I need some help. What ways are there to change the pitch of a wav file without changing the sampling rate? I have complete access to the data contained in said file as numbers so if there is an equation or method that would work please let me know. To recap

Change the pitch
Don't touch the sampling rate
I can mathematically change the samples freely


texts: Is there an OBJECT BASED Image Segmentation method, with pure non-segmented training images?

texts:
I wonder if there is a method works as follows:
1) The algorithm is trained with some images of a particular object (I train my Image with 1000 banana images)
2) I use the training data to segment a test image (semantic segmentation)

The critical point here is, the training images DO NOT contain labeling information (foreground - background). I expect the algorithm to learn the shape of the object with some common edge orientations or features etc.

I study over Graphical Models for Segmentation and the most similar method is "OBJCUT" proposed by the following paper:

http://eprints.pascal-network.org/archive/00006963/01/kumar10.pdf

But, they train the system with shape exemplars (ground truth images contain segmentations done by humans).

My objects to segment are non-articulated objects.


texts: Structure of function that describes sine signal

texts:
I need to create a vector of sine signal.

So I'm trying to figure out what the structure of the function that describe sine signal ?

For example the function   $\sin 2x$  meets the requirements?

If the answer is no, then what is the reason?


texts: Are Fast Wavelet Transform coefficients constant if the input signal frequency coefficients are constant?

texts:
I'm currently studying the Fast Wavelet Transform. As I currently understand, the Fast Wavelet Transform is implemented as a QMF filter bank where the frequency resolution decreases as the signal is low-pass filtered and sub-sampled. Pictured below:

However, I'm curious as to why this "works".
Let's say I have a true 750Hz signal that is sampled at 2kHz. The amplitude of this 750Hz signal is 2Vpp.
The Fast Wavelet Transform says that I can high-pass filter with a pass region of 500Hz - 1000Hz, decimate by a factor of 2, and it will output coefficients that correspond to the signal frequency. Because I know the signal is constant frequency 750Hz, I would hope for constant coefficients of magnitude 1 at Level 1.
But I'm confused. Because we're subsampling at 1KHz (half of the original sample rate), we're essentially aliasing the 750Hz signal into the 0Hz - 500Hz range. It seems like the coefficient value will vary depending on where we're sampling. However, because we know the true signal is a constant 750Hz signal, it would be desired that the wavelet coefficient is also constant.
How exactly do the coefficients of the Fast Wavelet Transform correspond to the true frequency content of a signal?
Edit: Are there any good resources for understanding the general properties of Wavelet coefficients? e.g. should a wavelet coefficient remain constant if the signal has constant frequency content?


texts:
It may be a really simple question, but I'm not sure about this one:

Given a complex white Gaussian noise process with iid real and imaginary parts and a double sided power spectral density of $N_0/2$. What is the PSD of this random process?

Is it just:
$$S_n(f) = \frac{N_0}2 + j \frac{N_0} 2\qquad ?$$
Or am I wrong?


texts: Least Squares Approximation for FIR Filter Design

texts:
Following this paper , I am trying to make a least-squares algorithm in MATLAB, but for type I (I know about firls()).

N = 41;
M = (N-1)/2;

wp=0.2;
ws = 0.4;
K = 2;
fp = wp/2;
fs = ws/2;
q= [fp+K*(1-fs), (fp*sinc(fp*[1:2*M])-K*fs*sinc(fs*[1:2*M]))];
Q1 = toeplitz(q([0:M]+1));
Q2 = hankel(q([0:M]+1),q([M:2*M]+1));
Q = (Q1 + Q2)/2;
b = fp*sinc(fp*[0:M]');
a = Q\b;
h = [a(M+1:-1:2)/2; a(1); a(2:M)/2]

How to find the design FIR filter if the desired  response :

How to find the group delay in this case response in matlab


texts:
Assume $M=3$ and $L=2$ for the multirate system shown below. Furthermore, assume an arbitrary fourier transform $X(e^{j\Omega})$. This one lies between $-\pi$ and $\pi$. The fourier transform after the filter $H_{LP}(z)$, which is an ideal low-pass filter with a cutoff-frequency of $\Omega=\frac{\pi}{3}$, lies between $-\frac{\pi}{3}$ and $\frac{\pi}{3}$, since we sample up by a factor of three. If I then sample down by two, the resulting fourier transform lies between $-\frac{2\pi}{3}$ and $\frac{2\pi}{3}$. Now the spectra of the input and output will clearly look different, therefore also the signal in time. Is it right to say that with this specific multirate system there will be information loss, although no aliasing effects occur?

Note: In case my problem description is confusing, I can provide an example.


texts: discrete wavelet transform matrix for vectorized image

texts:
Yesterday I asked about how to extract 2D DFT matrix for a vectorized image. Today my question is how can I extract 2D DWT matrix for a vectorized image.
Fourier transform have this property that rows of the image are transformed first, than columns are transformed. Is there a similiar property for discrete wavelet transform such that 1D DWT matrix can be utilized for 2D DWT?


texts:
The coarse acquisition C/A GPS signal has a symbol rate or chip rate of 1.023 Mcps with code period of 1 ms. The code bandwidth is +/-1.023MHz (to 1st sinc nulls). I'd like to record this signal using a USRP.

What's the minimum complex sample rate that I'd need to use so I could acquire this signal in software? 1 sample/chip = 1.023e6? 2 samples/chip = 2.046e6?

Since the signal is coming from a satellite, there will be an unknown carrier frequency offset due to Doppler. From what I've read, for GPS the Doppler is +/- 10 KHz. What's the maximum residual frequency offset allowable when doing the despreading of the CA code? 500 Hz? Or, what's the SNR loss due as a function of frequency missmatch due to residual Doppler?


texts: Phase noise modeling using gnuradio IIR filter

texts:
I need some advice regarding phase noise modeling using gnuradio IIR filter.  I'm trying to synthesize phase noise (multiplicative) with a given noise mask. In one paper, the authors were able to synthesize the noise process for their noise mask (DVB-S2x) by filtering a WSS signal (AWGN with unity variance) with a combination of filters derived by the least square approach.

In my case, I have to use a noise mask specified by the CCSDS standards, as shown below. I'm trying to generate phase noise by means of autoregressive moving-average (ARMA) modeling. I found a couple of parametric modeling functions. I have tried the [b,a] = invfreqz(h,w,n,m) function based on the phase mask. The function generates poles and zeros and expects the impulse response vector (h), frequency vector in radians (w), number of zeros (n) and number of poles (m). I have generated the frequency vector by choosing "appropriate" frequencies from the mask and converting them to radians. The respective amplitudes were converted from dBc/hz to linear scale.

fs = 2e6; %sample rate
offset_rad = 2*pi*[1e1 1e2 1e3 1.4e3 1e4 1e5 1.1e5 1e6]/(fs);
offset_hz = offset_rad *fs /(2*pi);
lf_db = -1*[25 55 85 95 95 95 95 110];
figure
semilogx((offset_hz),lf_db)
hold on
lf_linear = 10.^(lf_db/20);
[b,a] = invfreqz(lf_linear,offset_rad,5,6)
h = freqz(b,a,offset_hz,fs);
w = offset_hz;
%freqz(b,a)
semilogx((w),20*log10(abs(h)),'r')
legend('ideal','generated')
%xlim([0 1e6])
%10*log10(w*fs/(2*pi));
figure
zplane(b,a)

The obtained response is fairly close to the expected one. The problem I'm facing now is IIR filter instability in GNU Radio. The IIR block expects feed-forward and feedback taps. These are the "b" and "a" filter taps from Matlab (direct copy/paste). The IIR filter block also has a boolean parameter called oldstyle, which is set to false. After running the flowgraph, I was able to see phase noise spectrum for less than a second before "it moves up the graph" and disappear. Inserting a number sink revealed that the output of the IIR filter goes to infinity very quickly. Why does this happen while all the poles are inside the unit circle? Any suggestions on how to generate phase noise in a more elegant way are highly welcomed too.


texts: Can anyone help with cross correlation error in Simulink?

texts:
I am sorry if I say any thing silly. I have a vector x below and trying to implement correlation with itself (lets say x2) in Simulink, but getting error. I am getting all zeros in result of my output vector. Kindly can some one please help me out in it. Thanks in advance.

ts=1*10^-9;
t=0*10^-9:ts:150*10^-9;
fc=1*10^6;
yi=154*10^-9:ts:2000*10^-9;
y=yi*0;
xi=2*sin(2*pi*fc*t);
x=[xi y];
x2=x;
z=xcorr(x2,x);


texts: How do I improve the SNR of my signal given a noise profile (additive noise)?

texts:
From an iPhone microphone, I have a background recording of 5s(fs = 44.1 Khz), and immediately after that, I have a voice recording of a person saying 'aaa...' for 10s(fs = 44.1 Khz). How can I use the background recording (How do I characterize the noise from it?) to improve the SNR in the voice recording? The noise can be colored too, also I am modeling it as a additive noise.

Can someone please point me to any work related / in a similar direction (theory, implementation etc.,)


texts: RMS in frequency domain using Plancherel theorem

texts:
I have an acceleration measurement in x, y and z direction, i.e. three vectors each of length 3000.
$x = [x_{1}, x_{2}, x_{3},...,x_{3000}]$
$y = [y_{1}, y_{2}, y_{3},...,y_{3000}]$
$z = [z_{1}, z_{2}, z_{3},...,z_{3000}]$
Currently I'm taking the RMS of x, y and z then I'm calculating the FFT. This produces a vector of length 3000.
$RMS = [\sqrt{x_{1}^2 + y_{1}^2 + z_{1}^2}, \sqrt{x_{2}^2 + y_{2}^2 + z_{2}^2},...,\sqrt{x_{3000}^2 + y_{3000}^2 + z_{3000}^2}]$
I've now read that it is also possible to apply the RMS in frequency domain then using the Plancherel theorem which states that:
$\sum_{n=0}^{N-1}|x_{n}|^2 = \frac{1}{N}\sum_{k=0}^{N-1}|X_{k}|^2$
if $X_{k}$ is the FFT of $x_{n}$.
As far as I understand it this would produce a single output for my three vectors of length 3000. Is that correct? Is there a possiblity to compute the RMS in frequency domain and getting a vector of length 3000?


texts: How to prove that my combiner vector is $\mathbf H \vec f $?

texts:
When my received signal is $y=\mathbf H \vec f x+\vec n$,they we need a combiner vector in the receiver,how do i prove that combiner vector $\vec w$ is $\mathbf H \vec f $?

After the received signal multiply the combiner,the formula will become this

$\vec w^H y=\vec w^H\mathbf H \vec f x+\vec w^H\vec n$

$\mathbf H$ is a $2 \times 7$ matrix,and it is channel

$\vec f$ is a $7 \times 1$ vector,and it is beamforming

$x$ is a $1 \times 1$ value,and it is signal


texts:
I am working on an audio system where I need to implement an equal power crossfade between two signals over a specified number of samples. Initially I implemented a linear fade but this creates a dip in the middle of the crossfade which I don't want. I tried a few other methods that I found on Google but they either make a dip in another part of the curve or even make the sound louder.

If I have float samples ranging from -1.0 to 1.0 and I want to crossfade between them over the course of x samples, keeping the perceived volume the same over the entire crossfade, what would be the correct formula for calculating the curve? Preferably I want to have a value like 0.001 represent the minimum level that the fade in starts at and the fade out goes to. 1.0 should be the maximum. Thanks in advance for any advice!


texts:
I am a DSP newbie and I am really confused after reading Chapter two of

"An Introduction to Digital Signal Processing by J. H. Karl".

On page 32, first paragraph, the book defines a minimum phase couplet as one with zeros OUTSIDE the unit circle of the z-plane. But I have read everywhere else that a minimum phase system has all zeros and poles INSIDE the unit circle.

Is this a contradiction? When do we have a  minimum phase couplet?


texts: Relation between constellation SNR and baseband SNR

texts:
Let's say I have an I/Q modulation, ie. QPSK.

I can simulate an the effects of an AWGN channel by adding to my QPSK symbols the realizations of a complex normal.

$\mathcal{CN}(0,\frac{1}{\mathrm{SNR}}) = \mathcal{N}(0,\frac{1}{2 \mathrm{SNR}}) + j \quad \mathcal{N}(0,\frac{1}{2 \mathrm{SNR}})$

with
$\mathrm{SNR} = 10^\frac{\mathrm{SNR_{db}}}{10}$

I should also be able to generate a baseband waveform with my QPSK symbols, ie using a 15kHz carrier, add gaussian noise to my baseband signal and then make a soft demodulation.

Using the same variance on both channel models, I observe a lower SNR on the soft decoded symbols than if I directly apply the noise realizations over the I/Q constellation.

Why are both different? Shouldn't I expect even better SNR when applying AWGN to the baseband waveform, as the noise isn't band limited?

What's the relation between the noise variance applied to the constellation and the noise variance applied to the baseband signal?

Edit:

Defining the quadrature carriers as:

$\phi_1=\sqrt{\frac{2}{T}}\cos(2 \pi f_c t),\quad 0 \le t \le T$

$\phi_2=\sqrt{\frac{2}{T}}\sin(2 \pi f_c t),\quad 0 \le t \le T$

With $T$ the symbol interval and $f_c$ the  carrier frequency.

The baseband waveform for an I/Q symbol $s_i=[s_{i1} s_{i2}]$ is:

$s(t)=s_{i1} \phi_1 + s_{i2} \phi_2$

The received waveform at the output of the AWGN channel is

$x(t)=s(t)+ w(t)$

With $w(t)$ a continuous-time wide-sense-stationary white noise process.

Assuming perfect coherence, the demodulation performed by a correlator is:

$\hat{s}_{i1}=\int_0^T x(t) \phi_1 dt = s_{i1}+\mathrm{w_{i1}}$

$\hat{s}_{i2}=\int_0^T x(t) \phi_2 dt = s_{i2}+\mathrm{w_{i2}}$

I'm looking for the relation between $w(t)$ and the I/Q noise $[\mathrm{w_{i1}},  \mathrm{w_{i2}}]$.

My waveform simulation will sample $x(t)$.


texts: What is the reason for filtering out the negative frequencies of a signal?

texts:
I am reading this tutorial. Quoting the lines from the topic "Analytic signals and Hilbert transform filters":

  the corresponding analytic signal $ z(t)=x(t) + j {\cal H}_t\{x\}$ has the property that all ``negative frequencies'' of $ x(t)$ have been "filtered out''

I have a few questions about this:

what is the purpose of filtering out negative frequencies?
What is the interpretation of filtering out negative frequency?
What are negative frequencies in physical world? I read that
    negative frequencies doesn't exist and are introduced to simplify
    the mathematics behind theory.
If negative frequencies don't exist, then why do we all plot two
sided spectrum?
If negative frequencies exist, why do we only consider the first
half of an N-point FFT?


texts: How Does the RMS of White Noise Change with Sampling Frequency?

texts:
There is an analog system which includes the continuous-time linear equalizer (CTLE). With some .noise analysis the power-spectral density (PSD) of the noise in that system is provided. So let's not care about the details how PSD is provided, the only thing that matters is that I have access to the noise PSD of the circuit. The noise PSD is not necessarily white and it can have an arbitrary shape, but it dies out after a certain frequency (let's say most energy is contained within 100 GHz).
By having the noise PSD, I integrate it over a large enough frequency range (in this example from 0 to 100 GHz) to get the (almost) total power of the noise, let's name it $\sigma^2$.

Now, in my MATLAB simulation I have a clean sampled signal $s(n \,\,\Delta t)$ with sampling frequency $f_s$, where $f_s = 1 / \Delta t$ (To add more details, $s(t)$ is a continuous-time signal that has a bandwidth of 30GHz and I am oversampling $s(t)$ by taking 32 samples per symbol time. So the sampling frequency in this example would be $ f_s $ = 32 x 30GHz = 960 GHz . Maybe these details don't matter).

I want to add discrete white Gaussian noise to $s(n \,\,\Delta t)$ such that this discrete Gaussian noise has the same power (i.e. variance) as the continuous (non-white) noise which was calculated to be $\sigma^2$. I know I can do something like $N = \sigma$ x randn( size(s) ) in MATLAB to generate noise with power $\sigma^2$, but my main question is that what is the dependency of the discrete noise in this example to the sampling rate $f_s$. If I change $f_s$ how the randn() function or the $\sigma$ should change.

Cause I see in some pages it is mentioned that:

"Given a continuous White Noise (Wide Sense) with variance $\sigma^2$ and you want sample it at rate of $f_s$ you should generate discrete noise samples with variance of $f_s \, \sigma^2$",

but I don't know what is the mathematical reasoning behind it and if the same applies in my case or not.

If you have any good read on this matter, that would be really appreciated.


texts: dB SPL readings differ from external sound level meter?

texts:
I'm trying to write a program that takes values from a MEMS microphone and calculates the dB of the incoming audio data.

The microphone i'm using is this: https://www.mouser.es/datasheet/2/218/ph0645lm4h-b-datasheet-rev-c-1525723.pdf using the I2S protocol, 24 bits signed integers.

The problem i have is that my dB SPL calculations are reading at 80 dB SPL at silence whereas my sound level meter reads around 30/40 dB(A).

Here is how i'm calculating the db SPL of the microphone input:

#define DC_OFFSET 419430

// this int32_t contains a 24 bit signed integer
int32_t sample = GetSample();

// Remove the DC offset.
sample = sample + DC_OFFSET;

sample = sqrt(pow(sample, 2));
double sampleNormalized = ((double)(sample) / 8388607);

// Convert the normalized samples to dB FS
double dbFS = 20.0 * log10(sampleNormalized);
double dbSPL = dbFS + 120.0;

Why do my readings read from 80 dB SPL whereas my sound level meter reads from 30/40 dB(A)?

UPDATE

Now my GetSample() function returns a signed 18 bit values, a few test samples results in these values, before DB offset adjustment:

78528
77376
73472
80384
79872
83648

DC offset defined as (+-100 5%), max value of 18 bits is 131071, therefore DC offset is (131071 * 0.05) = 6553.0.

When i pass it through this algorithm to convert it to dBSPL, the values now come out even higher than before, but the value responds how it should.

double ref = 131071.0;
double dc_offset = ref * 0.05;

// read 18 bit signed value
int32_t sample = GetSample(); // 18 bits wide.

// remove the dc offset.
double adjustedSample = sample - dc_offset;
adjustedSample = fabs(adjustedSample);

// normalize. the sample to the range 0.0 to 1.0.
double sampleNormalized = ((double)(adjustedSample) / ref);

// map the range from 0.0->1.0 to -1.0->1.0.
sampleNormalized -= 0.5;

sampleNormalized = abs(sampleNormalized);

// convert the normalized samples to dB FS
double dbFS = 20.0 * log10(sampleNormalized);
double dbSPL = dbFS + 120.0;

I'm still calculating db SPL values that are 88 dB SPL at quiet background, but the value actually responds correctly, i will apply the A filter as @Max posted and return the results afterwards.


texts:
I feel my question is somewhat related to this, but I'm not sure how to do what I want.

I get transmissions of ~1000 bits at a time, which may be corrupted by noise.
To try and infer data from partially corrupted data, I wish to try and correlate the sequence with expected patterns.

For instance I expect to find a certain subsequence in the stream, followed by a known amount of bits which hold no known structure, and then another set of expected bits. To illustrate:

1010111---n bits unexpected structure - 111000111

In the above, the two specific sequences are expected, and sometimes even the amount of unknown bits in between (n) is apriori known.

The above example is uncorrupted, but if one several bits flip, I wish to scan the stream and look for the best matching pattern.

I already built "masks" of expected patterns and I need to somehow correlate the received data against them, but I'm new to the subject and I'm not sure how to continue from here.

Since I'm new to the subject any: suggestions, references, or proposed schemes will be appreciated.


texts:
An audio sampler (in "music" terms because I'm not experienced in dsp) is a software which loads an audio file and allows you to play it at different pitches by slowing it down and speeding it up.

The audio file should be playable at any new frequency, and the frequency at which I'm playing it can change at any time, even continuously for "pitch slides".

My understanding is that if I am pitching it up, I should use a lowpass filter first so I don't create frequencies above nyquist that will alias. After that I need to use an interpolation function.

What kind of filter I should use? IIR or FIR? What type? Please also keep in mind that for audio the phase changes and the ringing in the time domain can be undesired because they change the way the transients sound.


texts: Invertibility of Room Impulse Response: Reproducing Research Paper

texts:
I have been trying to reproduce this paper¹. Few things which are unclear to me. The paper talks about finding whether a given Room Impulse Response(RIR) is invertible or not based on Nyquist plot.

How can we plot the Nyquist plot from an RIR(using python/other
opensource tools)?
I understand a Z transform as having both poles and zeroes, but how
do we interpret an RIR in Z plane?
How do they do inverse filtering as mentioned in the
paper(convolving with inverse impulse response)?
How to interpret Nyquist plot of RIR in terms of its invertibility?

¹
Stephen T. Neely and J. B. Allen: Invertibility of a room impulse response in: The Journal of the Acoustical Society of America 66, 165 (1979)


texts: The inverse of an orthogonal matrix is its transpose

texts:
In the following statement I don't understand the case for $\ i = j$:

Let $\mathbf A$ be an $\ m \times \ n$ orthogonal matrix where $\ a_i$ is the $\ i^{th}$ column vector. The $\ ij^{th} $ element of $\mathbf A^{T}\mathbf A$ is

$$
\left(\mathbf A^T \mathbf A\right)_{ij}  = a_i^{T} a_j =
\begin{cases}
1 &\text{if} &i = j\\
0 &\text{if} &\text{otherwise}\\
\end{cases}
$$

Therefore, because $\mathbf A^{T} \mathbf A = \mathbf I $, it follows that $\mathbf A^{-1} = \mathbf A^{T}$.


texts:
In an OFDM transmitter, after IFFT & CP block, a DAC block is used.Since the output of IFFT is complex, and the input to DAC should be an integer, how can we integrate the two blocks? I implemented OFDM (upto CP block) in MATLAB but I could not integrate DAC with it.

Regards
Ehsan


texts:
I was wondering if anyone has had any success in implementing OQPSK in gnuradio.

I have a working theory of expanding on the Constellation Modulator to include a custom OQPSK function.

I have already reviewed a few post, but what has been done seems too complicated.

I am wanting to learn so any reading material, etc, would be appreciated.

Feel free to comment!


texts: Impulse response $h[n]$ from input and output signal using cross correlation on Matlab

texts: I have two signals on Matlab, one being input and the other output. I want to know if taking the cross correlation of these signals would give me the impulse response?And if so, will taking the Fourier transform of that resulting signal give me the frequency response?


texts:
I have a huge data for which I want to apply wavelet and plot scalogram.
the below code works for small data. However, I have a huge signal, for which I am not able to plot scalogram.

Can anyone suggest possible way to

wname = 'morl';
scales=1:64;
coefs = cwt(SignalXX1,scales,wname);%'scalCNT');
TAB_Sca2Frq = scal2frq(scales,wname,1/fs);
z=wscalogram('image',coefs,TAB_Sca2Frq,SignalXX1);


texts: When does a Window function's sidelobe level become out of band attenuation?

texts:
Another question I'm struggling to answer by myself.

For a Window function, say Blackman (response pictured), what is the out of band attenuation? If I'm designing a filter with out of band attenuation of 40dB, is it the level of the first sidelobe that I'm interested in? Or is the out of band attenuation dictated by something else (like 3rd sidelobe level for example?)

Thanks!


texts: Image Reconstruction:Phase vs. Magnitude

texts:
Figure 1.(c) shows the Test image reconstructed from MAGNITUDE spectrum only. We can say that the intensity values of LOW frequency pixels are comparatively more than HIGH frequency pixels.

Figure 1.(d) shows the Test image reconstructed from PHASE spectrum only. We can say that intensity values of HIGH frequency (edges,lines) pixels are comparatively more than LOW frequency pixels.

Why this magical contradiction of intensity change (or exchange) is present between Test image reconstructed from MAGNITUDE spectrum only and Test image reconstructed from PHASE spectrum only, which when combined together form the original Test image?

clc;
clear all;
close all;
i1=imread('C:\Users\Admin\Desktop\rough\Capture1.png');
i1=rgb2gray(i1);

f1=fftn(i1);
mag1=abs(f1);
s=log(1+fftshift(f1));
phase1=angle(f1);

r1=ifftshift(ifftn(mag1));
r2=ifftn(exp(1i*phase1));
figure,imshow(i1);
figure,imshow(s,[]);
figure,imshow(uint8(r1));
figure,imshow(r2,[]);
r2=histeq(r2);
r3=histeq(uint8(r2));
figure,imshow(r2);
figure,imshow(r3);


texts: Determining periodic pattern when there is pattern which period is multiple of former

texts:
I have periodic signal, say, it has year period, yearly pattern. I want to check if also has quarterly pattern.
Event if there is no visible quarterly pattern, Fourier Transform of this signal has peaks on all frequencies which are multiples of 1/year frequency. They're gradually decreasing. I understand why they're present but don't now how I can remove this effect.
How do I tell if signal has quarterly pattern?


texts: What is the difference between clustering and quantization?

texts:
The Lloyd-Max quantizer is a scalar quantizer which can be seen as a special case
of a vector quantizer (VQ) designed with the Linde Buzo
Gray (LBG) algorithm.
In k-means clustering, we are given a set of n data points in d-dimensional space  and an integer k and the problem is to determine a set of k points in $R^d$, called centers, so as to minimize the mean squared distance from each data point to its nearest center. A popular heuristic for k-means clustering is Lloyd's (1982) algorithm.
Confusion :

What is the difference between clustering and quantization. I understand the mechanism of clustering technique which is a unsupervised method of grouping data points whereby we map data points into indices of cluster centers which is closest to it. There are different algorithm to perform clustering  - Is clustering another way of doing quantization? Quantization are of 2 kinds -scalar and vector. k-means algorithm is applied for vector quantization. Again k-means algorithm is also applied in vector clustering. So, is vector quantization = vector clustering?
Thank you.


texts: I want to get impulse response of causal system

texts:
I want to get impulse response of this causal system on first photo.

So I got frequency response on second photo from this system.

and now I am trying to transform this frequency response to impulse response on third photo.

but I can't invert it because there is a square of s on numerator.

How can I invert it to impulse response?


texts: If I build a low pass gaussian filter in the frequency domain should I also normalize it with the sum of the filter coefficients?

texts:
I saw two examples of applying a Gaussian low pass filter to a Matlab image: one- when the filter was built in the space domain and was normalized with the sum of the filter coefficients, and another- where the filter was built in the frequency domain and no normalization was made.

Is it consistent that no normalization is needed in the frequency domain and why is it different in this respect from the spacial domain?


texts: Detecting complex signal within complex AWGN using Neyman Pearson statistic

texts:
I need to construct a Neyman Pearson statistic for detecting a complex signal with known amplitude and unknown phase in additive complex Gaussian noise. My hypotheses are:
\begin{align}
\mathcal{H}_0 &: R_i=W_i,& R_i \sim \mathcal{CN} (0,1) \\
\mathcal{H}_1 &: R_i=W_i + A,&R_i \sim \mathcal{CN}(A,1)
\end{align}
alternatively,
\begin{align}
\mathcal{H}_0 &: R_i=W_i,& W_i \sim \mathcal{CN} (0,1) \\
\mathcal{H}_1 &: R_i -A=W_i,&W_i \sim \mathcal{CN}(0,1)
\end{align}
where $A = A_0e^{j\theta_0}$ is the signal, and $i=1,2,...N$. The signal A occurs over a known region of the data, say $m_1$ to $m_2$, where $ 0 \leq m_1 \lt N, m_1 \lt m_2 \leq N$ and  as such the statistical model is to be applied only over that region.
My first steps are marginalizing the nuisance parameters $\theta_0$, and calculating the statistic, assuming $\theta_0 \sim \mathcal{U}(0,2\pi)$:
$$ \frac{f(\mathbf r|H_1)}{f(\mathbf r|H_0)} = \frac{\displaystyle \prod_{i=1}^{N}\int\limits_0^{2\pi} \frac{1}{2\pi} \frac{1}{\pi} \exp{\left(-\lVert r_i - A_0e^{j\theta_0}\rVert^2\right)}\,\mathrm d\theta_0 }{\displaystyle\prod_{i=1}^N \frac{1}{\pi} \exp{\left(-\lVert r_i\rVert^2\right)}}\text,$$
with $\mathbf r=(r_1,r_2,\ldots,r_N)$ being the observation.
Is my construction correct? I can't see how to proceed from here that will allow me to construct a statistic to compare against a threshold. Any suggestions?


texts: I wanted to calculate response of integrator  of sinusoidal input at steady state via these two methods as mention in image  but these two methods give two different answers at steady state, so where I'm going wrong? What concept I'm missing?


texts: Damped spring mass system - parameter estimation

texts:
I've tried to calculate the parameters of a damped spring mass system of the form

$m~ y''(t)+d~y'(t)+c~y(t)=F(t)$

but I have some problems determining the mass m of the system.

The damped spring mass system is given by a MATLAB function

y_out=mass_spring_system(F,Td),

where F is the input force vector and Td the sampling period. The output of the function is already corrupted by noise.

First, I considered the steady state solution for an input force of

$ F(t) = B\cdot \cos(\omega t)$

$\Rightarrow y_p(t) = \frac{B}{\sqrt{(c-m \omega)^2+d^2\omega^2}}\cos(\omega t-\theta)$.

I get the largest response for $\omega=\omega_0=\sqrt{\frac{c}{m}}$, which is the natural resonance frequency.
I estimated this parameter the following using the Power Spectral Density:

NDFT = 2^10;
h = hann(NDFT);
[Pxx,W] = pwelch(F_n,h,NDFT/2,NDFT);
[Pyy,W] = pwelch(y_mean,h,NDFT/2,NDFT);

Ha = sqrt(Pyy./Pxx);
Ha_dB = 20*log10(abs(Ha));

where F_n is white Gaussian noise and y_mean is the ensemble average for 500 realizations. The maximum value should correspond to resonance at the natural frequency $\omega_0$, right?

The second step is the determination of the homogeneous solution:

$y(t) = C\cdot e^{-\gamma t} \cos(\omega_d t + \varphi)$

For the homogenous solution I created an input force vector of

F = [30*ones(1000,1);zeros(1000,1)];

I estimated $\omega_d$ by taking the FFT after averaging the output $y(t)$.
The parameter $\gamma$ I get from determining the envelope of the damped sinusoidal signal using Hilbert transform, taking the log and performing a least square fit using polyfit.

Using the formulas
$\gamma = \frac{d}{2m}$
$\omega_d = \sqrt{\omega_0^2-\gamma^2}$
$\omega_0 = \sqrt{\frac{c}{m}}$

I can calculate the parameters $\frac{d}{m}$ and $\frac{c}{m}$. But how do I get a value for the mass m?


texts:
My question is about the action of a bandpass filter. First, I have a temporal signal. I convert it in a spectral signal via a FFT. Now, I want to apply a bandpass filter on that spectral signal. But my question is : how to do it ?

My temporal signal is a real signal so, I understand how to multiplicate by a coefficient. But my spectral signal is complex so, what do I have to do ? Multiplicate the real part and the imaginary part by the coefficient ?

Thank you for your answers,

Dr_Click


texts: How to take the FFT in a limited range of frequencies?

texts:
When using the fft function in Matlab to take the Fourier transform, one can specify the number of FFT points, nfft, and the command is written fft(x,nfft). In this case, what is the frequency vector?

How can the Fourier transform be calculated with the same nfft, but over a narrower range of frequencies? How would the frequency vector be calculated?


texts:
Lets say I want to double the number of points in an array f. I had the idea to do this:

F=fft(f);N=length(f);
FF=[F(1:N/2) zeros(1,N) F(N/2+1:N)];
f=ifft(FF);

But the result is not exactly correct. Why is this?


texts: On what basis, N of N-point-fft is chosen? What can go wrong if I don't select it to be  power of 2? How is this N related to sampling frequency?


texts: Resampling: how many samples to zero-stuff or down sample?

texts:
I am a little confused about how resampling works with one filter. I get why you choose the lower nyquist cutoff when you combine both interpolation and decimation for resampling by a rational fraction (you don't want to alias and you only get one shot to take care of it when you combine interpolation+decimation).

I will use a practical example: $44.1 \, \textrm{kHz}$ to $48\, \textrm{kHz}$.

$(44100 / 48000) = (441 / 480) = (147 / 160)$

Combining $147$ decimation factor by $160$ interpolation, you can use one anti-alias:
at $44.1/2= 22050\,\textrm{Hz} $ cutoff for a lowpass FIR. However, you still need to zero stuff, right? But how much? 160 zeroes?

Going the other way around $48 \, \textrm{kHz}$ to $44.1\, \textrm{kHz}$, how many samples do you downsample? Certainly not 147 zeroes?

Is $160 - 147 = 13$ zeroes the answer?

In the multi-stage solution, here's how this would work (I think): Stuff $160$ zeroes, anti-alias lowpass at $24\, \textrm{kHz}$. Downsample the output by $147$, antialias lowpass at $22050 \, \textrm{Hz}$. But wait, you're not computing outputs for a lot of those zeroes, so why not just combine the two? Does this leave me with $13$ zeroes for resampling up and $13$ samples to skip over when resampling down before applying the lower of the two cutoffs?

Next thing to wrap my head around: polyphase filterbanks.

EDIT: source of confusion is in last comment to answer, interpolate is upsampling followed by anti-image filter while decimation is anti-alias followed by downsampling. cascading these two obviates the need for two separate filters. if this is the case then would this mean that resampling by a rational factor < 1 entails 2 filters always by the nature of the decimation/interpolation process, while a resample factor > 1 entails that 1 filter may be used instead?

DOUBLE EDIT: thinking about it, you can still use 1 filter, but you must interpolate first, then anti-alias, then decimate.

sorry (definitely not sorry) this is how my brain debugs its malformed thoughts on DSP (:


texts: What Is the Correct Way to Apply Richarson Lucy Deconvolution to Luminance Data?

texts:
My question concerns the Richardson–Lucy deconvolution algorithm, which is described in Richardson's original paper. I am interested in applying it in the context of a raw image converter for digital photography. The RawTherapee discussion that motivated this post can be found here.

By using a Poisson noise model, the original paper assumes the algorithm is being applied to the linear data recovered from the sensor, which in particular is not gamma corrected any does not have any additional nonlinear tone mapping applied. A Poisson random variable is no longer Poisson after a nonlinear transformation, so the model no longer applies to such data.

My understanding is that if we simply demosaic the data from the Bayer array on the sensor, then apply deconvolution to the resulting RGB data, we are in good shape because it is linear.

However, it seems traditional to apply the deconvolution only to the luminance data. (I'm not totally sure why, or whether it is advisable, so please consider this a sub-question.) RawTherapee does this, for example, by converting the image to a LAB color space and deconvolving the L channel.

But there's the issue. It doesn't seem like L is linear in the number of photons captured. I found some data that suggests reflectance percentages of (12.5, 25, 50, 100) correspond to L values of (42.4, 57.0, 75.8, 100). Possibly this conversion depends on the illuminant being used; I am not well-versed enough in color science to say. However, the point is that the mapping is definitely non-linear.

So, is there a practical way to apply the Richardson–Lucy algorithm to a "linearized" version of the luminance data? Is this even desirable? I am not sure what the best practices are here, so any references you could provide would be greatly appreciated. I could not find any mention of this issue in the references I looked at.


texts: What's the Difference Between UQI and SSIM Measures for Image Similarity / Quality?

texts: I'm studying image quality assessment (IQA) problem and I'm currently learning about quality measures UQI and SSIM, I know UQI is a predecessor of SSIM, but they look exactly the same for me...


texts: Realtime image-stitching using a CMOS camera array

texts:
I am very new to DSP so please excuse any brevity...

I have to design a 4x4 CMOS camera array and stitch the images together for further processing... I have tried doing that with opencv but the frame rate count drops drastically... So i want to do the Stitching part using some Image Processor Kit which can take multiple CMOS camera input and stitch them and give me a single image as output and the fps should be 60fps ideally atleast 48 fps...

The problem is that I have a very basic level knowledge of DSP so I would need some suggetions regarding which Image Processor or Kit I should use to make a quick prototype... Also what kind of input does the Processor need from the camera(ttl, usb, etc)...

Thanks in advance


texts: Is it possible to simplify the convolution integral if the functions are non-zero in disjoint areas?

texts:
Possible Duplicate:
  When convolving two functions that are constants in a region and 0 everywhere else, where does the integration start?

I have a function $f(x,y)$ and $h(x,y)$.  $f(x,y)$ has a value of $\frac{1}{3}$ when $x$ is between $\frac{1}{3}$ and $\frac{5}{9}$, and $y$ is between $0$ and $1$.  The function has a value of $0$ everywhere else.  Meanwhile, $h(x,y)$ has a value of 1 when both $x$ and $y$ are between $-\frac{1}{19}$ and $\frac{1}{19}$.

These two functions are non-zero in completely disjoint regions.  Is there a way to leverage this property to simplify the convolution integral?


texts: Capacity loss of fast fading regarding AWGN

texts:
The theoric capacity of fast fading channel is:
$$
C_{\rm fading}=E\left\{\log\left(1+\lvert h \rvert^2\textrm{SNR}\right)\right\}, \quad\text{with $E$ is expectance operator.}$$

At high SNR:
$$\begin{align}
C_{\rm fading}&=E\left\{\log\left(\lvert h \rvert^2\textrm{SNR}\right)\right\}\\
&=\log(\textrm{SNR})&&+E\left\{\log\left(\lvert h \rvert^2\right)\right\}\\
&=C_{\rm AWGN} &&+ E\left\{\log\left(\lvert h \rvert^2\right)\right\}
\end{align}$$

The difference regarding AWGN is $\Delta=E\left\{\log\left(\lvert h \rvert^2\right)\right\}$ which is stated as $-0.83\textrm{ bits/s/Hz}$ for standard Rayleigh fading model, or $2.5\textrm{ dB}$ requirement.

Could anyone explain how they come to this result $-0.83\textrm{ bits}, 2.5\textrm{ dB}$ ?
Standard Rayleigh fading is $C\mathcal N\left(0,1\right)$ then $E\{\lvert h \rvert^2\} = 1$. Using Jensen's inequality $$\begin{align}\Delta&=E\left\{\log\left(\lvert h \rvert^2\right)\right\}\\ &\ge \log(E\{\lvert h \rvert^2\}) \\&= \log(1) \\&= 0 \\&> -0.83\end{align}$$ and this $\Delta$ should not be $-0.83$? I think I made some mistakes but I don't know what it is.

The (stupid) mistake I have committed is that log is concave, so the inquality is in wrong side.
A numerical evaluation of the integral proposed by Marcus confirmed the value 0.83 bits.
Thank Marcus.


texts: Redistributing Color in a RGB Image According to a Gaussian Distribution

texts:
I haven't done this stuff in a while.

If I have an image $I$ I can equalize the histogram of the image using some opencv procedure, it's already defined. Equalizing an histogram means essentially to construct a function such that applied to each pixel $p \in I$ the new value it's like it had been sampled from a uniform distribution, assuming gray scale image for such explanation.

I don't remember if what I want to do has a specific name, but I kind of remember it was possible, but suppose I have a gaussian distribution (instead of a uniform one), say $\mathcal{N}(0,\Sigma)$, and I want to re-distribute the colors according to such distribution.

Is there a way to do this in opencv?


texts: Filtering Passband Signals using Complex Baseband Filtering

texts:
So I'm given a pass band filter with specific transfer function $H_p(f)$, I want to implement this via baseband processing.

I already know how to take the input signal $u(t)$ and process it such that I get the I/Q components $u_I(t), u_Q(t)$.

I also have figured out that the following convolution relation holds for the filter output $y$'s I/Q components $$y_I=(u_I * h_I) - (u_Q*h_Q)$$ and $$y_Q = (u_I*h_Q) + (u_Q*h_I)$$ where $h_Q$ and $h_I$ are the filters I/Q components.

I have run into an issue where I can't figure out a simple way to determine $h_Q$ and $h_I$.

I know for the complex envelope of the filter I have $$h(t) = h_I(t) + j\cdot h_Q(t)$$

I can take the fourier transform of this to obtain $H(f) = H_I(f) + j\cdot H_Q(f)$.

So essentially my issue is with how to obtain $H_I$ and $H_Q$.

My book basically gives the answer but I must be missing something since I can't figure out the reasoning behind it. They state that $$H_I(f) = (H(f) + H^*(-f))/2$$ and $$j\cdot H_Q(f) = (H(f)-H^*(-f))/2$$.

I'm guessing its some property to do with the fourier that I've forgotten but I'm hoping someone can explain the reasoning for it. Thanks!


texts: Impulse response of a linear time-variant system

texts:
Is the impulse response of a LINEAR TIME-VARIANT system unique?
I know that it's unique for the L.T.I case!!!
Advance thanks for the HELP....


texts: Mathematically speaking, is a "signal" a function or the set of outputs from a function

texts:
In engineering it is well known that we abuse the notation of a function $f: A \to B$ with the image of that function $f(t) = \left\{y \in B\mid f(t) = y, t \in A\right\}$

When we say: let $f(t) = \sin(\omega t)$ be a signal, do we mean a function or the set of output from that function i.e. is signal a function or the output?
Similarly, looking at the diagram below, would one say that $\epsilon(t)$ is a set of inputs (i.e. scalars) going into the controller, or a function?


texts: Are there characteristics of photographic film that would make it an analog medium?

texts:
There was an article a few years ago entitled "Analog is Not the Opposite of Digital".

Film is not analog, period. I used to shoot film on a Canon AE-1 from the ’70s. Now I have a digital SLR from Canon, and they’re obviously extremely different. But we have to be careful not to confuse ‘old’ and ‘new’, with two very specific terms like analog and digital.
The word digital, to most people, refers to a device that can capture, store, or display data in a binary fashion. Ones and zeros, on and off, digital is all about numbers. Digital shouldn’t be confused with binary, of course, as digital simply means concrete values. Any system that utilizes solid values (or digits) is digital, binary is simply the most common system. Digital cameras capture light with a sensor, that light is converted into data (numbers), so the use of the word ‘digital’ for your cell phone camera or DSLR is accurate.
Analog, however, is a very abused word. I would venture a guess that the significant amount of readers have used the word ‘analog’ to refer to film cameras. If the new, fancy robot cameras are ‘digital’ then our aging film cameras are ‘analog’, right?  Not at all. Older cameras capture light with film, which is basically plastic, gelatin, and silver halide. When you take a photo (perhaps of your dog drinking a beer), photons hit this material and produce a latent (invisible) image, that can later be brought into view by bathing the film in various chemicals. You could write hundreds of blog posts on film development alone, but the point is that film photography is a chemical process.

I understand the author's point about this misuse of the word "analog" in general, but it seems as though chemical film possesses many of the characteristics of an analog medium, in that it can capture a continuous spectrum of color values over a particular range.
Does this characteristic (or any other) of chemical film make it directly comparable to other analog media such as a record or a cassette, or does the fact that the method of storage of the light information is not a continuous wave like a phonograph record preclude it from being classified as such?
Since I know it will likely show up in a search for these topics, I had found Analog Photography appears to be a subfield of the art of photography, in which "progressively changing [the] recording medium" creates the image, regardless of whether chemical film or digital capture is used.  The technique itself seems to be the analog process in that case, and not any sort of statement about the photographic medium.


texts: Transfering the concept of time delay to image processing

texts:
I've just come across a paper that interrelates the covariance matrix of time discrete signals to their autocorrelation function (or a time-delay, respectively), i.e.
$$\mathbf{C} = E\{\mathbf{x}\left(t\right)\mathbf{x}^H\left(t-\tau\right)\}.$$

Coming from the image processing field, I have usually just worked for the $\tau = 0$ case, whereas $t$ was just the pixel index instead of a time stamp.

Now, I am wondering how I can transfer the concept of this time delay $\tau$ to image processing.

Side-note: Since I don't use box-car based image processing and determine the pixels for this covariance matrix estimation adaptively, what I have basically is a $\left(K\times N\right)$ data matrix containing $N$ pixels with $K$ "measurements" each.

How can I calculate the time-delayed covariance matrix described in the formula above in my case?


texts: Doubt about sampling rate and spectrum with USRP and gnuradio

texts:
I have started working recently with software defined radio (SDR) as well as gnuradio. I have worked with digital signal processing for a while and one feature at one tutorial, which I saw for calculating the spectrum from a bandpass signal, intrigates me.
The implementation (tutorial) is a frequency scanner which uses a sample rate (Fs) for sampling a selected band (B) with a center frequency (f0) from the USRP, then the streamed data is sent to a block which creates a vector with 1024 elements, and the vector is sent to the fft module.
The process described before it seems correct (a bandpass signal has its frequency spectrum calculated). However, what it bothers me is to see that: the tutorial states that the measured spectrum lay between f0-B/2 and f0+B/2, which would lead to the spectrum band equal to B. However, I believe from sampling theorem, this spectrum band is subject to undesired aliasing. Therefore, I believe the correct spectrum should lay between f0-B/4 and f0+B/4.
Am I correct with my belief about the spectrum band, or am I missing something?


texts: How can I Improve my 2D Fourier Transform Graph?

texts:
I'm very new to signal processing. I have been tasked with running Fourier analysis on some velocity waves that I have found in a simulation. My supervisor would like me to plot the angular velocity vs. the wavenumber along the y=x axis using numpy. I have made an attempt at this, although I am not convinced the graph is as good as it could be. Could I get some advice on how to represent the data better? I've included my questions at the end of this post.
Context of the Code:
I'm trying to analyse waves along the diagonal axis of a simulation. Specifically, I'm looking at the contribution to the x component of the velocity. My data size is around 500 points in the time and position domain respectively. The matrix I'm analysing; vx_diagonal matrix is a 2-dimensional matrix storing vx values with respect to time and position along the y=x axis. I.e. vx_diagonal_matrix[0] would be all the values along the x-axis at t = 0.
What I'm expecting, are lines in the angular frequency / wavenumber plot. I'm seeing these, but I believe the graph can be made better. I'm just not really sure what direction I should take to do that.
Code:
#Calculate the Power of the Fourier transform of the Vx Waves
fourierAmplitudes_diagonal_axis = abs(np.fft.fft2(vx_diagonal_matrix))/(len(time)*len(y_coordinate_data))

#Get rid of the extra negative values from the transform
fourierAmplitudes_diagonal_axis = fourierAmplitudes_diagonal_axis[0:len(fourierAmplitudes_diagonal_axis)//2, 0:len(fourierAmplitudes_diagonal_axis[0])//2]

#Generate the Frequency Axis
frequencyList = 2*np.pi*np.fft.fftfreq(len(time), dt)
frequencyList = frequencyList[0:frequencyList.size//2]
wmin = frequencyList[0]
wmax = frequencyList[len(frequencyList)-1]

#Generate the Wavenumber Axis
k_diagonalList = np.fft.fftfreq(len(diagonal_coordinate_data), dr) #dr is the spacing between the diagonal data points
k_diagonalList = k_diagonalList[0:k_diagonalList.size//2]
k_diagonalmin = k_diagonalList[0]
k_diagonalmax = k_diagonalList[len(k_diagonalList)-1]

#Plot the Graph
fig, ax = plt.subplots()
IM = ax.imshow(fourierAmplitudes_diagonal_axis, origin = "lower", extent = (k_diagonalmin,k_diagonalmax,wmin,wmax), aspect = "auto")
fig.colorbar(IM)
plt.xlabel("Wavenumber")
plt.ylabel("Angular Frequency")

#Zoom into the non-zero result region. [The rest of the graph is blank].
zoom_max_k_value = k_diagonalList[int((0.5/6.5)*len(k_diagonalList))]
plt.axis([k_diagonalmin,zoom_max_k_value,wmin,wmax])

Here is what this produces:

My Questions:
From my understanding; the diagonal lines are what I'm looking for. However there are a number of problems that I am concerned about, and I was wondering if anybody with more experience could give me more advice about.
1 - The yellow peak is clearly a concern. From looking around, this could possibly be a DC-Offset? I'm not really sure. Would anyone have any suggestions of what to do with this?
2 - The diagonal lines that are present, are not exactly smooth. This is over a data set of 500 data time points and 500 position data points. I would imagine it would not be as, blocky as this. Does anybody know what may have caused this?
3 - I'm unconvinced I've plotted the wavenumber and frequency axis correctly. Does what I have done seem correct?
4 - I've had to zoom into the region to visualise the data better. Here's what it looks like unzoomed. Is it normal for the rest of the graph to just be zero? That strikes me as odd.

I appreciate the time looking at this post. Let me know if you wish for any more information from my end.


texts: Python scipy.fftpack.rfft frequency bin mapping

texts:
I'm trying to get the correct FFT bin index based on the given frequency. The audio is being sampled at 44.1k Hz and the FFT size is 1024. Given the signal is real (capture from PyAudio, decoded through numpy.fromstring, windowed by scipy.signal.hann), I then perform FFT through scipy.fftpack.rfft, and compute the decibel of the result, in whole, magnitude = 20 * scipy.log10(abs(rfft(audio_sample)))

Based on this, and this, I originally had my mapping from the FFT bin index, k, to any frequency, F, as:

F = k*Fs/N  for k = 0 ... N/2-1 where Fs is the sampling rate, and N is the FFT bin size, in this case, 1024. And the reverse as:

k = F*N/Fs  for F = 0Hz ... Fs/2-Fs/N

However, realizing that the rfft's result is no symmetric like fft, and provides the result, in an N size array. I now have some questions in regarding the mapping and the function. Documentation unfortunately did not provide much information as I'm novice in this area.

My questions:

To me, the result of rfft on an audio sample can be used directly from the first bin to the last bin, as no symmetry occurs in the output, is that correct?
Given the lack of symmetry from the above, the frequency resolution appears to have increased, is this interpretation correct?
Because of using rfft, my mapping function from bin index k to frequency F is now F = k*Fs/(2N)  for k = 0 ... N-1 is this correct?
Conversely, the reverse mapping function from frequency F to bin index k now becomes k = 2*F*N/Fs  for F = 0Hz ... Fs/2-(Fs/2/N), what about the correctness of this?

My general confusion arises from how rfft is related to fft, and how the mapping can be done correctly while using rfft. I believe my mapping is offset by a small amount, and that is crucial in my application. Please point out the mistake or advise on the matter if possible, thank you very much.


texts: How do I quantify the amount of information being lost by filtering a digital signal?

texts:
Let's say I'm doing feature extraction from some sensor reading using a DTW algorithm. I want to extract and catalogue a set of features.

Before I send my signal over to the DTW extractor, I have several tools to remove what I consider to be noise from the signal. I can use bandpass/notch filters or DFTs.

When doing this interactively, it's easy for me as a human to spot which part of the spectrum is interesting and which part of it is low entropy noise.

For instance, the following signal:

has three signals superimposed.

a low frequency square pulse (400 samples long - 2.5mHz)
a mid range sine wave (33 samples long - 30mHz)
and a high frequency noise.

To me, this is trivially visible. I say visible, but this is also audible  if I were to listen to the signal.

The ground noise happens to be the RF interference that a brushless fan makes when it sporadically comes on. So I could have the noise there or not there depending a completely random variable.

It's not clear whether the middle frequency sine wave is simply a superposition on the square wave, or if it is being generated by the same fundamental underlying phenomenon. If it is not the same phenomenon, there's a possibility that their respective phases will drift, and the rising edge of the square signal will only sometimes align with the rising edge of the sine wave. This means that for what ends up being a square wave signal, I could have dozens and dozens of different looking DTW features, and this is something I'd rather avoid.

Having described the above scenario, my questions are as follows:

what metric can I use (or even theory or branch of analysis if the answer isn't a simple one) to determine the amount of filtering necessary? For any given analysis window of say 500 samples, how do I determine programatically how to preprocess my signal before trying extract features?
given a sample window with a particular spectral distribution, does an optimal cutoff for a bandpass filter exist? (or a notch frequency)
if so, does it always exist?
if an optimal value exists, how do I obtain it?
what is the relationship between this and the amount of information in the signal?

Notes: the vagueness of my question is part of my problem. I don't know where to start looking for answers.

To be clear: I'm not attempting to get more information than there is. I'm wanting to degrade my signal as little as possible while removing obvious low entropy noise.

I've looked into things like adaptive filters, but these seem to often time have a separate measurement to aid in the feedback loop.


texts: How do I go from LPC coefficients to a filter polynomial?

texts:
I've recently begun experimenting with LPC, and while I understand that it works, I'm still slightly lost on why it works.

Specifically, I understand that LPC involves finding coefficients $a_1, a_2, \ldots a_p$ such that for a signal $h$,

\begin{align}
h(n) &\approx \sum_{i=1}^{p} a_i h(n-i) \\ &= a_1 h(n-1)+a_2 h(n-2)+\ldots+a_p h(n-p)\end{align}

In other words, we want to find $a_1, a_2, \ldots a_p$ which will let us linearly approximate the next sample given the previous $p$.

However, when LPC is actually used, these coefficients $a_1, a_2, \ldots a_p$ are usually treated as the coefficients of a polynomial:

$$x^p + \sum_{i=1}^{p} a_i x^{p-i} = x^p + a_1 x^{p-1} + a_2 x^{p-2} + \ldots + a_p$$

Matlab's LPC documentation refers to this as "the prediction filter polynomial", and it's very useful for e.g. finding the formants of a signal (which are simply the zeroes of this polynomial).

But I don't understand the connection between these two uses. Why is it useful to turn the linear prediction coefficients into a polynomial like this? And why does it work?


texts: Frequency bin of the positive and negative frequency

texts:
I am using this MATLAB tutorial for Frequency-Domain Linear Regression. There is one part of code they provide, where it is necessary to determine the "frequency bin of the positive and negative frequency", for fft and ifft:

freqbin = 72/12;
freqbins = [freqbin 72-freqbin]+1;
tsfit = zeros(72,1);
tsfit(freqbins) = tsdft(freqbins);
tsfit = ifft(tsfit);
mu = mean(ts);
tsfit = mu+tsfit;

Length of time series is $72$, and $12$ months is considered as one cycle.

How can there be only one frequency bin for positive frequency?
How do we know that it is exactly $72/12 + 1$ ($+1$ is because first bin is for zero frequency)? Is this some formula or what?


texts: Comparing two small audio files (e.g., hand clap) Android

texts:
Dears,
I am new in signal processing and I am dealing with an Android application that should do some audio files (having the same "small" duration) comparison in real time.
My main focus is about comparing two small audio files (about 10 seconds each one) which are related to a hand clap for example (or other sounds related to one short event like ball falling).
In other terms, I should save a reference audio file (or a significant representation you may suggest) and then compare (in real time) a newly recorded audio file with the one already saved.
For instance, if my audio of reference is a hand clap so I decide that the 2 audios are similar if the second one includes a clap even if there are some noises or the 2 sounds intensity is not the same.
I already checked musicg lib for doing this kind of comparison however, I have to convert my recorded audio files to .wav which may add some slowness in my real time application.
I am wondering if there is another alternative that fits better small audio files (about 10 seconds) having the same duration ?
Thanks in advance for your suggestions and help.


texts: Why the carrier frequency does not appear after AM modulation?

texts:
Okay, this might be a very basic question, but I'm completely new to the field.

I have a $40$ Hz sine wave ($1000$ samples, $fs = 1$ kHz), which I have amplitude modulated to a $4$ Hz sine wave ($1000$ samples, $fs = 1$ kHz). In the time-domain, the result looks like I expected. However, when I calculate the FFT, there are $2$ sidebands but nothing at $40$ Hz anymore. If I increase the FFT resolution with zero-padding I can see something there, but very modest.

Is this expected behaviour? Can someone explain what is happening here and why. Thank you!

EDIT: to clarify, I get two peaks which are at 16 Hz and 24 Hz.


texts: Correct comparison of FFT noise with FFT signal

texts:
I have a sensor which I am using to track a strong but noisy varying signal. When the signal is present there is a lot of energy in the FFT bins, when there is no signal the is much less energy - but still not zero.

I am wondering what is the correct way to use the no-signal FFT data (at startup I know there is no signal) to determine whether there is actually a signal at a later time.

Currently I am creating a noise reference by taking an FFT of the background noise and averaging it across all the bins. I am then using that reference via simple comparison to determine whether a bin contains a signal or not.

This works ok, but not perfectly and I am wondering whether calculating the Power Spectral Density of the noise and then comparing that bin-by-bin with the PSD of the signal would give better results. I would plan to scale the PSD as per https://holometer.fnal.gov/GH_FFT.pdf to eliminate the effects of the Hanning window that I am using.

I can't use ensemble averages as time resolution is important, although I could generate an ensemble average of the noise at startup if that would help.

Update: I found this paper: https://pdfs.semanticscholar.org/e59e/303e8f8fca708fe557bc36babb254ffa07f8.pdf which seems almost exactly what I want but too complex. It seems to suggest finding a noise floor by taking the minimum PSD from the spectra and then doing a direct comparison with the PSD of the signal - but clearly its got to be more compicated than that. Could I use this to calculate SNR of the peak bin and configure the minimum acceptable value?


texts: Shannon formula, why square the channel gain?

texts:
In Shannon formula below:
$$R = W\log_2\left(1 +\frac{P_t\cdot H^2}{W N_0}\right)$$
why square the channel gain? When calculating the channel capacity, why do we need to square the channel gain?

Notation: $W$ is the bandwidth, $P_t$ is the transmit power, $H$ is the channel gain and $N_0$ is the thermal noise power density.


texts:
My goal is to build a real-time filter for EEG data. The data are being recorded in 512Hz and I wish to filter it (in a real-time, from 0.5Hz to 30Hz in the most optimal way).

Unfortunately I'm not expert in DSP stuff, so I'm not sure which type of filter I should use, an IIR or FIR? What is the recommended filter type butter, cheb1, or cheb2?

All I understand is there are some properties to set when designing a filter. Among of them are the transition band, the order of the filter, a group delay, Gibbs-effect, linear phase, etc.

Currently, I'm really puzzled about these lists of properties. I also believe that some of the properties mentioned above are dependent on each other. I couldn't find any recipe that is "easy to use " on the web of how to build such a filter and what considerations should be made when building such filters.


texts: is there a way to apply a median filter diagonally?

texts: For example, some studies show how to smooth a spectrogram with a 2D median filter vertically and horizontally (along time or frequency). My question is how would I apply 2D median filter along the diagonals of the spectrogram instead? This would be useful for a certain study I am doing, so I need to understand how that would be coded in MATLAB.


texts: frequency/phase offset estimation in satellite DOA beamforming applications

texts:
In satellite QPSK demodulation, we first estimate and remove the phase offset introduced by the channel. Then, we demodulate the baseband information bits (I & Q channels).

In satellite DOA beamforming, Channel can also introduce phase / frequency offset. When we detect the desired signal from a certain direction in the frequency domain (spectrum), how could we detect this channel introduced phase / frequency offset?


texts: What is the standard way to use z-score with EEG?

texts:
I can't understand how to use normalization technique for EEG.

Is it used for results of FFT or for raw EEG signal?
Are there different methods with common name "$z$-transform"?
And what z-score used in this article?


texts: What are the various techniques for detecting Walls in the Building(Architectural) floor plan Images?

texts: I have to process building floor plan images to extract Walls from the structure. It is trivial in case of binary images; but the images in our case are colored, and each have different color of walls. This restricts us from using any thresholding operation as we cannot assume the Walls will always be the darkest. How possibly can I come up with an algorithm that can work  on all sort of Images?


texts:
I found this formula in a paper : http://ieeexplore.ieee.org/xpl/articleDetails.jsp?reload=true&arnumber=5655742

could someone explain the meaning of "j is frequency sample in the highlighted part.(since frequency limits have defined earlier)