texts: Questions about reproducing results of a paper on Communication Systems and Wireless Communication

texts:
I'm trying to reproduce the results of a paper on Device-to-Device Communication, channel assignment and mode selection (G. Yu, L. Xu, D. Feng, R. Yin, G. Y. Li and Y. Jiang, "Joint Mode Selection and Resource Allocation for Device-to-Device Communications," Nov. 2014.).

I have a couple of questions on the noise spectral density, how to generate the channel taps, and finally the throughput of the entire system.

It is mentioned the following in the paper:

  "All links are assumed to experience independent block fading.
  Hence, the instantaneous channel gain of the interference
  link between CU m and the receiver of D2D pair k can be
  expressed as
  $$h_{k,m} = Gβ_{k,m}d^{−\alpha}_{k, m},$$
  where G is the path loss constant, α is the path loss exponent,
  $β_{k,m}$ denotes the channel fading component, and $d_{k,m}$ is the
  distance between CU m and the receiver of D2D pair k.
  Similarly, we can express the channel gain between CU m and
  the eNB as $h^C_{m,B}$, the channel gain between D2D pair k as $h^D_k$,
  and the channel gain between the transmitter of D2D pair k and
  the eNB as $h^D_{k,B}$."

Towards the last section, we have the following relevant simulation parameters:

Noise Spectral Density: -174 dBm / Hz.
Path Loss Model (for Cellular Links): $128.1 + 37.6 \times \log_{10}(d[km])$.
Log Normal Shadowing: 10 dB.
Uplink Bandwidth: 3 MHz.

Now my questions are related to the channel model:

How do I convert from the Path Loss Model combined with the Shadowing Effect of 10 dB to the first mentioned equation ($h_{k,m} = Gβ_{k,m}d^{−\alpha}_{k, m}$) ?
When I convert the noise spectral density from dBm/Hz to dBm, should
I divide the value by the number of Uplink channels, which is 20 ?
And if so, should I do that in the dB or Watts units ?


texts: Approximate using DFT phase shifting property

texts:
I have a discrete signal $x(n)$ having $N$ samples with DFT $X(n)$. Here $N$ is large say $N=600$. Let the samples of $x(n)$ be, $x(n) = \big[x(0), x(1), x(2), ..., x(N-1)\big]$. But if suppose instead of sampling at the zeroth instant, we began at $k^{th}$ instant, then the new DFT using phase shifting property will be:
$$x_k(n) = \big[x(k), x(k+1), ..., x(N-1), x(0), x(1), ...,x(k-1)\big] \rightarrow  e^{\frac{j2 \pi kn}{N}}X(n)$$
Now consider $k$ to be small say, $k=2$ but instead the obtained signal is:
$$x_2(n)=\big[x(2), x(3), ..., x(N-1), x(0), x(1)\big]$$
My question is about what happens when you change $k$ entries to something else, like this:
$$\tilde{x}_2(n)=\big[x(2), x(3), ..., x(N-1), y(0), y(1)\big]$$
where $y(0) \ne x(0)$ and $y(1) \ne x(1)$.
Thus given that $N$ is large and $k$ is small can we approximate DFT of new sequence w.r.t $X(n)$? Can we say that still the DFT of new signal will approximately the same as $e^{\frac{j2 \pi kn}{N}}X(n)$? If not can we have some error bound? You may further assume that:
$$\text{max}|y(0) - x(0)| = \text{max}|y(1)-x(1)| =c$$
where $c$ is some real valued constant.


texts: Why are the observation features of an HMM-based recognition/synthesis system modeled by a Gaussian distribution?

texts:
Why are the observation features (namely MFCCs) of an HMM-based recognition/synthesis system modeled by a Gaussian distribution?

Even the state duration is modeled by a Gaussian in this paper:

K. Tokuda et al., "Speech synthesis based on hidden Markov models", Proceedings of the IEEE, May 2013


texts: Manipulating RGB color by copying channel value

texts:
I have a picture and trying to swap the RGB channel into RRB, it's like I need to change all G pixels value into R. I did that by copying the R pixel value to G and store it into a new image file (Mat).
Here's my code :

ori.at<Vec3b>(i, j)[0] = image.at<Vec3b>(i, j)[0]; // B-->B
ori.at<Vec3b>(i, j)[1] = image.at<Vec3b>(i, j)[2]; // G-->R
ori.at<Vec3b>(i, j)[2] = image.at<Vec3b>(i, j)[2]; // R-->R

however the result is a little bit confusing me,

Could someone give me a brief explanation how is that the Red and Yellow pixels in the wall end up have the almost same result color (aquamarine) ? I checked the pixel value at some points and the value is almost the same.

What I expected was, if I copied R pixels --> G, part with original R values shouldn't be affected, right ? And because Yellow is a mixture of Red + Green, if I try to copy Red pixels into Green, shouldn't it be Red + Red then ?

[UPDATE]

I have checked sample color values (RGB) for Red wall and Yellow wall.

Red wall :              Yellow wall :
R --> 200 195 185       R --> 237 226 220
G -->   2   1   1       G --> 197 179 168
B -->   1   1   0       B -->   4  29  23

Copying the value R --> G :
Red wall :              Yellow wall :
R --> 200 195 185       R --> 237 226 220
G --> 200 195 185       G --> 237 226 220
B -->   1   1   0       B -->   4  29  23

will result in the Red wall and Yellow wall have almost the same color. But when I checked the result color, it's not aquamarine. Did I make wrong calculation ?


texts: how can we represent a Big Data set in Graph Signals. (Big Data in Graph signal Processing)

texts:
Hi

I have read this article https://pdfs.semanticscholar.org/65d6/1afd9c35b0a75d9de77c2a4a2428af0f7f7b.pdf about Big Data analysis in Graph Signal Processing

I have a couple of questions :

first of all, consider a big dataset which is given to us, it has 150 rows and 365 columns which are sensor measurements and days of one year respectively($D$).

In Example Application section it tries to represent the big dataset in the product of two graphs (the article introduces 3 types of products but I think the dataset have been tested on the Cartesian product).
It breake the data set into two graphs (time series($T$) with 365 nodes and the sensor network($A$) with 150 nodes).

Because there was no tag named Graph signal processing, I could not mention it.

the article has represented this two graph in FIG1 a,b.

Now my question :

if I correctly understand!

if there was no time we just need to decompose the sensor matrix ($A$) and after that, we could take Fourier transform, now during the time we just have more samples from our sensors so why we should eig-decompose our product($T \cdot A$) instead of sensor network($A$).
I mean how increasing the number of sensor data from our sensors, represents itself in graph product.

or why matrix $D$ can be broken to ${T, A}$?(and definitely ${T,A}$ can be processed parallelly)

I would appreciate your solutions.


texts:
I am currently reading up on some speech recognition theory.
As it speech signals usually are non stationary, but if a speech is seen in a very small time frame can it be deemed stationary, as would only consist of the sound of phoneme.

As explained: http://recognize-speech.com/preprocessing#[object HTMLHeadingElement]

But then this comes up:

  As it can be seen in the figure above there can be an overlapping of
  the windows. The frame length and window length is dependent on the
  scope of application and algorithms used. In speech processing the
  value for the frame length  typically varies between 5 to 25ms and for
  the window length  between 20 and 25ms [ref.]. Smaller overlapping
  means larger time shift in the signal, therefore lower processor
  demand, but the difference of parameter values (e.g. feature vectors)
  of neighbouring frames can be higher. Whereas larger overlapping can
  result in a smoother change of the parameter values of the frames,
  although higher processing power is needed.

My question is how the overlap lower/increase the processor demand?


texts: Fast ICA : recover the sources with their mean after prewhitening

texts:
I want to apply the FastICA algorithm on a certain dataset. I believe the sources the Fast ICA can recover are able to explain other data as well (in similar way to factor analysis). However, I also believe the mean of the data contain valuable information. I was wondering if, after the prewhitening which removes the mean from the data, it is possible tk somehow "bring back" the means to the sources extracted, such that the mean of the sources is meaningful as well.
Any reference is well accepted. Thank you in advance


texts: Let $|A|$ denote the power of the analog signal $A$. Is that true to assume $|A+B|=|A|+|B|$, for two arbitrary signals $|A|$ and $|B|$? That is the power of the summation of two signals is -at least roughly- equal to the summation of the powers of the two signals. What if $A$ and $B$ are random signals?


texts: How replicas are formed in Frequency domain when a signal is sampled in Time Domain?

texts:
I know that sampling in one domian (time or frequency) gives raise to replicas in another domain (frequency / time).
How replicas are formed?
What is this Time domain periodicity and frequency domain periodicity here in sampling?


texts: Is there any connection between frequency response of a system and spectrum of a signal?

texts:
For instance, an AR process can be generated through a white noise $w[n]$ passing through an all-poles system (IIR),
$$y[n] = \sum_{k = 1}^{p}a_ky[n-k] + w[n]$$
and we can analyze the IIR system, instead of obtaining the spectrum of the signal, to get some properties of the signal and therefore apply it to some applications, e.g., prediction.

My question is whether all the signals, stationary or non-stationary, can be treated as it is generated via some system driven by white noise, and we can analyze the system, which hopefully can completely characterize the signals.


texts: How can one find symbol duration if bit rate is 100 kbs (as BPSK is used, symbol rate = bit rate)?


texts: Detecting characters in high contrast/oversaturated regions

texts: I'm currently tagging insects with reflective tags and filming them in night-mode on my camera which has an infrared light. I have some markings on the tags which I need to catch on film. Unfortunately, when I film about 2 metres back, the the contrast between the reflective tag and the marking is too great, resulting in the white reflection being the only thing visible. Does anyone have any ideas for how to fix this problem? For instance by using filters?


texts: The quadratic performance surface of an adaptive filter is a paraboloid. Its minimum can be found wherever the gradient is zero. However, since there are two types of paraboloids (elliptical and hyperbolic), is there a way to tell if the minimum detected is a global minimum or just a saddle point?


texts: how can i do a 2^20 FFT if the hardware imposes a limit of 2^10 FFT size?

texts: I need to implement an fft of size 2^20, but the hardware has imposed a limit of 2^10 fft size.Then how can i implement the FFT. please point me to some literature where i can go through the process


texts: Why is the size of the output of the FFT of a signal is same as the size of the input signal?

texts:
I am kind of new to the DSP domain. I was trying to get the frequencies associated with a signal by performing FFT over it. I used numpy.fft.fft for this.

Input signal frequency = 1000 Hz

Noise signal frequency = 200 Hz

The input signal, noise signal and combined_signal all have 2400 samples. When I apply FFT on the combined_signal, i get 2400 size array as output. According to my understanding the FFT operation returns the frequencies present in the signal (in this case, 1000 Hz and 200 Hz). Also, when I plot the output of the FFT, I get some weird plot. Instead of giving two peaks at 1000 Hz and 200 Hz it is giving some random frequencies. Also, why is it giving are two peeks at the beginning and 2 at the end ?

I have used the following jupyter notebook code.

import numpy as np
from matplotlib import pyplot as plt
import matplotlib
%matplotlib inline
matplotlib.rcParams['figure.figsize'] = (16, 9)

freq_noise = 200
freq_signal = 1000

sampling_rate = 48000.0
seconds = 0.05
num_samples = int(sampling_rate * seconds)

signal = [np.sin(2 * np.pi * freq_signal * x1/sampling_rate) for x1 in range(num_samples)]
noise = [np.sin(2 * np.pi * freq_noise * x1/sampling_rate) for x1 in range(num_samples)]

# convert to numpy array
signal = np.array(signal) * 1000      # Increase the amplitude of the signal
noise = np.array(noise) * 500         # Lower amplitude for the noise

combined_signal = signal + noise

print "The length of the signal is : ", len(signal)
print "The length of the noise is : ", len(noise)
print "The length of the combined signal  is : ", len(combined_signal)

plt.subplot(4, 1, 1)
plt.plot(signal)
plt.xlabel('Sample number')
plt.ylabel('Amplitude')
plt.title('Original Signal')

plt.subplot(4, 1, 2)
plt.plot(noise)
plt.xlabel('Sample number')
plt.ylabel('Amplitude')
plt.title('Noise Signal')

plt.subplot(4, 1, 3)
plt.plot(combined_signal)
plt.xlabel('Sample number')
plt.ylabel('Amplitude')
plt.title('Combined Signal')

combined_fft = np.fft.fft(combined_signal)
combined_freq = np.abs(combined_fft)

plt.subplot(4, 1, 4)
plt.plot(combined_freq)
plt.xlabel('Frequency')
plt.ylabel('Amplitude')
plt.title('FFT of Combined signal')

plt.tight_layout()
plt.show()

The output of the above code is:

The length of the signal is :  2400
The length of the noise is :  2400
The length of the combined signal  is :  2400


texts:
I am trying to measure (estimate) THD of an audio path using a consumer sound card. The idea is that I do two measurements, one without the system being measured in the path and one with it:

MEASUREMENT 1:

Sound card DAC ---> Sound card ADC

MEASUREMENT 2:

Sound card DAC ---> Measured system ---> Sound card ADC

There is the same sine wave generated in both measurements. In each measurement, the ADC data is transformed using FFT (with Hann windowing function) and the THD is calculated by the standard THD(F/R) equation (https://en.wikipedia.org/wiki/Total_harmonic_distortion) from the amplitudes from the peaks of harmonics on the FFT.

The idea is that by doing two measurements, one with and one without the measured system, is that I could eliminate the TDH introduced by the non-ideal DAC, ADC and sine wave generated and (only) measure TDH of the measured system.

Since the THD is not a linear phenomenon, I assume this is not a correct way to estimate THD of the measured system:

THD_measured_system = THD_measurement_2 - THD_measurement_1

So my question is what would be the correct approach in eliminating the THD of the measurement system and if this is even a possible approach?

EDIT:

Would checking the difference between amplitudes of the harmonics between measurements 1 and 2 and calculating the THD from this difference work?


texts: How to apply the gain factor G to SOS coefficients when designing a IIR filter on sptool (MATLAB)?

texts:
After designing an IIR filter on sptool, I export the coefficients (SOS) and gain (G) into a .mat file. When filtering, we can use sosfilt(SOS, x), which results in the correct output but with the wrong scaling. I guess this is what the gain variable G is for, but the documentation of sosfilt() never mentions this.
The impulse response of the filter has values of several powers of ten (~10^16), whereas in sptool it had very small coefficients (< 1).
How is one supposed to use the SOS and G variables when exported from sptool, to correctly scale the weights of the exported IIR filter?


texts:
is negative sine, or zero, which is trouble; the imaginary DFT basis (sine) is likewise zero. Is there a way to meaningfully define the derivative of cosine when sampled at $f_s/2$ (such that it's not all-zero)? That is, what'd such $(f'(t))[n]$ look like?:
$$
f(t)=\cos(2\pi 64t);\ n=[0, 1,..., 127] \\
\left(\frac{d}{dt} f(t)\right) [n]=\ ?
$$

The context is 'trigonometric differentiation', or derivative via frequency-domain interpolation as:
$$
f'[n] = \text{iDFT}(\text{DFT}(f) \cdot 1j \cdot \xi), \\ \xi =[0, 1, ..., 64, -63, -62, ..., -1]
$$
(if $N=129$, do $64, -64$ instead). Whenever the sequence is of even length, we end up with a complex $f'[n]$, even if $f[n]$ is real, as the real DFT coefficient at Nyquist turns imaginary. The method is thus unsound for even-length sequences; we also cannot seem to cleanly upsample then subsample via zero-padding, since to subset the original's spectrum we must pad by integer-multiple length of original signal, guaranteeing even length.
Let $x[n] = \cos(2\pi 63n) + \cos(2\pi 64n), N=128$ be the test signal for workaround attempts.


texts: Good book for MATLAB image processing functions explanation

texts: What are some books that can explain what is actually going on in the MATLAB functions for image processing? Some examples would be the function im2bw, imtophat, etc.


texts:
Conceptually when we want to represent a peroidic series, for example a pulse train, we find the Fourier coefficients and get a representation in the time domain.

However what is conceptually wrong with using say an infinite sum of time shifted rect functions to represent it?

Sorry I had a formatting problem, so i am putting this in the bottom post...

My method is as such:

Assuming we have a periodic pulse $x(t)$ such that $x(t) = 1$ for $0 < t< T$ and $0$ for $T < t < T_p$; thus $x(t)$ has a period of $T_p$.

Finding the fourier coefficients Ck via:

$$Ck = \frac{1}{T_p} \int_{-\infty}^{\infty} x(t) * e^\tfrac{-j2\pi kt}{T_p}$$

over 1 period, and thus we can represent x(t) as:

$$x(t) = \sum_k C_k e^\tfrac{j2\pi kt}{T_p}$$

and doing the Fourier transform we will get this form:

$$X(f) = \sum_k C_k \delta(f-kf_p)$$

which is discrete.

However if we consider $x(t)$ to be of this form:

$$x(t) = \sum_n \text{rect} \left(\frac{t-nT_p}{T} \right) $$

Applying fourier transform of $x(t)$ to get (in the form):

$$X(f) = \sum \text{sinc} * e^{-j*W}$$ -- (2)

where sinc() is due to the FT of rect and the e^(-j*W) comes out due to the time shifting property of FT.

Comparing X(f) in (1) and (2), we see that 1 is discrete and the other continuous.

However they come from the same x(t), so isn't this a contradiction?

Sorry for the long post.


texts: Looking for the Concept About All In One Curve Fitting

texts:
I know that there are some Technic for finding curve fitting like polyfit() or pinv(). so we can get a some polynomial equation. it's OK good.

I can get successively a polynomial equation for each curve by using polyfit() or pinv()..it easy. no problem.

But my problems is here.

Actually I want to know a concept to make one equation which from all polynomial equations.

for example,
if I get 3 equations(curve #1's equation ,curve #2's equation ,curve #3's equation) then I want to make all in one equation.\

Is this possible? what kinds of Technic is exist?

  UPDATE:

  I've got 4 polynomials

  1) 3.1734x^3−5.8952x^2+3.8194x−0.0377

  2) 2.1489x^3−3.9026x^2+2.8169x−0.0160

  3) 1.0796x^3−1.9013x^2+1.8500x+0.0008

  4) x

  But I'm not sure how to apply your way about extracting the corresponding parameters.  Would you please help me a small bit more ?


texts:
Im thinking whether any convolutional operation can output a unit impulse, an example to further explain:
where a convolution between system $h[n]$ and unknown system $g[n]$ would output $\delta[n]$.

conv(h[n], g[n]) = δ[n]

where

$$h[n] = \delta[n] + \delta[n-1]/2 + \delta[n-2]/4 + \delta[n-3]/8 + \delta[n-4]/16$$

In regards to Z-transform - it can be concluded that

$$H[z] = 1 + 1/2z^{-1} + 1/4z^{-2} + 1/8z^{-3} + 1/16z^{-4}$$

Moreover, it is know that in z domain conv becomes multiplication therefore:

$$H[z] G[z] = 1$$

And in conclusion $G[z] = 1/H[z]$

$$G[z] = \frac{1}{1 + 1/2z^{-1} + 1/4z^{-2} + 1/8z^{-3} + 1/16z^{-4}}$$

If the above is somehow correct, then when trying to achieve the inverse z-transform we will get:

$$g[n] = \delta[n] + 2/\delta[n-1] + 4/\delta[n-2] + 8/\delta[n-3] + 16/\delta[n-4]$$

We get $1/\delta[n]$, which does not bode well for anybody.

Im I in the wrong? Can it be possible?


texts: Is is possible to use a digital filter (eg. Butterworth) for smoothing sensor data that is sampled at a non-constant rate? The time between two consecutive samples is between 12 to 40 ms.


texts:
For apply least-squares linear-phase FIR filter design,with frequency domain specification is not symmetrical.

The pass-band error function,

$$E(\mathbf{h})_p=\int_{\omega_{p_1}}^{\omega_{p_2}}| \mathbf{c}^T(\omega)\cdot \mathbf{h}-D(\omega)|^2d\omega \tag{1}$$

The stop-band error function,

$$E(\mathbf{h})_s=\int_{\omega_{s_1}}^{\omega_{s_2}}|\mathbf{c}^T(\omega)\cdot \mathbf{h}|^2d\omega\tag{2}$$

h: unknown complex coefficients.

Total error function,

$$E(\mathbf{h})_t=w_1E(\mathbf{h})_p+w_2E(\mathbf{h})_s\tag{3}$$

$w_1,w_2$: weighting constants .

The equation $(1)$ : which represents a filter approximation by determine the error between actual and desired response in pass-band , why the method is different in stop-band (equation $(2)$)  ?

-In matlab ,does it the freqz function achieve the transpose operation ?


texts: Filtering noise from FFT, where noise is known to be shifted delta function

texts:
I have an input sequence $x=\{x_1, x_2 , ... x_n\}$ of reals, where $n=2^m$ for some $m$. I wish to calculate FFT of $x$.

$X=FFT(x)$

However, before  I calculate the FFT, the signal $x$ gets corrupted with noise $\eta$, so $\hat{x}=x+\eta$, and calculated FFT, $\hat{X}$ is FFT of $\hat{x}$, rather than $x$.

$\hat{X}=FFT(\hat{x})$

I wish to  recover actual or approximate FFT $X$ from $\hat{X}$ under following assumption:

$x$ is band limited signal; FFT spectrum dies down quickly.
The noise $\eta$ is a shifted delta function, with unknown shift and magnitude, i.e. $\eta(n)= R\delta (n-l)$, where $R,\ l$ are real and integer respectively.

How will I do it in a compute efficient way?


texts:
I have a sample test with an answer but don't understand how they got to the answer:

  $x(t)$ has info only between $2 < |\omega| < 4$
  $X^F(\omega) = 0$
  for other frequencies.

  All noise was cleaned from the signal and its CTFT is defined:
  $X^F(\omega) = \begin{cases} |\omega|-2  & \text{for 2 < |$\omega$| < 4}\\ 0 & \text{otherwise} \end{cases}$

  The clean signal was sampled
  at $\frac 1 T = \frac {12}\pi$
  Draw the DTFT of $X^f(\theta) \text{ for} -3\pi < \theta < 3\pi$

  Use a clear coordinate system.

The answer has x-axis span of -9 to 9 and y-axis of 0 to 4.

How was this answer calculated?

I thought I would get a simple repeating copy of the original CTFT, but the original CTFT with $2 < |w| < 4$ is a line beginning at w=2 (on the x-axis) with the value zero, going up until the point  $w=4$ where it has the value of 2...

Please help!!  Thanks


texts: Given input and output time series, obtaining filter design

texts:
Say I have 5 minutes of input and output audio. One method I know is to do FFT on windows from input and output. Then divide FFT output into bins and find average energy in various bands for all the windows. After this, take ratio of energy bands from output windows with corresponding energy bands from input windows. And then take average of all these ratios. This will give an approximation to FFT(Y)/FFT(X) i.e. transfer function of the filter.

Are there any other methods to find transfer function given input and output time series?


texts:
I keep reading that a constant group delay over a frequency band is desirable for a filter. But a constant group delay means a linear decrease in the phase delay, which would distort the signal.

Let's assume that I have a group delay of 5 [time units] over a frequency band from 100 Hz to 110 Hz.
That means that the phase delay is $-5\omega + c$.

So:

the 100 Hz component will be delayed by 500 + c rad
the 101 Hz component will be delayed by 505 + c rad
the 102 Hz component will be delayed by 510 + c rad
etc

Considering that different frequencies components have different phase delays, the signal will be distorted.

So why is a flat group delay response desirable for a filter?


texts: Unclear time-to-frequency integration step

texts:
From here; $\hat f=\mathcal{F}(f)$, bar = complex conjugate:

Time-shift property: $x(t-b) \Leftrightarrow e^{-j\omega b}{\bf X} (\omega)$, so why is it $+$ (red)?
What at all is happening? Looks like convolution theorem, taking $f$ and $\psi$'s time-domain multiplication into frequency domain, except don't we need a second integral for conv.?


texts: Relationship between amplitude and sampling rate

texts:
I'm loading a signal with librosa in python. With the original sample rate of 22050 Hz, i get the following waveplot:

When i choose to resample my signal with the sample rate of 512 Hz, i get the following waveplot:

My question is, how come the amplitude is less with 512 Hz than with 22050 Hz? I assume that since the interval at which we record the amplitude is larger at 512Hz, the values should be larger too.


texts: Cascading a filter with itself to form a higher-order filter?

texts:
I've been given an audio signal with some high frequency noise. The professor gave me the system function of a lowpass filter that gives me control over the cutoff frequency:

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

where $\alpha=\frac{1-\sin \omega_c}{\cos \omega_c}$ and $\omega_c$ is the cutoff frequency in radians/sec.

The professor claims that if the filter is cascaded with itself $L$ times, it becomes a filter of order $L$.

My questions are the following:

What does increasing the order of the filter do?
Why would I want to increase the order of the filter?
If I cascade the filter with itself $L$ times, won't I lose the ability to control the cutoff frequency? How could I correct for this?

edit:

Here's what I get for the frequency response of the filter when the cutoff frequency is set to 1000 Hz. Sampling frequency is 16000 Hz.


texts: Apply band pass filter (BPF) before calculating the power of a signal (PSD)?

texts:
From my understanding PSD tries to give a good estimation on the power each frequency attributes to the overall signal power.

If I am only interested in a frequency range $(\omega_0 \leq \omega \leq \omega_1)$, could I put the signal through a band pass filter for that same frequency range then calculate the signal power using this formula:

$$P_x = \frac{1}{N_1 - N_0 + 1}\sum_{n = N_0}^{n=N_1} \left|x(n)\right|^2  $$

Would this technique give me what I am looking for, (the power of the signal from $N_0$ to $N_1$ for the specific frequency range)?

How does this compare to PSD estimation using any of the popular techniques?

Pros, Cons?

Thank you!


texts: filter and resample or resample and smooth?

texts:
Currently doing some signal analysis in python for a major project in my physics degree which is due really soon. I need some help!

Say I have two signals, f(t) and g(t) which are recorded over the same amount of time.

f(t) and g(t) both feature different amounts of noise.

f(t) and g(t) are recorded at different sampling rates so I have two different quantities of measurements.

For analysis I imagine I will need:

a) Interpolation (or resampling?) both signals to get the same sampling rate.

b) Smoothing/filtering to reduce noise.

a then b or b then a?

I would like to ideally have two signals with equal sampling frequencies. I would also like all values of f(t) and g(t) to occupy the same points on the t (x) axis. In addition to this I would like to smooth the data as much as possible with a minimal loss of information.


texts: How to create matched "minimum phase" for a system of parallel FIRs?

texts:
Problem statement

I have a collection of magnitude (only) responses I'd like to turn into FIR filter kernels that are

matched in phase
"minimal" in phase, with respect to the complete system
have little to no "pre-ring"
maintain the design magnitude to a high degree

In the ideal case, I'm hoping for a solution that is analytic rather than parametric.

Meeting the first requirement can be expressed as:

$$
\phi_s(f) = arg\{ H_0(f) \} = arg\{ H_1(f) \} = ... = arg\{ H_n(f) \}
$$

And, is very easy to meet in a linear phase case.

Accumulating phase

In attempting to meet the second condition, I've tried something akin to a "cascading" of phase given minimum phase versions of each filter. I've done this by

returning the "Hilbert" minimum phase for each magnitude response
accumulating these responses, and then defining this as the system phase
reassigning system phase to each response

Here's an expression for the accumulation I used to find the system phase:

$$
\phi_s(f) = \sum_{i=0}^n - \mathscr H \{ \ln( |H_n(f)| ) \}
$$

Unfortunately, for my responses, this results in a poor outcome, returning kernels with way too much smearing / lag at low frequencies in the time domain. All this smearing is far from "minimal" in any sense!

I'm suspicious this may be a result of the two definitions of minimum phase discussed earlier.

Minimum phase of L-2 norm of magnitudes

The most promising results I've had after a bit of hacking is to take the L-2 norm of my magnitude set and use this as a minimum phase target:

$$
|H_{L2}(f)| = \biggl(\sum_{i=0}^n |H_n(f)|^2 \biggr)^{1/2}
$$

$$
\phi_{L2}(f) = - \mathscr H \{ \ln (|H_{L2}(f)|) \}
$$

Setting $\phi_s(f) = \phi_{L2}(f)$ returns kernels with coefficients packed towards the beginning and end of the table. Rotating the table, a linear phase offset, by $\max\bigl(\phi_{L2}(f)\bigr)$, turns the end into a "pre-ring".

While not beautiful, this does get me closer to what I think I'm looking for.

Some more details on my magnitudes

$|H_0(f)|$ is allpass; the remaining are highpass with every increasing roll-off slopes and cut-off frequencies.

Questions

Is there a nice off the shelf know solution? In my mind, the problem is very similar to a matched phase cross-over.

Am I asking too much? In the ideal case I'd like to have no "pre-ring". Maybe that's not possible?

Am I asking the right question(s)? In other words, is "minimum phase" what I'm looking for. What I really want is for my kernels to be matched in phase, and respond "as soon as possible".


texts:
I would like to use Matlab's lsqnonlin functinon with regularization (e.g. Tikhonov regularization). However, to my big surprise there is no such option it seems!
I like lsqnonlin because I can specify a cost function with an  anonymous function in a convenient way.
Is a way to add regularization or another function that does the same thing with regularization?
Also please note that my function takes as input a matrix and convolves with a small kernel k:
fun = @(x)conv2(x,k,'same')-y;
lsqnonlinoptions=optimset('Algorithm','Levenberg-Marquardt','MaxIter',maxIterArg);
x_lsq =lsqnonlin(fun,x0,[],[],lsqnonlinoptions);

It seems that lsqnonlin doesn't expect the user to explicitely write the least squares cost i.e. only what I wrote and it takes care of the squared L2 norm. that is why I don't think it makes sense to just add a +lambda*norm(x) to my func (the latter being the regularization term).


texts: How does sampling a frequency result in a duplicated spectrum?

texts:
This is probably a simple question, but I have a hard time understanding how c. converts to d. in the following picture.

Further reading as an answer is acceptable. After all, I am only asking this because the PDF I’m reading doesn't really explain why c. results in d. The more I look at it, the more I am confused.


texts: In my simulation of an equation in Simulink, I need to have variable factors the same as shown in following picture but Simulink doesn't accept that. how I can solve this problem?


texts: How to generate amplitude and Phase values (I,Q) for 8 QAM?

texts:
I am trying to simulate 8QAM with (4,4) amplitude levels. I want to get a constellation like this:

What are the amplitude levels should I select for I and Q?


texts: Non-causal FIR filter in the feedback loop

texts: I have a feedback loop with a transfer function $H(z) = \sum_{i=0}^{L-1} h(i) z^{-i}$. Is there a way to make this FIR filter non-causal? If it was a feed-forward filter, we could simply do so by buffering data. But can it be done in a feedback loop?


texts:
I am building a scanning-type system.  I have samples in $\textrm{dB}$ and I want to determine when to allow them through.  I am trying to decide between allowing them through when they reach a percentage over their normal value (or of the noise floor), or if I should use a static $\textrm{dB}$ value.  An example of the first way would be to let a sample through when it is $50\%$ over the normal value.  An example of the second would be to let it through when it is one $\textrm{dB}$ over the floor.

Which way is standard, and what are standard values I might use?  I understand it depends on the application, but I'm just trying to get a rough idea.


texts: Practicality of the assumption of i.i.d. for Rayleigh channels

texts:
I would like to understand how far the assumption of having i.i.d. is accurate/valid (from the practical point of view) when an OFDM system is working in Rayleigh channels. Is this means the channel has to encounter flat and slow fading? If not, under what conditions, the assumption of i.i.d. can be considered practically acceptable?

Any hints?


texts: How to invert FFT these magnitude equations to get FIR filters?

texts:
I need filter designs that will reproduce or approximate each of the following:

1) ${\rm magnitude} = e^{-cf^{2}}$

2) ${\rm magnitude} = e^{-\frac{a+b\sqrt{gf}}{df}}$

Where: $a$, $b$, $c$, $d$ and $g$ are constants $>0$, and $f$ is frequency in Hz.

Or if it's easier or an equivalent problem, the multiple of the two filters as one combined filter would even be better:

3) ${\rm magnitude} = e^{-cf^{2} -\frac{a+b\sqrt{gf}}{df}}$

I have never done an inverse FFT or designed a FIR filter though so I'm looking for any basic guidance on how this might be done.

This is for physical modeling so the closer I can get the better. Accuracy is only important with audio frequencies (20-20,000 Hz).

I really need something that will work one way or another. Any guidance is appreciated.


texts:
Given a signal that's been FM encoded on a carrier wave, what are the individual steps I need to take to demodulate it, and recover the original signal?

I don't have a background in signal processing, I have a background in programming. I'm trying to solve a real-world problem: I'm trying to demodulate a raw analog RF signal captured from a Laserdisc at 40 MHz. This RF signal contains a composite video signal centered at 8.1 MHz, two analog audio channels at 2.81 MHz and 2.3 MHz, and one digital data channel encoded at 1.75 MHz and below.For a handy reference image, here's a visual representation:

More specifics on the frequency ranges and other related info can be found in a concise summary here, and in the official documentation describing the final encoding format here.

Now, I've captured a signal like this straight from a Laserdisc. Under a spectral view, it looks like the following:

Here I can clearly see the composite video signal encoded, in this case showing the end of a field, the vblank period, and the start of a new field. I can see the two analog audio channels encoded below, and at the bottom, a band of low-frequency noise which contains the digital channel data. I can also see the harmonic sidebands coming through, giving lower intensity "echoes" of the data at lower and higher frequencies than the encoded frequency.

When the composite video is decoded to a waveform, I expect to see a signal like this:

The problem I'm having right now is, I don't really understand the math involved in demodulating the FM signals from the combined RF signal here.

I understand I can use some basic band-pass filtering to eliminate frequencies outside the bands I'm interested in, which would allow me to separate the audio, video, and digital signals from each other, but after that it's the actual process of de-modulating the FM-encoded data back into their data streams (i.e., turning FM encoded data back to a basic waveform), that I'm not getting right now.

For an answer here, I'm after something that explains this in layman's terms, so that if I was a green student walking into a classroom on signals processing for the first time (which I pretty much am here) I can understand what you're talking about. I still need it to be specific though, because ultimately, I need to understand this well enough so that I can write real code that performs these operations. I'm after a pure digital solution to this problem.


texts: Where can I find an authoritative (peer-reviewed or textbook) reference to sampling-induced beating?

texts:
I presume we are all here well aware about foldback aliasing when sampling signals above the Nyquist frequency; i.e. half the sampling rate.

By contrast, the phenomenon of beating occurs when sampling a signal with a frequency slightly below the Nyquist frequency. You can see and hear some beats on this educational site.

Here is an example of a 60Hz signal sampled at 128Hz resulting in a 4Hz beat tone with a 8Hz beat. The beat tone frequency is calculated as follows: $f_{t}=\frac{f_{s}}{2}-f_{m}=\frac{128}{2}-60=4\,Hz$. The beat has a frequency twice this tone frequency.

In preparation of a paper, I am looking for a peer-reviewed or textbook reference about beating. I browsed through almost 20 DSP textbooks and found nothing.

On the web, I only found a Proceedings of the 2003 American Society for Engineering Education Annual Conference & Exposition session paper by Kostic which included the nice explanatory figure shown below. I am looking for something more authoritative though; peer-reviewed or textbook.


texts: Calculating smoothed derivative of a signal by using difference with larger step=convolving with rectangular window

texts:
I have a signal sampled at $\Delta t: fi(ti=i\Delta t)$ where i = 0..n-1. I want to find the first derivative of the signal: f'(t).

My first thought was to estimate this by a central difference:

$f'(t_i) =\frac{f(t_{i+1})−f(t_{i−1})}{2\Delta t}$

However the signal may have a lot of high frequency noise that may cause quick fluctuations in f'.
I guess the proper thing might be to smooth the signal by convolving with a window function e.g., Hann and then find the derivative from the differences.

A colleague suggested a faster way of finding a smoothed estimate of the derivative:
use a central difference over 2n samples, where n >> 1:

$f'(t_i) =\frac{f(t_{i+n})−f(t_{i−n})}{2n\Delta t}$

This would of course be computationally faster than first convolving with a window function but is it a good solution?

If we form the sum:

$S=2\Delta t[f'(t_{i-n+1})+f'(t_{i-n+2})+..+f'(t_{i+n-1})]$

and expand each derivative by the central difference with step $\Delta t$:

$S=f(t_{i-n+2})-f(t_{i-n})+f(t_{i-n+3})-f(t_{i-n+2})+..+f(t_{i+n})-f(t_{i+n-2})$

all terms except two cancels out:

$S=f(t_{i+n})-f(t_{i-n})=2n\Delta tf'(t_i)$

Therefore:

$f'(t_i)=\frac{1}{n}[f'(t_{i-n+1})+f'(t_{i-n+2})+..+f'(t_{i+n-1})]$

So taking the central difference over 2n samples is equivalent with first convolving by a rectangular window of size 2n - 2 and then taking a central difference over +/- 1 sample.

How "bad" is it to smooth with a rectangular window?

If we take the FFT this will cause "ringing", but we do not need to take the FFT.

Thanks in advance for any answers!


texts: Synchronization Accuracy Required For Communications

texts:
One of the biggest problems in phase modulated communication is discovering when your signal has arrived. High accuracy in this area is important, as even a quarter wavelength error can cause the estimated data bit to be completely wrong. One solution to this problem is to use a preamble code, which is a waveform whose autocorrelation has a very tight mainlobe and low sidelobes.

My question is, are preamble codes actually used in practice? If so, which ones? How is it that such a high accuracy of quarter wavelength is achieved in the real world using a preamble code? If this bound isn't correct, what bound do we really need to achieve synchronization using say BPSK?


texts: Correlation between 2 signals of uneven dimensions

texts:
As a part of my work, I am trying to correlate the audio signal in a video with the pixels of each frame. The steps I follow are:

Audio sampling rate and frame rate of the video are known. So, division of the audio sampling rate with the frame rate gives us the number of audio samples associated with each frame of the video.
That being done, I want to correlate the audio signal (the collection of audio signal amplitudes with 2 channels) with the 2-D pixel intensity of each frame, to basically decipher whether which part of the frame has maximum correspondence with the audio.
But, obviously the 2 signals ( audio signal with 2 channels and the 2-D frame pixel intensity values ) are having different dimensions. In this scenario, how do I compute the correlation between them?

I did encounter the idea of performing canonical correlation, but then again, the unequal dimensions of the 2 signals makes it a problem.


texts: Why do eigenvectors specify the direction of fastest and slowest changes in the image

texts:
I'm reading Szeliski's Computer Vision:Algorithms and Applications. On page 212-213, corner detector is discussed. I can understand upto equation 4.8. From then on, I don't understand how come eigenvectors of A specified the direction of slowest and fastest changes, and why does the shorter eigenvector specified the direction of fastest change?

The book can be found here http://szeliski.org/Book/drafts/SzeliskiBook_20100903_draft.pdf

Any help would be appreciated.


texts:
I'm working on an OFDM-QPSK based system and I'm having some trouble to match the BER of my system with the theoretical curve. My OFDM implementation fills all the subcarriers and have t_symbol s duration and t_cp cyclic prefix duration. I'm trying this to validate my system implementation before go on with my system development. I'm using the following python code for apply noise:

def awgn(self,data,noise_power):
    Es = sum(abs(data)**2.)
    Eb = Es/(self.OFDM_SIZE*2.)
    EbN0_dB = noise_power + 10*log10(self.t_symbol/(self.t_symbol - self.t_cp))
    EbN0 = (10.)**(EbN0_dB/10.)
    N0 = Eb/EbN0
    noise = sqrt(N0/2.)*(random.standard_normal((data.size,)) + random.standard_normal((data.size,))*1j)
    return (data + noise)

I just realize that the normalization done here is not correct, since this signal is post OFDM modulation. But the problem remains.

The equation for the theoretical curve used is:

BER_t   = 0.5*erfc(n_sqrt((10.**(EbN0_vt/10.))/2.))

I had used similar code before on Matlab and it worked fine, but now It's not.

I guess it's not working now because my system don't have the constellation with (1+1j),(1-1j)... I have to use values that are less than one, since my system is implemented in fixed point. Now I'm using 0.75 as an amplitude but it's arbitrary chosen. I tried to correct the energy symbol to get a mean energy equal to one but using this take my signal out of range. I'm thinking about correct the EbN0 before apply on the theoretical BER equation,  but I'm just wondering how.

Can anyone give me the correct way to go?

My central question

Trying to improve my question:
What I need to know is how to apply awgn on ofdm signal when it's subcarrier are modulated with QPSK and how to match it with the theoretical BER plot. My QPSK symbol have an amplitude of 0.75. The above code is to show how I'm doing it now.

update

I run the code changing the FFT subsystem to the scipy.fftpack and the amplitude of the QPSK symbol to 1 instead of 0.75. And the curve match. I also had corrected the theoretical curve to:

BER_t   = 0.5*erfc(n_sqrt((10.**(EbN0_vec/10.))))

So I still think it's a matter of correct the EbN0 in the above equation.

update2

After several testing of the above functions I found thet they are correct. The source of the problem is that I'm dealing with fixed point implementation of IFFT/FFT and I have to handle the channel model with this situation in mind. I'm trying to figure out how to handle this but without a solution yet.


texts: Compensating effects of a system with a known transfer function

texts: Suppose we have a system which we want to know the exact transient times. In ideal case, we can extract the transient times, but in practice it will be affected by another system with a known transfer function. We know  the system with a transfer function  changes the transient times and also the shape of output. How can we obviate its effects on the (main) system?


texts: Remove audible click sound from high frequency hopped sound

texts:
I have a total sound duration of $T$.

I break the duration to several sub-durations and in each sub-duration, I generate several discrete sinusoidal values by the formula $\sin(2\pi ft)$ with a sampling rate of $f_s$.

The frequency used in consecutive sub-durations is different.

Now when I generate sound using these values (using matlab) I get audible click sounds between the point where two sub-durations meet; I think because of discontinuity. How can I remove these audible click sounds?

Edit
the code I am using (roughly)

duration = 0.25;
sampleRate = 44100;
numberOfSamples = (duration * sampleRate);
frequency = 19000:500:21500;
freqHopSamples = [];

for j = 1:length(frequency)
    startIndex = ((j-1)*(numberSubSamples) + 1);
    endIndex = j*numberSubSamples;
    t = (startIndex:endIndex)/sampleRate;
    freqHopSamples(startIndex:endIndex) = sin(2*pi*frequency(j)*t);
end

wavwrite(freqHopSamples, sampleRate, 'sound.wav');


texts: Best method to do piano key pitch detection?

texts:
For example, when playing piano, at the same time, print out the key notes by analyzing the signal of piano sound.

How to do some pre-processing to remove the noise? When calculate with FFT, noise will influence the result?


texts: Finding poles of an abstract transfer function

texts:
When finding the poles of something like the following transfer function, would I be able to write $z=\sqrt[L]{\mu}$ since square roots aren't technically defined on the complex plane?

$$Y(z) = \frac{a_L z^L + a_{L-1}z^{L-1} + ... +a_0}{z^L - \mu}$$

where $a_j,\mu \in\mathbb{C}$, $|\mu|<1$ and L is positive.


texts:
I am trying to understand the function of a Wiener filter. I get that the Wiener filter minimizes the mean square error between the estimated random process and the desired process.

Does this mean that when you have a signal with low and high frequencies and the noise is only in in the high frequencies that the resulting signal only contains low frequencies? What would the transfer function look like in this case?


texts:
I have a chirp signal. In Matlab, the command is $s=chirp(t,f0,t1,f1)$. I make the following sum with Matlab: $y=rand(1,size(s,2))+s;$

Here is my question:
What is the theoritical equation for $y$? I know the expression of $s$.


texts: Removal of shadow from Matlab image for recognition of seeds of plants

texts:
I am working recognition of plant seeds with Matlab, and I am looking for code to remove the shadow in grayscale and binary.

The binarized version of the image is:


texts:
I am using an on-chip SAR adc where for some reason I can only get accurate signal reconstruction and filtering when I sample at almost 100 times my desired sampling frequency. Specifically, the minimum working sampling rate right now is 25KHz but the desired sampling rate is 500Hz for a signal with a BW of 0.5-100Hz, therefore I am still slightly oversampling this signal but the 25KHz is clearly overkill.

I want to downsample this signal to that desired sampling frequency to aid in the design and implementation of digital filters in that passband of .5-100Hz. I know a key thing to look out for is aliasing and because of that you should use a digital anti-aliasing filter but is it pretty much the same concept as an analog anti-aliasing filter where I simply lowpass filter the signal bandwidth that I want to reconstruct after sampling? Also are there any other pitfalls/design tips that I should know to be able to implement downsampling in a DSP? Thanks in advance.


texts: What Is the Relationship Between a Kalman Filter and Polynomial Regression?

texts: What is the relationship, if any, between Kalman filtering and (repeated, if necessary) least squares polynomial regression?


texts: How are the time scales of an exponential decaying autocorrelation related to the duration of peaks?

texts:
I'm reading a biology article (full text here) that analyzed the variation of a cell parameter over time at different stimulus intensities (calcium concentration). The frequency of peaks increases at increasing calcium concentrations, while their duration remains constant, fitting a gamma distribution. Here are the bursts duration distributions for two different calcium concentrations:

By eye those bursts were classified in isolated bursts (with duration $\tau_{burst}$) and clustered bursts (with duration $\tau_{cluster}$). Here are typical traces of two different cells:

The authors computed the autocorrelation from the traces and then estimated $\tau_{burst}$ and $\tau_{cluster}$. Moreover, they estimated the relative weight of the two types of burst at varying calcium concentrations.

Their method was as follow: (1) they calculated the mean autocorrelation function (for every calcium concentration); (2) they described the curve with an exponential decay function; (3) they estimated the weight of the two types of bursts.

As an example, below are shown two mean autocorrelation functions: one of the tracks obtained at 50mM calcium concentration (blue), the other at 200mM (red).

The first is well fitted by a single decaying exponential with time constant $\tau=70s$; the second is fitted by the sum of two exponentials: $c(t)=a_1e^{-t/\tau_1}+a_2e^{-t/\tau_2}$ with $\tau_1=70s$ and $\tau_2=720s$. They concluded that $\tau=\tau_1=\tau_{burst}$ and $\tau_2=\tau_{cluster}$.

In the supplemental document (here) they say:

  The autocorrelation function was computed from the thresholded traces and provided an unbiased estimate of the size of the isolated bursts and clusters, which correspond to the two time scale of the exponential decay in the fit.

My question is: how are the time scales of the exponential decay of the autocorrelations related to the duration of bursts and their relative frequencies?


texts: What is the (approximate) function for amplitude of a plucked string over time? Does it differ between string types?

texts:
Impetus: I'm currently working on my first synthesizer, after many years of playing with them.  I've just added an ADSR envelope to modify the gain, but for some reason it sounds wrong to me.  I've asked around forums whether or not I'm supposed to apply the envelope to the amplitude linearly, and it would appear that's the normal way to do it.  So my guess here is that experience with acoustic instruments has conditioned me to hear drop off differently.  In order to prove/disprove this theory, I would like an answer to the following...

Question: What is the (approximate) function for amplitude of a plucked string over time?  Is it linear, or curved?  And, if it's curved, is that curvature a product of the material, tension, or other variables?

Detailed Explanation: In this context, what I mean by "function for amplitude" is the value by which a simple wave would be multiplied over time to create a reasonably approximate rate of fade out.  I do realize the actual function of amplitude is probably a very complex physics problem, and while I think it would be interesting to cite some sources in that vein, it is not my intention to create a true to physics simulation of a string.


texts: Why do we do oversampling in raised cosine?

texts:
I would like to do know why do we need to do oversampling in a raise cosine filter. There is already a question like this but the answer is very confusing and not sufficient.

I need to know what this oversampling factor does?
The higher is better or the lower?


texts: How do state of the art echo cancellation algorithms deal with variable propagation times of sound due to multiple speakers or reverberation?

texts: I understand how echo cancellation works with a single speaker and no reverberation (using adaptive filtering and freezing the coefficients during double-talk). However, in cases with more than one speaker source or reverberation causing different propagation of the same signal, would this not completely change the combined signal and cause issues with the ability to cancel the echo? Thanks in advance!


texts:
I have a signal matrix which is a $256\times 192$, and I want to calculate the SNR considering that my $259\times 192$ matrix  is an average of a $256\times 192\times 330$ matrix, where $330$ is the number of frames, after a reshape I have obtained a $256\times 192$ matrix. So my code is:

signal = 256x192x256;

signal_average = mean(signal,3);
noise_estimation = signal_average - repmat(mean(signal_average,2),1);
signal_power = mean(abs(signal_average).^2);
noise_power = mean(abs(noise_estimation).^2,2);
SNR = 10*log10(signal_power./noise_power);

Is it correct to use this approach in order to obtain a matrix with different SNR for each frame?


texts: Correlation : How is correlation analogous to dot product?

texts:
Till now I know correlation tells about similarity. I was watching a video lecture on image similarity in which I came to know that correlation is analogous to dot product. And hence correlation of two images is maximum when these images are similar as happens in dot product of two aligned(similar) vectors.

Dot product of two vectors $a=[a_1, a_2, a_3, .....]$ and $b=[b_1, b_2, b_3, .....]$ is given as $$a\cdot b=\sum_{i}a_ib_i$$ While cross correlation of two discrete real signals $x(n)$ and $y(n)$ is given as $$R_{xy}(n)=\sum_mx(m)y(n+m)$$
So now I want to know:

How correlation is analogous to dot product?
When is correlation maximum?
Why?


texts: I was curious if there is a way to extract metadata from wav files. For instance, can you get bpm, key, instrument, etc. from a .wav file on upload?


texts: How to understand multiple peaks in FFT analysis?

texts:
In FFT analysis, it is often seen that there are several peak frequencies with integer times of frequencies. For example, in the figure below, there are peaks at $0.2, 0.4, 0.6, 0.8$.

The FFT figure below is transformed from the time history of force coefficient $C_{Y2}=F_y/(\rho D U_m^2/2)$ and displacement $Y_2$ on $Y$ direction of a cylinder immersed in a disturbed fluid and attached to a spring.

How to understand this pattern?
If you are familiar with vibration, how to explain these vibration data?


texts:
In OFDM system, after IFFT a CP is added at the beginning of the OFDM symbol by copying it from the end of the symbol.My question is that with CP why we use cyclic convolution instead of linear convolution. I would appreciate if someone please could explain the reason of cyclic convolution here.

Many thanks!


texts: Matlab produces two unknown spikes in custom FFT

texts:
I'm making a relatively small Cooley–Tukey FFT in Matlab and I'm noticing unusual spikes in the result compared with Matlab's own FFT.

The figure below shows the signal flow of my program. It's a standard Cooley-Tukey scheme.

My results when computing the FFT of a 16 Hz sinusoid are shown below

And Matlab's own FFT is shown below

Clearly I'm missing something. What could be the cause of the two large spikes in the middle? I conjecture that it has something to do with how the even and odd parts are combined, for instance if there is a discontinuity there. I'm really not sure though.

Any help is greatly appreciated!

The code I'm using is as follows

clear all

% Generate input data sequence and plot
N=128;
f1=16;
num_cycles=2;
fs=f1*N/num_cycles;
x_time=0:1/fs:num_cycles/f1-1/fs;
x=sin(x_time*2*pi*f1);
plot(x_time,x);

% split inputs into even and odd samples and compute fft of each division
X_o=x(1:2:N);
X_e=x(2:2:N);

fft_x_o=fft(X_o);
fft_x_e=fft(X_e);

% Generate base twiddle factor
W32=exp(-1i*2*pi/32);

% Combine fft even and odd with twiddle factors to produce final output
for k=0:N-1
    if k<N/2
        X(k+1)=fft_x_e(k+1)+(W32^k)*fft_x_o(k+1);
    else
        X(k+1)=fft_x_e(k+1-N/2)+(W32^k)*fft_x_o(k+1-N/2);
    end
end

% plot butterfly fft and matlab fft
FFT_xaxis=0:fs/N:fs-fs/N;
figure
plot(FFT_xaxis,abs(X))
title('Butterfly FFT')
xlabel('Frequency')
ylabel('Magnitude')
matlab_fft=fft(x);
figure
plot(FFT_xaxis,abs(matlab_fft))
title('Matlab FFT')
xlabel('Frequency')
ylabel('Magnitude')


texts: How does transformation matrix and camera matrix change when indexing pixels from -1 to 1 instead of 0 to N-1

texts:
I'm trying to warp frames for view synthesis. Specifically, given a frame1 and poses of frame1 and frame2, I'm trying to warp frame1 to frame2's pose/view. For this, I'm indexing pixel locations as $[0,M-1]$ and $[0,N-1]$ for y and x directions respectively. My code is working perfectly (except occluded regions, which lead to holes).
I'm looking at this View Synthesis paper. The code released by authors instead index pixels in the range $[-1,1]$. Directly using transformation and camera matrices is giving invalid transformed coordinates. Can you please explain how I should change the transformation and camera matrices so that it works for pixel indices in the range $[-1,1]$?
More Details:
I've implemented the warping in python using the below equations. Given a pixel location $p_1=[x,y]^T$
$$\hat{p}_1 = [x,y,1]^T$$
$$ P_1 = Z K^{-1} \hat{p}_1$$
$$ \hat{P}_1 = [P_1^T,1]^T $$
$$ \hat{P}_2 = T \hat{P}_1 $$
$$ P_2 = \hat{P}_2[0:3] $$
$$ \hat{p}_2 = K P_2 $$
$$ p_2 = \hat{p}_2[0:2] / \hat{p}_2[2] $$
$p_1$ is the point in frame1 and $p_2$ is the corresponding point in frame2. $P_1$ is 3D world point in the view of frame1 and $P_2$ is the same 3D point in the view of frame2. $T$ is the transformation between views of frame1 and frame2. $Z$ is true depth of point $p_1$ in view of frame1. K is the camera intrinsic matrix. hat ($\hat{p})$ is used for projective coordinates.
I compute the $4 \times 4$ transformation matrix $T$ as follows:
$$T = \begin{bmatrix}
R & t \\
0 & 1 \\
\end{bmatrix}$$
Where $R$ is the rotation matrix and $t$ is the translation.
The $3 \times 3$ camera matrix $K$ is given by
$$K = \begin{bmatrix}
\frac{W/2}{tan(hfov/2)} & 0 & W/2 \\
0 & \frac{H/2}{tan(vfov/2)} & H/2 \\
0 & 0 & 1 \\
\end{bmatrix}$$
Where $H=W=240$ are the height and width of the frames, $hfov=60,\ vfov=45$
In short, what changes should I make to matrices $T$ & $K$?


texts: Why is PSNR used for image quality metrics instead of SNR?

texts:
I have been trying to determine why is PSNR used instead of SNR for image quality metrics?

I have read through several books and most just say that PSNR is used instead of SNR, but don't give a reason as to why.

I also read these posts:

Signal quality metric

and

Difference between SNR and PSNR

But I still feel like they don't explain clearly why PSNR is used instead of SNR.

I know that PSNR measures how the peak signal intensity compares to the mean squared error whereas SNR measures signal power compared to noise, but this doesn't help me reach an intuitive conclusion as to why to use PSNR instead of SNR.

Can someone please explain this?


texts: Biquad filter from Matlab to C using "tf2sos"

texts:
This is a biquad implementation I use:

float my_biquad(float in, float *coeffs, int scalefactor, float *delay)
{
    float out = coeffs[2] * in + delay[0];
    out = out * (float)pow( 2, scalefactor );
    delay[0] = delay[1] + coeffs[1]*in - coeffs[4]*out;
    delay[1] = coeffs[0]*in - coeffs[3]*out;

    // etc...

    return out;
}

It cascades several biquads:

for(int j=0;j<num_biquads;j++)
{
    tmp = biquad(...);
}

Each biquad has its coefficients and a scalefactor. Scalefactor is integer.

On the other hand, I have a Matlab filter:

[B_coeffs, A_coeffs ] = butter(4, 100/(16000/2), 'high');

I need to implement the Matlab filter using the above biquad.

I tried using tf2sos which gave me two biquad sections:

[ss, gn] = tf2sos(B_coeffs, A_coeffs);

but feeding ss coefficients to biquad in C gave different results than Matlab version. Gain gn is left unused, and I don't even have scalefactors.

How do I go from Matlab to this particular version of biquad in C? How to create proper coeffs and scalefactors?


texts: Which of DFT, DCT and DWT transforms is more robust to noise and geometric transformations in image processing

texts:
I am currently studying about image watermarking and I have been testing these 3 domains to embed the watermark at.

My question is, how do these domains hold up to noise addition or geometric transformations, regardless of watermark embeding.

To be more precise, I want to find the relationship between an original image and a transformed one (with the transformations being the ones I mentioned in the title, where geometric transformations include rigid and non-rigid ones) in these domains, and which one holds best at minor and/or major changes.

Also, I would like to learn about the generic pros and cons of these domains, compared to one another, aswell as why it is said to be better to use one instead of the others in watermarking.

Sources would be apreciated. Thank you!


texts:
I am using Audio Precision 2700 to test some test vectors.
I have connected resistor banks as class - d power amp's load, and tap from this resistor using a XLR cable to the AP.
When my power amp is powered up, but no song is playing(song is paused, not muted) i checked fft spectrum of the input channel , i get the following spectrum

You can see the three peaks in frequency above 4k. (Likely, the 3 peaks in frontier band upto 200 hz). Which causes this harmonics?
The test vector is playing at 48khz.
Note : There's a passive filter after the power amp's output having cut off at 30khz.
Please point me some insight. Does the passive filter aliases the input signal, so that i get undesired signal at folded frequency?[This is my initial doubt, because this is what i believe should be the cause, if it's not then some one please tell me why this couldn't be the problem]
Thanks folks in advance...


texts: Does the Nyquist Sampling Theorem hold for a triangular wave?

texts:
Does the Nyquist Sampling theorem hold for a triangular wave produced by a function generator?

For example, a triangular wave with a frequency of 15 Hz sampled at 60 Hz, 180 Hz and 15000 Hz.


texts: Does Radio Communication Have to Account for the Doppler Effect?

texts: I know that radio communications take place in specific frequency bands, and that the doppler effect can cause shifts in frequency when there is relative motion between the source and destination objects. Do moving objects such as planes and rockets have to account for this, or is this movement negligible compared to the speed of light? How about communications with high velocity spacecraft and probes on other planets?


texts: Why isn't the Hartley Transform more widely used?

texts: I was doing a bit of research on FFT and I came across this paper. According to it, computing the Fourier Transform using FHT is significantly faster than other methods, even for complex data. Yet, from what I gathered so far, the FHT method is rarely used. Why is that?


texts: Prefactors computing PSD of a signal with numpy.fft VS. scipy.signal.welch

texts:
The power spectral density St of a signal u may be computed as the product of the FFT of the signal, u_fft with its complex conjugate u_fft_c. In Python, this would be written as:

import numpy as np

u = # Some numpy array containing signal
u_fft = np.fft.rfft(u-np.nanmean(u))

St = np.multiply(u_fft, np.conj(u_fft))

However, the FFT definition in Numpy requires the multiplication of the result with a factor of 1/N, where N=u.size in order to have an energetically consistent transformation between u and its FFT. This leads to the corrected definition of the PSD using numpy's fft:

St = np.multiply(u_fft, np.conj(u_fft))
St = np.divide(St, u.size)

On the other hand, Scipy's function signal.welch computes the PSD directly from input u:

from spicy.signal import welch

freqs_st, St_welch = welch(u-np.nanmean(u),
     return_onesided=True, nperseg=seg_size, axis=0)

The resulting PSD, St_welch, is obtained by performing several FFTs in segments of the array u with size seg_size. Thus, my question is:

Should St_welch be multiplied by a factor of 1/seg_size to give an energetically consistent PSD? Should it be multiplied by 1/N? Should it not be multiplied at all?

PD: Comparison by performing both operations on a signal is not straightforward, since the Welch method also introduces smoothing of the signal and changes the display in the frequency domain.

Information on the necessity of the prefactor when using numpy.fft :

Journal article on the matter


texts: What kind of Modulation causes "drops" in FFT Power Spectrum

texts:
I am analyzing a signal with FFT and I found an effect, that I do not understand. I am expecting to have an amplitude modulated signal. As far as I know, if you compute the power of the FFT of an amplitude-modulated signal, you observe sidebands to the left and right of the carrier frequency.

E.g., if you have a sinusoidal carrier at frequency $w_c$, and you modulate the carrier's amplitude with frequency $w_m$, you see sidebands appearing at $w_c \pm w_m$.

However, in my signal's power spectrum, I see the following. Note, that the x-axis shows not the absolute frequency, but the delta from the carrier.

At the frequency, where I would expect peaks in the spectrum, I see drops. How can that be? What does that tell me about the signal?
If I look at the signal in the time domain. I can see, that amplitude modulation is actually happening!


texts:
I have a vector of zeros and ones with length ~2000. I need to detect and filter out the once, that appear periodically. There may be none to 3 such periodic signal's components. Typical periods are in range 10-90, but due to measuring noise, some spikes may be shifted by 1-3 steps in any direction.

I assume there is a technique that I can set a threshold on the energy of periodical spikes, and filter them out causally. I'm not asking for a ready solution (though it would be very nice), but I would be grateful for a link to a relevant and simple example, preferably in Matlab and using FFT.


texts:
I am trying to compute the phase spectrum of the signal
$$
s(t)=\frac{A}{\pi}\left[H\left(t+\frac{\tau}{2}\right)-H\left(t-\frac{\tau}{2}\right)\right]
$$

The Fourier transform is
$$
S(j\omega)=\frac{A\tau}{\pi}\cdot\frac{\sin\left(\omega\displaystyle\frac{\tau}{2}\right)}{\omega\displaystyle\frac{\tau}{2}}$$

My answer for this would be $\phi_s(\omega)=0$, since the imaginary part of $S(j\omega)$ is zero, but according to my textbook, the answer is

$$\phi_s(\omega)=\begin{cases}0&S(j\omega)>0\\\pm\pi&S(j\omega)<0\end{cases}$$

How can I get to this answer and interpret it correctly?


texts: Spectrum Visualization of TX AD9361 Simulink Model

texts:
I'm using the SIMULINK AD9361 simrf support models, more specifically i use an AD931 TX model to generate a 20MHz RF LTE signal around 2GHz Fc
the issue is when i use a spectrum analyzer scope i see the 20MHz LTE signal is around 0 DC ? it should be around 2GHz, what's the catch?
setup

Spectrum


texts: Why am I seeing a negative frequency on FFT graph in Gnu Radio Companion?

texts: I have a pretty simple flow here in GRC, but for some reason I'm seeing a negative frequency mirroring the 10khz signal that I'm generating. What's causing this negative frequency to show up in the FFT graph?


texts:
The OQPSK modulator in Matlab produces an output sequence as follows:
[0.707
0.707+0.707i
-0.707-0.707i
...and so on
]
This is the sequence and I am trying to view the frequency spectrum of this sequence. With the sampling rate being 2 MHz, I expect to see a peak centered at 0 with a bandwidth of 2 MHz. When I upsample the above sequence and multiply this sequence with a carrier of say 2.4 GHz, I should have 2 peaks centered at -2.4 GHz and 2.4 GHz and zeros at other frequencies. But somehow, the spectrum is distributed between -2.4 and 2.4 GHz. Am I making a mistake in plotting the frequency spectrum or is there something wrong in what I understood?


texts: Are there any methods of detecting, measuring, and/or estimating I/Q imbalance (from, say, an upstream quadrature mixer and/or ADC), without using or having a controlled test signal (other than whatever the RF antenna is currently picking up).


texts: Need a similarity metric that describes these two curve as highly similar

texts:
I have a large number of time series data and I need to be able to compare the similarity of the curves to a reference curve. The reference curve in question is shown in red in the figure, and I'd like to compute a similarity metric for the blue curve that tells me that these things are highly similar.

The reason that I want these things to be similar is that they are both characterized by a largely monotonic rise, followed by a near plateau at some upper value, and then followed by a largely monotonic drop back down to the initial level. I have computed some similarity measures such as the Frechet distance and DTW distances, using the following snippets of code

sim_score = similaritymeasures.frechet_dist(Z, Zref)
sim_score2, _ = fastdtw(Z, Zref, dist=euclidean)

The problem is is that the scores computed using these techniques do not tell me that these curves are similar. For example, the fastdtw score is shown in the title of the figure, and this score is actually quite large when compared to other comparison scores for curves which I would judge by eye to be much more dissimilar. For example, the following curve has roughly half the fastdtw distance, despite it not having as similar of a shape.

The issue seems to be that the temporal overlap is still being weighted much more heavily then more shape-related attributes. Now it's perhaps possible that a more nuanced understanding of what's possible with DTW would give me the outcome I need, perhaps using a more fully-featured DTW package like https://pypi.org/project/dtw-python/, but I'm unfamiliar with many of the options available in that package. Or maybe there's another method that would be better suited. Of course, I could brute force this by creating a piecewise continuous curve with a bunch of variable parameters describing a step-like shape, optimize them against the target curve, and take the best error match as the distance metric. If the optimization doesn't get trapped in a local minima, then this should work, but it's more involved and messier than I'd like to do if I can avoid it.

Would appreciate any other ideas you all might be able to suggest. Thanks in advance!


texts: Is it possible to do deconvolution with two data sets that have different sampling rates?

texts:
I have some terahertz spectroscopy time series data, a reference set with 2048 data points taken every 0.0521 picoseconds, and the sample data set with 544 data points taken every 0.0781 picoseconds. I'm using Matlab to take a FFT of both sets, with zero padding on the sample set up to 2048 points, and I'm supposed to do a deconvolution on those transformed sets on a specific range of frequencies, but because the sampling rate is different for the two, these correspond to different starting points in the data sets, and the lengths won't match up either. When plotting the two DFTs, they only line up properly when accounting for the different sampling rates when spacing the frequency axis.

Is there something I'm doing wrong, or how should I do the deconvolution?


texts: How to find analytic description of filtered signal

texts:
I am looking for an exact analytic description of a filtered signal. I have an electronic circuit whose input is a monoexponential decay. First (1) the signal gets filtered by a simple RC-Lowpass. Then (2) a second order Bessel Low Pass Filter is applied.

The respective transfer functions are:
\begin{align}
TF_{\rm RC}(s) &= \frac{R}{s+\frac{1}{R C}}\\
TF_{\rm Bessel}(s) &= \frac{a_1}{a_1+b_1s+s^2}
\end{align}

The input signal is:
$$I_{\rm in}(t) = \frac{V}{R_1 + R_2} + R_1\cdot V_{\rm cmd}\cdot e^{-t\frac{R_1 + R_2}{CR_1R_2}}\frac{1}{R_2(R_1 + R_2)}$$

My idea was to calculate the output in the $s$-domain:
$$ V_{\rm out}(s) = I_{\rm in}(s)\cdot TF_{RC}(s)\cdot TF_{Bessel}(s) $$
and then apply the inverse laplace transform to obtain $V_{\rm out}$(t).
The inverse transform is what troubles me. I tried to calculate it with MATLAB, but it doesn't return an analytic function. I guess this comes from the complex poles of the bessel filter. Could you give me a hint on how to calculate the output?

Edit:
$$ TF(s) =  \frac{6.092e12}{s^3 + 5.367e04*s^2 + 1.038e09*s + 6.074e12}$$

$$ I_{in}(t) = \frac{10e-3}{500e6 + 10e6} + 500e6*10e-3*e^{-t*\frac{500e6 + 10e6}{10e-12*500e6*10e6}}\frac{1}{10e6*(500e6 + 10e6)} $$

$$ TF(s)*I_{in}(s) = \frac{4.9e-21*(1.2e24*s + 2.5e26)}{s*(s + 1.0e4)*(s^3 + 5.4e4*s^2 + 1.0e9*s + 6.1e12)} $$

PFE yields:

$$ TF(s)*I_{in}(s) = \frac{1.1e-9 - 4.8e-11i}{s + 2.2e4 - 1.1e4i} + \frac{1.1e-9 + 4.8e-11i}{s + 2.2e4 + 1.1e4i} + \frac{6.6223e-7}{s + 10200.0} + \frac{-6.6441e-7}{s + 10235.0} + \frac{1.9666e-11}{s}  $$

Inverse Laplace transform of the upper term will yield a complex result. I know that I can simplify the first two terms by bringing them to a common denominator and thereby get rid of the complex numbers. The inverse transform is however still complex.

Edit: calculation with symbolic parameters yields:
$$ H(s) = \frac{a*s+b}{c*s^2+d*s+f} $$
$$ H(t) = a*dirac(1, t) + (2*b*sin((t*(- d^2 + 4*c*f)^{1/2})/(2*c))*exp(-\frac{d*t}{2*c}))\frac{1}{(4*c*f - d^2)^{1/2}} $$

I don't quite understand where the dirac comes from, but this is basically the same result as doing the calculation with complex numbers and the taking the real part.


texts: Variance of periodogram estimate of the power spectrum

texts:
I have been reading chapter 13.4. ("Power Spectrum Estimation Using the FFT") of the Numerical Recipies Book.

Some things related to the expectation value of the "periodogram estimate of the power spectrum" became not clear to me though.

Background

Suppose we have an equally-spaced $N$-point sample  $ c_0 ... c_{N-1} $ of the  discrete signal $ c(t)  $ and its Fourier transform $$ C_k = \sum_{j=0}^{N−1} =
c_j \exp{
2πijk/N}
 $$ $$k = 0, . . . , N − 1 $$
Then we can define the "periodogram estimate of the power spectrum" $P(f_k)$ as
$$ P(f_k) = \frac{1}{N^2}
[{|C_k|^2 + |C_{N−k}|^2]}
$$

with $k = 1, 2, . . . , \frac{N}{2}− 1$ and $f_k = \frac{k}{N\Delta}$ for $k = 0, 1, . . . ,\frac{N}{2}$ ($\Delta$ is the sampling interval).

So far everything is clear.

Now the book says

  "... the variance of the periodogram estimate at a frequency $f_k$ is
  always equal to the square of its expectation value at that frequency.
  In other words, the standard deviation is always 100 percent of the
  value, independent of N!"

My questions

What is the expectation value of the periodogram estimate at a frequency $f_k$? Why (can you show me how to derive it)?
The expectation value doesn't change with $N$, so changing the length of the input doesn't improve the power spectrum estimation. But if I have a non-integral number of periods in my sample, more samples (longer sampling time) means less weight on the end of the signal (where there is a disconnect) so less smearing of frequency components measured due to the disconnect between the beginning and the end of the signal. So in this case a higher $N$ means a better estimation? Where am I wrong?
If the SD is always 100 percent of the expectation value, doesn't that make this simple periodogram estimation (using just a single DFT of $N$ samples) pretty much useless? Should one always opt for Welch's or Bartlett's method e.g.?

Thank you very much in advance!


texts: Meaning of concentric circles on the spectrogram

texts: I've got this unusual spectrogram and I am really curious what these circles could mean. Has anyone ever seen anything like this? What could be the source of the signal?


texts: Calculating phase response of filter with latency

texts:
I have an LTI filter that I want to treat like a black box. It has a latency of 24 samples.

This is what I'm doing (which works for a filter with no latency):

Send unit impulse through my filter
Capture first 16,384 samples output
Perform FFT
For each of first 8192 complex values, take phase

Doing the above gives the red line shown below.

If I replace step 2 with:

Skip first 24 samples and then capture first 16,384 samples output

Then I get the blue line shown below.

The green line represents the actual phase of the signal.

Bode plot of phase vs log frequency (10Hz to 20kHz)


texts: What is the meaning of Mutual Information beyond the numerical calculation?

texts:
Beyond the raw equation for calculating mutual information, what does it mean in physical terms?
For example:
From Information Theory, we know that entropy is the smallest loss-less compression scheme that we can use on a alphabet with a specific probability distribution.

What would that mean in terms of Mutual Information?

Background:
I'm trying to calculate the mutual information of uni-gram words and determine which of two books did they come from.

essential $$ I(book; word)$$


texts: White noise and its phase after narrowband bandpass filtering

texts:
Suppose $w(t)$ is a white noise with uniform spectral density $N_0/2$.  It is bandpass filtered to get a narrowband white noise $n(t)$ for $f_c - B_T/2 \leqslant |f| \leqslant f_c + B_T/2$. We write this in polar form by splitting into inphase and quad-phase wrt $2\pi f_ct$. $n(t) = n_I \cos(2\pi f_ct) - n_Q \sin(2\pi f_ct)$ such that we get a phase $\phi_n = \tan^{-1} (n_Q/n_I) $.
Is this phase uniformly distributed on $[ 0, 2\pi ]$? If yes, why?


texts: Audio codec analysis in MATLAB: crosstalk

texts:
I am currently making an analysis of a codec system that can encode up to 2 wave files into 1 file and decode it. Here is a simple illustration:

L+R ---> ENC --->L'+R'
with L'=L+∂1 & R'=R+∂2

The codec is not lossless and we are analyzing if there is any leakage from one source to another, some kind of crosstalk, if you will.

So what I have to do is find out if ∂1 has traces of R & ∂2 has traces of L.

I feel a general approach would be to subtract L/Ldec to get ∂ and then compare that ∂ to R.

I have done some reading on correlation but it's all a little vague at the moment.

So, onto my question(s):

Is this possible with limited knowledge of DSP (reading up but it's a
lot & not that easy)?
If not, on what should I do some reading (From what I've gathered so far: correlation & Pearson's coefficient)?
Is there some sort of standardized test for this - I presume there is
but I can't find any!

Even though code snippets would be handy since I only just started to use MATLAB & still am somewhat inexperienced, the main point of this topic is understanding so I can implement it. If someone else implements it I won't learn!


texts: power spectrum of passband of ofdm is it right shape or not?

texts:
I tried to get power spectrum of ofdm passband signal.
Below is the whole matlab code for that with BPSK+ofdm.
and I got the power spectrum as below.
I wonder that is that right form or not?

Nfft=64; % number of sub carriers
Ncp=Nfft/4; % number of cyclic prefix
Nuc=52; % number of used carriers

d=2*(randn(1,Nuc)>0.5)-1; % bpsk mapped data

X=zeros(1,Nfft);
X([2:Nuc/2+1,(end-Nuc/2)+1:end])=d; % only used sub carrier
x=ifft(X,Nfft); % time domain
sbb=[x(end-Ncp+1:end) x]; % to add cyclic prefix

%% passband

fc=4; % carrier frequency
sps=16; % samples per symbol
Tsym=1; % symbol duration time
Ts=Tsym/sps; % sampling time
fs=1/Ts;

Tmax=(Nfft+Ncp)*sps*Ts-Ts;
t=0:Ts:Tmax;

%return

% I-ch up-sampling
I=real(sbb);
I=repmat(I,sps,1);
I=I(:).';
sI=I.*cos(2*pi*fc*t);

% Q-ch up-sampling
Q=imag(sbb);
Q=repmat(Q,sps,1);
Q=Q(:).';
sQ=Q.*-sin(2*pi*fc*t);

s=sI+sQ; % transimt signal

%% plot power spectrum

f=(0:(numel(s)-1))/numel(s)*fs;

figure
plot(f,10*log10(abs(fft(s))))
%plot(f/1e6,10*log10(abs(fft(s))))
xlabel('frequency[Hz]');
ylim([-40 20]);
title('full span');
grid on

figure
plot(f,10*log10(abs(fft(s))))
%plot(f/1e6,10*log10(abs(fft(s))))
xlabel('frequency[Hz]');
ylim([-40 20]);
xlim([3 5])
title('main band');
grid on


texts: Generation of a stream of data bits 0s and 1s at a certain bit or data rate

texts: This is a very basic question ,but basics important, : How is a stream of data bits (zeros and ones that express a certain information) and at a certain rate being generated? what is involved in terms of processing and possibly hardware?  a clock? a quantizer? a sampler? or what exactly?  I really would like to get a brief idea


texts:
I am having some trouble understanding how a digital signal gets converted back to an analog signal.

Why is the "sample and hold" effect a problem?


texts: How to calculate the luminance of an image?

texts:
Is there any method for calculating luminance of an image?

As from wiki got to know Luminance is an amount of energy perceived by the observer and LUX is different from brightness.

But still is there any way to calculate LUX value of an image using any of image processing technique.


texts: Number of discrete samples required for the longest wavelength

texts:
I am trying to understand the effect of the critical Nyquist frequency when applying the Goertzel algorithm for estimating the power spectrum of a discrete signal.  (Goertzel doesn't really matter, FFT applies as well).  Assume that I generate $N=512$ discrete values of the sine wave $x[n]$, where $n=1,2,\ldots,N$ using the function

$$x[n]=\sin(30 \cdot 2\pi n/N) + \sin(60 \cdot 2\pi n/N)  + \sin(120 \cdot 2\pi n/N) \\  + \sin(260 \cdot 2\pi n/N) + \sin(330 \cdot 2\pi n/N)$$

After applying the Goertzel algorithm, there will be 512 power spectrum values at $k=1,2,\ldots,512$. Looking at the results, it appears that in the upper half of the scale ($257 \leq k \leq 512$), the real (vs. imag) power spikes for frequencies 1/260 and 1/333 are negative.  In addition, in the lower half of the range ($1 \leq k \leq 256$), the power value for $k=60$ is zero, and there are false positive spikes at spurious frequency values.  However, if I double the length of the generated signal to 1024, and only apply the Goertzel algorithm to $k=1,2,\ldots,512$, all of the power spectrum values at $k=30, 60, 120, 260$, and $330$ are positive and there are no false positive values below 512 and the value at k=60 is non-zero and positive.  My understanding of the Nyquist requirement is that the power value at $k=1$ is unreliable.

So what I believe I have observed is that if you want to use the real (vs imag) power spectrum values at frequencies between 1/2 and 1/512, you need to provide a signal with twice (1024) the number of samples.  Another way of saying this is that if you have a discrete signal with 1024 samples, you can only determine the power spectrum at frequencies greater than 1/512, i.e., 1/511,1/510,...,1/2.

Is there a law or equation that states that the sample size $N$ needs to be at least twice as long as the longest wavelength to be assessed?  Or, is this really the meaning of Nyquist, i.e. $2N$.


texts:
Suppose I have a real (physical) dynamical system with some sensors and actuators, and I also have an idealized state-space model of this same system.  How, in general, can I adjust the model to match the real system?

My experience is with frequency-domain analysis.  My first reaction would be to measure swept-sine transfer functions, and to calculate parameters of my system from the location of resonances or poles in the measured data.

Is there a better, more general technique for system identification in the state-space representation?


texts: How to derive the $\log$ of fraction of two Gaussian distribution

texts:
I am designing a function $A(x)$ that map $x\to \{0,1\}$. The function $A$ can be expressed as
$$A(x)=\begin{cases}
1 \text{  if $p_1(x) \ge p_2(x)$}\\
0 \text{  otherwises}
\end{cases}
$$
where $p_1$ and $p_2$ are Gaussian distribution:

$p_1(x)=\frac{1}{\sqrt{2\pi}\sigma_1}\exp\left( {-\frac{(x-\mu_1)^2}{2\sigma_1^2}} \right)$;

$p_2(x)=\frac{1}{\sqrt{2\pi}\sigma_2}\exp\left( {-\frac{(x-\mu_2)^2}{2\sigma_2^2}} \right)$.

Could you help me to derive a short form of $A(x)$ ?

This is my solution. However, it maybe wrong

I derived the bellow eq.
$$A(x)=\begin{cases}
1 \text{  if $\log \frac{p_1(x)}{p_2(x)}>0$}\\
0 \text{  otherwises}
\end{cases}
$$
where $\log(e)=1$

We have
$$\log \frac{p_1(x)}{p_2(x)}=\log \frac{\sigma_2}{\sigma_1}+\left( {\frac{x-\mu_2}{\sqrt{2}\sigma_2}} \right)^2-\left( {\frac{x-\mu_1}{\sqrt{2}\sigma_1}} \right)^2$$

Hence, my final solution is
$$A(x)=\begin{cases}
1 \text{  if $\left (\log \frac{\sigma_2}{\sigma_1}+\left( {\frac{x-\mu_2}{\sqrt{2}\sigma_2}} \right)^2-\left( {\frac{x-\mu_1}{\sqrt{2}\sigma_1}} \right)^2 \right) \ge 0$}\\
0 \text{  otherwises}
\end{cases}
$$

Thanks in advance


texts: If a filter has a non-linear phase lag will the output signal be distorted?

texts:
EDIT: From here (https://www.mathworks.com/help/signal/ref/butter.html)
it appears that an Butterworth filter has a non-linear phase lag at the normalized frequency of 0.5 onwards.
Does that mean that if the signal has multiple frequencies beyond 0.5, its output after passing through this filter will be distorted by these frequencies?

NOTE: To help clarify things, this question had been modified from the original one below.

A butterworth low pass filter would have a large phase lag for signals passing through at higher frequencies.
[EDIT:
The assumption is, for some filters to work, you'd have to do some kind of transform on the input (into the s domain), multiply that by the filter's s-domain transfer function (such as ones from here ece.uic.edu/~jmorisak/blpf.html), then to see the time-domain filtered output, the resulting output has to undergo an inverse transform.]
I envision that this means that the output signal coming out of the filter, when converted back into time domain, would come out distorted due to this lag.
[EDIT: The assumption is, phase lag is the same as a time delay between the input and the output, and if the input has more than one frequency component, the different components would be delayed differently, and end up being combined with another signal output at a later point in time.]
Would there be similar distortions (i.e.: as if there were a phase lag) for an IIR time-domain digital filter?


texts:
I want to up-convert a baseband sinusoidal signal of frequency $4000 \, \rm{Hz}$ using a carrier frequency of $20\, \rm{KHz}$. I read that one can up-convert a baseband signal to a passband signal using the following.

But after convolution with a carrier cosine, I get a $16\, \rm{KHz}$ signal. What could be the problem?

I used following Python code to do convolution with carrier frequency: samples * cosineSignal.


texts: Histogram Equalization Ignoring zero in Python

texts:
I wrote a small python code for histogram equalization (ignoring zero values) for an image but it is taking too long to run. Code is written with reference to wikipedia article on histogram equalization.

def hist_eq_updated(img):
    img_flatten = img.flatten()
    img_flat_nnzero = img_flatten[np.where(img_flatten!=0)]
    image_histogram, bins = np.histogram(img_flat_nnzero,256)
    #print(image_histogram,bins)
    cdf = image_histogram.cumsum() # cumulative distribution function
    #print(cdf)

    eq_img_flatten = np.zeros_like(img_flatten)
    cdf_min = min(cdf)
    cdf_max = max(cdf)
    for i in range(len(img_flatten)):
         val = img_flatten[i]
         ind = 0
         flag=False
         for j in range(len(bins)):
             if val==0:
                 flag=True
             if (val>=bins[j]) and (val<=bins[j+1]) or flag==True:
                 break;
             ind = ind + 1
         if flag==True:
             eq_img_flatten[i] = 0
         else:
             eq_img_flatten[i] = (((cdf[ind] - cdf_min)/(cdf_max-cdf_min))*254) + 1
    return eq_img_flatten.reshape(img.shape)

I am new to image processing in python and not familiar with any other "pythonic" way to do same. Thanks in advance !!!

EDIT - Sorry, I forgot to mention that image size is ~4000 * ~4000, image type is float32.


texts:
Why is $\delta[an]=\frac 1 a \delta[n]$ in discrete time? Prove.

Hi, I want to prove it but I don’t know on what facts to rely. In continuous time we have integral and dt in that integral that turns to be dt/a for continuous time delta function. But in discrete time, we don’t have dt. So how to prove it then?


texts: When Does $ {L}_{1} $ Regularization Give a Sparse Solution?

texts:
I was maximising a likelihood function, which is convex. I know that the system has a K-sparse solution. I wanted to know the conditions (or some sufficient conditions) on the likelihood function under which l1-regularisation is guaranteed to give me the sparse solution.

In other words, I am trying to solve the following problem

$$\max_{x \in \mathbb{R}^n} f(x),$$

where $f$ is a concave function. However I know that the true solution is K-sparse. I wanted to know the conditions under which the following optimisation returns the correct solution

$$\max_{x \in \mathbb{R}^n} f(x) - \lambda ||x||_1.$$

I found online that if $f=||y-Ax||_2^2$, then we need $A$ to satisfy a property called the restricted eigenvalue property. I wanted to know if anything is known in the general case or even special cases where such guarantees were known.


texts:
In a problem I'm working on, I have a signal that is sampled 256 times for 3.5 microseconds.

After performing a Discrete Fourier Transform on the data, what is the frequency step between the calculated values of Cn? The prompt states that this is also called the bandwidth of each point in the power spectrum.

How can I go about solving this problem?

Thanks


texts:
Assume an arbitrary (discrete) signal that is periodic and known over a whole period. I need a way to select a characteristic point along the signal such that I can always retrieve it even when the signal is phase shifted.

Candidates could be zero crossings or extrema or inflections... provided they exist and are unique or can be discriminated from others.

For example, the position of the maximum would work well for a triangle wave, a little less accurately for sinusoid and would raise an ambiguity for $\sin x+\frac13\sin 3x$ as this one has two maxima of equal height.

The method should be robust to noise and work for relatively smooth signals.

Do you know any general feature/formula having the desired properties ?


texts: Applying DFT to OFDM signals above Nyquist rate

texts:
Consider an OFDM signal which has $N=8$ subcarriers and is band-limited to $[-B/2,B/2]$, $\sum_k x_k e^{j2\pi kf_0t}$ where $f_0=\frac{B}{N}$. For this case the Nyquist rate is $2\times B/2=B$, i.e., the sampling interval $T$ (distant between two consecutive samples) is $T=\frac{1}{B}$ and the $l$'th sample is obtained as

$x(lT)=\large\sum_k x_k e^{j2\pi kf_0 lT}=\sum_k x_k e^{j\frac{2\pi}{N}kl}$

which is the same as IDFT formulation and the frequency components $x_K$ can be obtained using these $x(lT)$ and DFT. $x_k$'s are the symbols in frequency subcarriers $[-4f_0,-3f_0,-2f_0,-f_0,0,f_0,2f_0,3f_0]$

Now consider a sampling frequency $B^{'}$ bigger than Nyquist ($B^{'}>B$). In other words, a sampling interval $T^{'}<T$. In this case, we will obtain $N$ new samples, let's call them $x(lT^{'})$'s:

$x(lT^{'})=\large\sum_k x_k e^{j2\pi kf_0 lT^{'}}=\sum_k x_k e^{j\frac{2\pi}{N}\frac{B}{B^{'}}kl}$.

My question is:

Is there a way to obatin $x_k$ (symbols in subcarriers $[-4f_0,-3f_0,-2f_0,-f_0,0,f_0,2f_0,3f_0]$) from these new set of samples $x(lT^{'})$.


texts:
Suppose we have a sequence of still images each of which has been contaminated by some particles(ex, dust/sand/smoke) making the images very poor in certain areas.
What approach would be best to teach image regeneration using multiple frames? The simplest technique is to simply find a way to detect what parts of the image are contaminated and uncontaminated and pull uncontaminated sections from each frame.
For spec based contamination you could probably use one of the classic denoising techniques but I'm thinking of having whole sections be contaminated making these single image local approaches not work great.
I am thinking some form of deep learning restoration(advice on architectures would be great), but I am open to more classic signal processing techniques.


texts: How to go about implementing a program to measure the "intensity"/"saturation" of a color in an image?

texts:
I have a task of developing an application to measure the intensity of a color in an image. I am researching about how to go about it.

  Intensity refers to the purity of a hue. Intensity is also known as
  Chroma or Saturation. The highest intensity or purity of a hue is the
  hue as it appears in the spectrum or on the color wheel.
  Source

The image will be RGB but I can convert it HSV or any other color space.

I am not asking for code, but just high level steps on how to do it.

I am thinking of using OpenCV library for Java.


texts:
For a typical conventional control loop (no delay, measurement transfer function is 1) the closed-loop transfer function is
$$H(s) = \frac{G_{comp}G_{process}}{1+G_{comp}G_{process}}$$
With an actuator delay of $T$, the transfer function is as follows.
$$H(s) = \frac{G_{comp}G_{process}e^{-Ts} }{1+G_{comp}G_{process}e^{-Ts} }$$
However, suppose that the delay is in the measurement, i.e. in the feedback loop, the transfer function would then be as follows.
$$H(s) = \frac{G_{comp}G_{process}}{1+G_{comp}G_{process}e^{-Ts}}$$
I have the following two questions.

How do I analyze the stability of the last case? As far as I understand, I cannot use the margin and allmargin Matlab functions as the delay is not in the open loop transfer function ($G_{comp}G_{process}$).
Now assume that I'm stuck with a total delay of $100 \; \mu s$. Say, I can split the delay, $50 \;\mu s$ in the measurement and $50 \; \mu s$ in the actuator. Can I improve the performance, somewhat, of my controller compared to having an actuator delay of $100 \; \mu s$?

Edit Corrected the 2nd transfer function, with the actuator delay only. Thanks to @TimWescott for finding my mistake.


texts:
Possible Duplicate:
  Calculate Phase angle between two signals i.e. digital Phase meter

I want to calculate the phase difference between the two channels(signals) of a stereo audio file.(uncompressed data i.e. PCM samples).

I want to determine whether whether two signals are in phase or not i.e. kind of phase meter ?
Not the exact phase difference.

Can you please suggest a algorithm to calculate the phase difference.


texts:
I have tried to solve, but do not know if the answer is correct or not.

A person has a periodic voltage input to a circuit. The input repeats itself every 0.02 seconds i.e. the fundamental period is 0.02 seconds.
The person measures the input with a voltmeter and finds that the input may be approximated by:

Vin(t) = V0*exp(a*t) ; 0 sec <= t < 0.02 sec ;  V0= 2 Volts a = -100sec^-1

Graphically the function may be represented by:

A)Write a mathematical expression containing an integral that may be used for calculation of the Fourier-coefficients of Vin(t).

B)Evaluate the integral analytically for the DC-term (often referred to as c0 ). Looking at the graph, how can you check if your answer is reasonable?

I know I need to use these formulas:

I've tried to calculate a0 and here I got approx 2,
according to this matlab code:

v0 = 2;
a = (-100^-1)
T = 0.02
a_0 = (1/T)*int((v0*exp(a*t)),t,0,T)

a0 calculated to 1.9998

I stuck with the an and bn, because i'm not sure if that is correct or not...

an = ((2/T)*int((v0*exp(a*t))*cos(n*w0*t),t,0,T))

an calculated to (20000*exp(-1/5000)*(exp(1/5000) - 1))/(100000000*pi^2 + 1)

bn = (2/T)*int((v0*exp(a*t))*sin(n*w0*t),t,0,T)

bn calculated to (200000000*pi*exp(-1/5000)*(exp(1/5000) - 1))/(100000000*pi^2 + 1)


texts: Digital filter design of time series for specified frequency response function

texts:
I am currently working with vibration measurements in structures. In the netherlands there is a guideline for verifying vibration measurements for damage to machinery. This is the so-called "SBR Trillingsrichtlijn". In this guideline a frequency weighting function is specified to modify the time series from vibration measurements:

with small f being the frequency in Hz, f_0 a constant frequency of 5.6 Hz and v_0 a predefined constant velocity of 1 mm/s. To illustrate the filter behaviour i show it graphically here as well (axes are linear):

In post-processing this is fine as i can perform an FFT on the time series data, apply this weighing function in the frequency domain (through multiplication with the FFT transformed time series) and then do an inverse FFT to get the adjusted time signal data.

My question is the following: the equipment we use to measure vibrations does this weighing internally using a digital filtering on the time series in real-time. To verify the data i would like to design a digital filter corresponding to the weighing function above, can anyone point me to a resource how to approach this? From what i have looked up so far on IIR / FIR filter design usaully the examples start from a transfer function in the s-domain or z-domain (which correspond to a laplace transform in discrete or continuous time) but now i only have a frequency response in the fourier / frequency domain.

EDIT: Based on suggestion of @Hilmar it appears a simple first order high pass function. i made a python implementation to check the method of FFT and the highpass filter. shown in the graph here. Strangely the frequency response of the highpass filter and SBR specified response match exactly but the implementation using the highpass has some artefacts whereby non existant frequencies are generated with non-zero amplitudes:

this can be seen in the green dotted line which is a mock signal filtered in python using the scipy package.