texts:
I have a doubt. We calculate metrics such as PSNR, SNR, and MSE to prove that the proposed algorithm is better than its counterparts. To do this we need Original image (Noise free), Noisy image (we add this noise) and then Filtered image. Now we compare original and filtered output to estimate these parameters.

My question is, how to estimate these parameters when we have only the noisy image? That is in case of an ultrasound image, speckle noise occurs in the original image itself. We don't have a noisy free images in case of ultrasound. If we do as usual (as stated above), we will get poor results (PSNR,SNR & MSE) even though the subjective quality of the image is improved considerably.

Well, one way to do is to add extra noise to this noisy (original) image and apply our filter and show that the added noise has been removed efficiently. I want you to suggest any other ways to prove that my proposed algorithm is performing better than its counterparts.

Example:

              Original Image                       Filtered Image


texts: Relating transfer functions with step responses

texts:
Relate the transfer function to its' corresponding step response.

First, I tried setting up the poles and zeros of the transfer functions. This helped a bit since I know that $G_A (s)$, $G_B(s)$ and $G_C (s)$ should all have oscillatory step responses (because of conjugate poles). The problem arose when trying to determine the differences between $G_A (s)$ and $G_B (s)$, and we thought that there has to be a better way of understanding this. Working with poles and zeros, we didn't really understand what we were doing.

So, I'm now thinking that I need to use the fact that the inverse laplace of the transfer function equals the impulse response, and the integral of the impulse response equals the step response. This is how we relate the step response with the transfer function, but we don't know how to work with it.

Edit:
For example, when I take the inverse laplace of $G_D (s)$, I get $\frac{1}{5}$, and the integral of this is $t/5$, which can be translated to the first curve above.

But when I take the inverse laplace of $G_F (s)$, I get $9te^{-3t}$, and the integral of this is $-9e^{-3t}(3t+1)$, but this is too much of a hassle, how can I relate this to one of the step responses in an easy way?

Maybe I should go back to using the poles and zeros, any tips would be greatly appreciated.


texts:
I want to write a program to filter a single radio station using SDR. Not for any application just a stepping stone to learn. I am going to use C so it is easier for me to see the steps and things going on rather than them being obfuscated with scipy, matlab fmdemod object or something. I have IQ data from my receiver and I know how to change the sampling rate and base-band frequency but once I have the correct data and I know the desired information is there, I need to do frequency demodulation. From reading I have seen that possible would be to do in software:

polar discriminator
digital limiter
phase angle estimate

then I would send this filtered data to some audio api. Does this seem like the right idea? Just looking for input or someone to tell me I am wrong lol and point me in the right direction. Thanks


texts: Extracting signal from data using Otsu thresholding in MATLAB?

texts:
I am very new to signal processing. I am trying to extract a signal from my data using otsu threholding but it's not fully accurate. In this case I know that there is going to be a single signal in the data. It looks something like this --

I get the histogram of the signal. It is a bimodal histogram. I give those counts to the [otsuthresh][2] function in MATLAB. It returns the threshold between 0 and 1. Then I take the histogram of the same signal but normalized with the cdf, which looks something like this,

I binarize the data at the threshold value percentile and I can extract the signal out of it. The extract signal looks something like this (ignore the part where extracted signal drops to zero),

My questions are as follows:-

How can I find the exact beginning of the signal (roughly at x = 200)? I tried increasing/decreasing the number of histograms but it didn't really help. Once I find the exact beginning of the signal, I can consider that to be zero and offset the rest of the signal by that amount.
What should I do in the case where I don't know how many signals are there in my data? Can you suggest me some resources considering the fact that I'm very new to signal processing?

For your info, the data is strain measurements of a beam, including the noise, measured with time.


texts: Different amplitude for different time resolution with FFT

texts:
In the following Matlab code, the amplitude of the spectral domain signal is printed
Temporal_Pulse_Width=1;
t=linspace(0,20,2000);
dt=t(2)-t(1);
Time_Period=round(Temporal_Pulse_Width/dt);
AA=zeros(length(t),1);
AA(Time_Period/2:3*Time_Period/2)=1;
A=20*log10(abs(fftshift(fft(AA/length(t)))));
frequencyAxis = ((0:length(t)-1) -ceil((length(t)-1)/2))/length(t)/(t(2)-t(1));
plot(frequencyAxis,A,'LineWidth',3);
max(A)

If we change the vector t as follows,
t=linspace(0,10,2000);

The result (max(A)) will also change. How can I normalize the fft without having this problem?


texts: U-nets : how exaclty is upsampling performed

texts:
In U-nets I would like a more straight-forward/detailed explanation in how the upsampling part ("right part of the U") is performed.

I read that it can be done by "transposed convolution layers" aka. "deconv layers". I would like a clear example (possibly a bit simplified) on how that is performed?

e.g. initialization, why/how? Detailed sizes/shapes input output etc?

e.g. from here https://medium.com/@keremturgutlu/semantic-segmentation-u-net-part-1-d8d6f6005066 :
"In transposed convolutions we have weights that we learn through back-propagation. In papers I’ve come across all of these upsampling methods for various cases and also in practice you may change your architecture and try all of them to see which works best for your own problem. I personally prefer transposed convolutions since we have more control over it but you may go for bilinear interpolation or nearest neighbor for simplicity as well."

I don't understand the stated "ou may go for bilinear interp etc."

why would we choose a fixed interp method over learned transposed conv2D filter?
How can that make sense?


texts: Finding amplitudes from FFT data of a 2D image

texts:
I'm working on experiments with 2D - wrinkled surfaces.
Basically I am getting a bunch of images like this.
This example is hand-made.

Edit: The test image was made by sampling a cosine function something like
$4*\cos{2\pi(3x+4y)}$

Is there a method to find the amplitudes of the first few dominant sinusoids in the image?

I've seen for simple 1-D examples how to get to the amplitudes. It seems to require that the coefficients of the expansion cancel out just right.

Here is a sample image from experiment:


texts: Having trouble transmitting/receiving a BPSK signal using two USRPs

texts:
I am relatively new to working with SDRs and GNURadio and I am attempting to generate a BPSK signal and transmit it through one USRP and receive and demodulate it on another USRP. I am using GNURadio and I am transmitting from an E310 to an N210 directly connected with a 30 dB attenuator.

Here are the flow graphs for my transmitter/receiver setup:

When I implemented this transmitter/receiver in a single flow graph with no SDRs, it worked fine. However, when I try and use it with the SDRs I get some strange results in the constellation at the receiver end.

From my simulations I expect to see the constellation points in two clusters on the x-axis. Running it through the clock sync and costas loop did not help. I have tried synchronizing the SDRs using an external reference signal but it didn't help. I also used an external synthesizer to send a BPSK signal to my receiver and got a similar constellation.

Any suggestions or ideas would be appreciated.


texts: Deep learning model for phone recognition - issues with dimension the model

texts:
I seem to have some problems understandind how the model described in this paper has been designed

This is what is written about the model dimension..

  ...In these experiments we used one convolution ply, one poolingply
  and two fully connected hidden layers on the top. The fullyconnected
  layers had 1000 units in each. The convolution andpooling parameters
  were: pooling size of 6, shift size of 2,filtersize of 8, 150 feature
  maps for FWS..

So according to ^ does the model consist of

Input

Convolution

Pooling

Input being the 150 feature maps (each with shape (8,3)

Covolution being 1d as kernel size is 8

and pooling is with size 6 and stride 2.

What was expected of output would be a shape of (1,"number of filters), but what i get is (14,"number of filters)

Which I understand why i get, but I don't understand how the paper suggest this can give an output shape of (1,"number of filters")

when using 100 filters I get these outputs from each layer

convolution1d give me (33,100)

pooling (14,100)..

Why i expect the output to be 1 instead of 14

The model is supposed to recognise phones, it takes in a 50 frames (150 including deltas) as input, these being a context frame, meaning that these are used as support to detect one single frame... That usually why context windows are used.


texts: What is an intuitive understanding of the difference between the minimum required bandwidth and actual bandwidth?

texts:
When using for example, polar signaling of full width rectangular pulses, the bandwidth is equal to the pulse rate $Rb$, which is also the bit rate.

$B=R_b$ for full width rectangular pulses.

But also we can describe the minimum required bandwidth for transmission of this data.

$B_{min}=\frac{R_b}{2}$

If I am sending data over a channel using polar signaling of the full width rectangular pulses, surely I cannot achieve a bandwidth of $B_{min}$?

How would I go about converting my data to achieve this theoretical minimum bandwidth? I'm looking for an intuitive understanding of minimum bandwidth versus the actual bandwidth of my data.


texts: Extract a signal that is known to repeat every 4 samples at a known phase i.e. exactly 500 Hz signal in 2000 Hz data

texts:
I have a situation that the measurement device is generating a known signal in the data at 1/4 of the rate. In this case that is a 500 Hz signal when sampling at 2000 Hz.

Now I can extract the 500 Hz signal using a Butterworth Band-Pass at 490–510 Hz which works, but it takes >1/8th of a second to extract the signal. In our setup we know the measurement and the signal are perfect in sync (same clock source etc) so a perfect 500 Hz filter would be ideal.

Is there a better extraction method? Ideally without the latency of the Butterworth or more computationally efficient or more accurate to extract the signal?

Here is an example of the 2 kHz data:

-34561, -32650, -3422, -5064, -33967, -32061, -2807, -4436, -33364, -31447, -2108, -3771, -32863, -31050, -1801, -3516, -32634, -30871, -1711, -3496, -32761, -31115, -2027, -3919, -33237, -31652, -2684, -4606, -33831, -32279, -3394, -5291, -34442, -32775, -3764, -5582, -34667, -32969, -3872, -5589, -34532, -32665, -3446, -5064, -34014, -32082, -2771, -4385, -33355, -31437, -2114, -3771, -32854, -31032, -1796, -3501, -32651, -30819, -1593, -3425, -32671, -31067, -1924, -3859, -33214, -31631, -2639, -4519, -33802, -32225, -3301, -5254, -34406, -32736, -3741, -5594, -34652, -32932, -3886, -5572, -34453, -32592, -3377, -5056, -33953, -32043

Edit:
 The aspect of interest is the amplitude of the 500hz signal, the amplitude can be converted to give a circuit impedance which we wish to efficiently measure. The butterworth takes a long time for a 1st order butterworth to converge to a steady amplitude measurement with fixed circuit impedance. A higher order could viable but this feels like overkill and as there is so many known about the signal so a more efficient approach would feel to be viable, just magnitude is the key. I tried some simple approached but these still had oscillation of amplitude at 50hz albeit the average was good, there may be formal methods I'd like to try.

Update: Here is a plot showing the butterworth vs the given solution shown in Green which does exactly what I had hoped with 4 samples latency! Previously I had looked at peak-signal but taking peak-average does look preferable anyway and what the solution is giving us. Thankyou for your help.


texts: (Rough) camera calibration of manipulated images

texts:
Coming from a photogrammetry background, I am quite familiar with camera calibration, and the common assumptions that are used for crowd-sourced images (i.e. center point in the middle of the image, pixel resolution from the EXIF data etc.).

But I have been wondering a maybe naive question lately: Can I calibrate a manipulated image (e.g. downloaded from twitter or facebook), that has possibly been cropped or maybe even resampled and has no EXIF data? The nice thing is: More than 4 ground control points in the desired object coordinate system can be asummed available.

Is there any papers available on this problem?


texts: Adjusting corner frequency to constrain maximum disturbance in a high-pass filter

texts:
I have a first-order high-pass filter with transfer function:
$$G(f)=\dfrac{G_0 jf}{jf + f_c}$$

where $G_0$ is the gain at high frequencies.

If I input a sine wave with frequency 1 KHz and I want a maximum disturbance of 0.1% in the amplitude, how can I know the maximum value of the corner frequency ($f_c$) that allows this error?

P.S.: Corner frequency is the frequency when the gain is at -2dB.


texts: How to measure algorithmically side lobes in real data

texts:
I want to measure algorithmically the following parameters in a given FFT data (real numbers).

peak side lobe levels
side lobe levels

side lobes nulls

The following image illustrates parameters I need.

How I can compute those values using Matlab or Python?


texts: How to Implement the moving average filter in time domain in MATLAB?

texts:
I have a sinusoidal x(t) = sin(500πt) signal  is corrupted by random noise. The corrupted signal is sampled with sampling frequency fs = 5 kHz and passed through a 5-point moving average filter to reduce the noise. I have written the following MATLAB Code for the same:-
fs=5e3; t=0:1/fs:0.02;
x=sin(2*pi*250*t);           % uncorrupted signal

r=0.1*randn(1,length(x));    % random noise
XR=x+r;                      % signal + noise (filter input)

y=0;
% Implement the moving average filter in time domain
for n=5:length(XR),          % M = 5
    y(n) =                   % missing code using the command 'sum'
end

figure(12);
subplot(211); plot(t,XR,t,x,'r'); legend('Corrupted','Original');
xlabel('time (sec)'); ylabel('Amplitude'); grid;

subplot(212); plot(t,y,'k',t,x,'r'); legend('After filtering','Original');
xlabel('time (sec)'); ylabel('Amplitude');grid;

But, I am not sure how to implement the moving average filter in time domain using sum command to get the output (y[n]) in the for loop as mentioned above in the code. Any suggestions would be highly appreciated! Thanks a lot for the help


texts:
I’m trying to figure out how to go about differentiating between a single-voice speech and a multiple-voice conversation in which voices overlap. In other words – the moment of overlap is what I’m interested in. I want to be able to listen to a vocal conversation between two people, and detect when voices collide/overlap.

I guess that spectrum-wise there’s a major difference between a speech-based conversation with no overlaps compared with a conversation with speech overlaps, so perhaps one direction to look at is detecting sudden spectral changes?


texts: Why is not Fourier Transform Good for Non-linear Processes

texts:
I was reading through slides about Hilbert Huang Transform. In slide 14, which talks about the motivations of a new method instead of Fourier Transform (FT), the author provides those two reasons in addition to other reasons

Physical processes are mostly nonstationary.
Physical processes are mostly nonlinear.

I understand the first reason, where the FT cannot locate the frequencies in time. However, I do not understand why the FT does not work for the nonlinear processes. Based on my understanding, the FT is a method to decompose the signal into different and simple components. Why does condition of non-linearity play a role here?

My first thought for reasoning that assumption, was that we cannot represent a non-linearity in the signal by a linear combinations of sinusoids. However, I think sinusoids themselves will handle the non-linearity.


texts:
I have the following phase noise plot (-100dBc/Hz @ 1MHz offset and -130dBc/Hz noise floor):

The integrated phase jitter for this phase noise is around $4.6^\circ$. According to my understanding the spread of the constellation due to phase noise is equal to the integrated phase jitter.

But I see that the spread is much higher than the jitter. Why is it so? I am adding phase noise as follows:

modData = pskModulator(data);  % Modulate data
tx_out = modData.*exp(1i*phas_noise_signal);  % Add noise


texts: Identify the Type of Image Distortion (On Lena Image)

texts:
I want to understand the distortions in each of the lena images shown in the image below.

I see that the first one is a blur. How can undo the blur. The second one is colored salt and pepper noise. This is different from the gray salt and pepper. How can remove that from the image. And for the last one, I don't know that kind of distortion is in that image.


texts: Initial Rest Condition for LCCDE causal LTI systems

texts:
I am self studying Alan Opennheim's course Signals and Systems. I am a math major and have no background in EE.
I understand that for a linear constant-coefficient difference equation (LCCDE) system to be linear its auxiliary conditions must be 0.
I also understand that for the system to be causal every output $y(t)$ corresponding to input $x(t)$ s.t. $x(t) = 0$ for $t<t_0$ must satisfy $y(t)=0$ for all $t<t_0$.
What i am having trouble understanding is why this imposes the initial condition $y(t_0)=0$ for the response to the unit step function. Couldn't we have set the initial condition to be $y(t_1)=0$ for $t_1>t_0$ and kept the system causal and LTI?
In the 6th lecture in the video series Alan Oppenheim finds the unit step function response to the causal LTI system described by $y'(t) +ay(t) = x(t)$. He imposes the initial condition $y(0)=0$, as i mentioned above i can't understand why this is imposed by the properties of the system. He then goes on to find the unit impulse response, and finds that it is $e^{-at}u(t)$ where $u(t)$ is the unit step. Clearly this function does not satisfy the initial condition which were imposed on the unit step response, I am having trouble understanding why this is valid.
Any help is appreciated, Thanks!


texts: Estimating Velocity from Fuzzy Position Data with Known Uncertainty

texts:
I was wondering if someone could maybe clarify or direct me to the best answer for this question.  I want to estimate velocities from position data.  However, the position data is fuzzy but I have an idea of its uncertainty.  My initial sense is that one is supposed to use some iteration of a Kalman filter, but from my limited understanding of it, you need to know at least one parameter with some certainty, i.e., the position data.  So to summarize, my questions are:

Is it possible to estimate velocity values given fuzzy position data?
If so, what are some available techniques for doing so?


texts:
Consider the filter which equation can be represented by $y[n] = \frac{1}{3}x[n] + \frac{1}{3}x[n-1] + \frac{1}{3}x[n-2]$, in $x[n]$ and $y[n]$ are sequence of input and output of the system respectively.

Mark the correct option about the sequence ${h[n]}$ filter response with respect the unit pulse ($\delta_0[n] = 1$) in $n=0$.

A) h[1] = 2/3

B) h[1] + h[5] = 1/3

C) h[0] = 0

D) h[2] = 3

E) h[2] + h[3] = 1

My attempt:

$Y(z) = \frac{1}{3}\cdot X(z)(1+z^{-1} + z^{-2})$

$1+z^{-1} + z^{-2}$ its a Finite Filter Response right?

$X(z) = 1$

$Y(z) = H(z) \Rightarrow H(z) = \frac{1}{3}(1+z^{-1} +z^{-2})$

I stuck here because and really cant realize to envolve the discrete domain with this expression on Z domain.

Any hints?


texts: Evaluating the spectral density of generated noise through the autocovariance

texts:
I'm working on generating noise signals $X(t)$ (with $t \in \left[0,T\right]$ with step size $\delta t$) with a prescribed power spectral density $S_{XX}(f)$ and I'm figuring out how well I am generating said noise; I need to know if it has the desired spectral density (in case someone notices, I posted about similar issues a while back. I've since then worked on a different topic and I am now returning to this).

However, estimating the spectral density of a signal is a tedious process. It involves a lot of estimation through periodograms and window functions and such, and it is therefore not so easy to quantify numerically if the spectrum not only has the desired trends but actually implements the desired parameters (in practice I will simply use a spectrum analyzer, but that is a different story).

The auto covariance function $C_{XX}(t') = \left<X(t) X(t+t') \right>$ on the other hand is very easy to calculate and subsequently fit. Moreover, the auto covariance function and the spectral density are related by a Fourier transform: $C_{XX}(t') = \int_0^\infty{S_{XX}(f)\cos{2\pi f}}$ (note that I'm using the cosine transform as this is customary in these settings, due to the evenness of the spectral density).

So my question is, how much do I learn about the spectral density by fitting the auto covariance? Say I am working with the Ornstein Uhlenbeck process, which has auto covariance $C_{XX}(t') = \frac{c \tau}{2}e^{-t'/\tau}$ and spectral density $S_{XX}(f) = \frac{2 c \tau^2}{1+(2\pi f \tau)^2}$. If I fit my numerically generated $X(t)$ to $C_{XX}$ to extract $c$ and $\tau$ and am able to do so with high accuracy, do I then know that my spectral density is also as desired?

Intuitively I would say yes, Fourier pairs contain the exact same information. The reason I am not sure if this is the case is due to the fact that I am unsure of how sampling rates and such come into this. Of course one cannot expect a spectral density to be accurate up to the GHz range if one has an $X(t)$ point every millisecond. But is this equally true for the auto covariance? In other words, if my autocovariance is accurate up to some specific lag time (which means I have a high enough sampling rate I suppose), is the spectral density then also accurate up to some specific frequency? And can I quantify this in some way?

What I mean by that last sentence is, can I tell up to (or down to) what frequencies I can trust my spectral density, given my auto covariance fit, $T$ and $\delta t$?


texts: Is there a penis detection demo similar to face detection?

texts:
Tutorial for face detection:

Is there a script / tutorial / demo for penis detection? These guys ran into some issues:

Fairly serious question, future of internet memes is at stake. Breast / nipple detection could be handy too.


texts: For a standard second order transfer function, what is the equivalent time domain significance of $\zeta>0.707$?

texts: In frequency domain $\zeta>0.707$ implies no resonant peak in frequency response. But I am unable to correlate what is the exact time domain effect of this? If  $\zeta=0.707$ acts like an 'important damping value' in the frequency domain, shouldn't its importance be reflected in time domain as well. But in time domain $\zeta=1$ is where system behavior changes.


texts: How to unify the color and brightness measure of two given cameras?

texts:
I have two cameras, and I'm up to perform some matching stuff on videos captured by these two cameras since I'm a CV guy.
The problem is that two cameras have very different performance in terms of color and brightness. The same object in the same illumination appears different in two cameras.
Is there any academic work on measuring something like the "rgb response curve" of a camera? I want to unify the rgb response of these two cameras to guarantee the same object "looks" the same in two cameras.
Note Below is just an intuitive illustration of my problem. Same object appears different in different cameras.


texts: Recovering a signal that was upsampled by non-integer factor

texts:
I have a data sequence $G[n]$ that is upsampled at the transmitter of a communication system by a non-integer factor ($1.5$ to be precise). At the receiver, I am supposed to downsample the data and recover the sequence $G[n]$.

In MATLAB, I am implementing both the transmit and receive signals. When I upsample the signal $G[n]$ by integer factors at the transmitter and then downsample the received signal by the same integer factors, I am able to successfully recover the exact $G[n]$, and continue further processing.

However, when I upsample $G[n]$ by a non-integer factor $1.5$, and then downsample the received data by the same factor $1.5$, I do not get the exact sequence $G[n]$ that was transmitted. In MATLAB, for the transmitted sequence, I simply used commands upsample with factor $3$ and downsample with factor $2$. At the receiver I upsampled by factor $2$ and then downsampled by factor $3$.

Is there a way to recover the exact sequence $G[n]$ at the receiver even if it was upsampled by a non-integer factor at the transmitter?

~ryan


texts: Wiener Filter for Image Noise Reduction (Image Denoising)

texts:
I'm trying to get my head round the operation of the Wiener filter for the purpose of image noise reduction. In my case I'll have used another noise reduction filter first and will then use the result of this as an approximation of the noise characteristics for the Wiener filter.

Regarding information on the Wiener filter, I found the following Matlab code & description useful:

http://www.mathworks.co.uk/help/toolbox/images/f11-12251.html#f11-14272

and a few other good links such as

http://blogs.mathworks.com/steve/2007/11/02/image-deblurring-wiener-filter/

So from a Matlab perspective I can see how to use the inbuilt Matlab function, but I'd like to gain a more fundamental understanding rather than just use the function call, yet at the same time I'd prefer to find something more digestible than the Wikipedia entry on Wiener filtering.

Anyone care to offer a brief explanation on Wiener filtering?


texts: How we do frequency shaping and amplitude shaping in hearing aid

texts:
My project name is Digital Hearing Aid with frequency shaping and amplitude shaping using adaptive noise cancellation DSKC6713 DSP kit .
Almost I have done this noise cancellation using code studio ver 3.1( specifically used for this kit)
Now I get stuck on frequency shaping and amplitude shaping ....after so many searches I still dont get any coding or algo ...can anyone help me..??


texts:
I read in a paper that a sound signal can be delayed by D degrees. How can such a thing be possible? Isn't sound signals expressed only in terms of the time axis?

Thanks in advance


texts: How to model transfer functions using Python (without using higher level packages) and plot the magnitude and phase response?

texts:
I'm trying to model a transfer function in Python and thought I could do it by simply plotting the transfer function at many frequencies. This seemed to work for a 2nd order LPF. See the below figure.

A bit of sample code would be like:
s = np.arange(10 * 2 * np.pi, 1000 * 2 * np.pi, 10) # Creates vector of frequencies in rad/s
lpf = np.zeros(len(s))  # Initializes lowpass filter vector

for k in range(len(s)):
   lpf[k] = lpfilt_func(s[k])  # Calls LPF function with TF = cutoff_freq ** 2 / (s + cutoff_freq) ** 2

To get the frequency response, instead of using s as a real number, I changed s to a complex number using the built-in function complex(0, s) and calculated the magnitude and phase as follows:
for k in range(len(s)):
   lpf[k] = lpfilt_func(complex(0, s[k]))
   gain_lpf[k] = np.absolute(lpfilt_func(complex(0, s[k])))
   phase_lpf[k] = np.angle(lpfilt_func(complex(0, s[k]))) * 180 / np.pi  # Phase in deg

I think this is essentially like setting $s = j\omega$ and then calculating the magnitude as $|LPF| = \sqrt{A^2 + B^2}$ and the phase response as $\angle{LPF} = \angle{(B/A)}$.
Is this a valid method to calculate the frequency response of a transfer function? Is this conversion of s from a real number to a complex number valid?


texts: Non-symmetrical lemniscate curve parameterization

texts:
I'm trying to fit a function to data points. The data generally resembles a butterfly/lemniscate shape, see drawing.

The problem is that the shape in my data can be rotated, skewed and/or non-symmetrical.

I've been looking at Bernoulli's, Devil's curve, Watt's curve, however, these are, as far as I can see, symmetrical.

Does anyone know of a plane curve that is able to represent the example shapes? Preferably in Cartesian coordinates.


texts: confusion understanding FFT butterfly diagram?where are 12 , two input butterflies?

texts:
https://www.alwayslearn.com/DFT%20and%20FFT%20Tutorial/DFTandFFT_FFT_Butterfly_8_Input.html

I was reading above link and then  i came across a confusing sentence, i highlighted it in attached snapshot. It says there are 12 , 2 input butterflies ,while i am only able to see only 4, 2 input butterflies at left most of snap in stage1  but unable to see any 2 input butterfly in stage 2 or stage 3


texts: FSK analysis: only ever find one frequency

texts:
I've trying to write this software as part of my amateur radio and SDR hobby for analyzing various two-level FSK signals.  I'm very new to DSP but have been having a lot of fun so far learning more about this subject.  My long-term idea is to try and automatically detect the mark and space frequencies, as well as the baud rate, in arbitrary signals (maybe this is incredibly naive -- I don't know).

Right now I'm trying to get this software working using a known signal, known mark and space frequencies, and known baud rate.  I've tested it with samples from sources like the Emergency Alert System, Bell 103, and even a pure $2600\textrm{ Hz}$ sine wave just to make sure that my frequency detection algorithm (very basic zero-crossing detection) works correctly.

I get back a $2600\textrm{ Hz}$ frequency when I test against the pure $2600\textrm{ Hz}$ sine wave sample so I'm at least somewhat confident about my ability to detect one particular frequency.  However, I can never extract the mark and space frequencies when I use two-level FSK samples.

What I'm using right now is a sliding window that starts at the beginning of a sample and moves along the sample, changing both the starting position of the window and the overall length of the window, in an effort to try and find a good fit.  But no matter what, I never see the mark or space frequency (actually, sometimes I will get a somewhat near read on one of the two frequencies but never both).

The project is pretty small if you want to have a look for yourself and give me some pointers (no pun intended).  Apologies for the amount of code in the Cli project right now.  I've put the code that I'm pretty confident about into their own classes and behind interfaces but my "work-in-progress" code is currently just living in the Main() method of the command-line interface.

Here is my current "window algorithm" in a nutshell:

Start with window positioned at $0\textrm{ ms}$ into the sample
Start with window length of a fairly minimal value
Count zero crossings across the current window to determine the frequency of the signal in that window
Output various information about the window and the frequency it contains
Increase the window length and repeat until the window length hits a maximum (maybe $1-10\textrm{ ms}$ in $0.25-1.0\textrm{ ms}$ increments -- all easily configured)
Once the maximum window length is reached, reset to the minimum window length, increment the window position a little, and repeat the varying window length frequency samples again
Repeat this until the maximum window position is reached

And here is a screenshot of the software running, where it was able to accurately detect one of the two FSK frequencies but not the other:

EDIT 1: Turns out I had a bad sample.  When I looked at a frequency spectrum plot I only saw one frequency peak, whereas a 2-level FSK spectrum should have two peaks. I found another sample and now my window algorithm is working much better.

EDIT 2: Still using the same algorithm but now with the new sample and a few very minor tweaks to the software am successfully decoding at least the first non-trivial amount of this NWR/NOAA SAME signal.  You can very clearly see the 1010 1011 preamble repeating at the beginning.

EDIT 3: Refactored a lot of the code and cleaned-up the Cli project so that it is actually just a Command-Line Interface.


texts:
I understand that the DFT (Discrete Fourier Transform) process produces the same number of outputs as there are inputs, and this is clear from the basic DFT expression

$$X[k]=\sum_{n=0}^{N-1}x[n]e^{-i2\pi kn/N} \tag{1}$$

However, when computing the FFT (Fast Fourier Transform), the final output is the sum of two N/2 transforms (the even and odd DFTs).

$$X[k]=X_{even}[k]+W_{N}^{k}\cdot X_{odd}[k] \tag{2}$$

Suppose $N=16$. If I try to evaluate Equation (2) for $N=12$, it won't work, since $X_{even}[12]$ doesn't exist (since $X_{even}$ only has $N/2=8$ inputs in the first place).

This is normally solved, as in the case of Fat32's solution to a previous question of mine, by concatenating the outputs of each $N/2$ DFT with themselves to produce an output that is N long, which would result in an $X_{even}[12]$ which does exist. Alternatively, one could add a condition that if $k>N/2$, then substrat $N/2$ from each $X_{even}$ and $X_{odd}$ index.

I was wondering what the mathematical foundations of each of these work around are, since they aren't explicitly included in expression 2? Is it something that is basically common sense? Or does it really mean that only half of the outputs are actually valid?

I also have a hard time interpreting the output of the FFT since the second half of the outputs are mirrored and actually represent negative frequencies. So doesn't that mean that in order to really present an accurate output I must chop that second half off and move it to the negative frequency domain?

There seem to be several ways to present the output (some leave the negative frequencies mirrored and in the positive domain, and some do order them properly). Any advice about the best way to interpret and present the outputs of the FFT process is greatly appreciated!


texts:
Let $R$ and $\Theta$  be zero mean independent Gaussian random variables with variances $\sigma_R^2$ and $\sigma_\Theta^2$ respectively. Hence their covariance matrix would be diagonal give by
$$C_{R,\Theta}=\begin{bmatrix}\sigma_R^2 & 0 \\ 0 & \sigma_\Theta^2\end{bmatrix}$$
Now, let's define transformations given by $X = R\cos (\Theta)$ and $Y = R \sin (\Theta)$.

The question is, what would be $C_{X,Y}$, the covariance matrix of transformed random variable $X$ and $Y$?

Thank you in advance.


texts: DFT Frequency domain analysis and interpolation

texts:
I have a 2 part question, one may be related to why I'm not understanding the other.

A while back, I remember some professor saying that for the $N$ point DFT frequency domain, the values for $k$ around $N/2$ (halfway point for the spectral domain - those around $\pi$) are considered high frequencies, and the end points $(0, N-1)$ (those around 0, $2\pi$) are considered low frequencies.  I am having trouble seeing this mathematically given the definition for computing the coefficients of the DFT.

$$X[k]=\sum_{n=0}^{N-1}\left(x_ne^{\frac{-i2\pi{}kn}{N}}\right)$$

To me, it seems that the frequencies grow as $k$ increases right up until $k=N-1$, where the coefficients repeat due to periodicity.  Perhaps I am mis-remembering and this may be an incorrect assumption?

This leads on to my second question.  I have seen a few examples floating around online, where to interpolate in the time domain, one may pad in the frequency domain (running the IFFT after of course).  The examples I have seen seem to decide to pad the coefficients around the middle, i.e. for $N=4$ and a desired sampling of $2N$, [1 2 3 4] become [1 2 0 0 0 0 3 4].  Is there a reason why the middle is targeted and not the end like [1 2 3 4 0 0 0 0]?

If my first question is correct, then this simply pads the higher frequencies as 0.  I am aware that these zeroes add nothing to the frequency information as $X[k]=0$ cancels any effect during the inverse DFT.  But wouldn't moving any of the coefficients change some frequency information, i.e. for [1 2 3 4] becoming [1 2 0 0 0 0 3 4]

$N=4, k=2, X[k]=3$, becomes $2N=8, k=6, X[k]=3$

$$e^{\frac{+i2\pi{}2n}{N}}\neq{}e^{\frac{+i2\pi{}6n}{2N}}$$


texts: How do I find the convolution of the following signals?

texts:
$h(t)$ is the inputted to convolve the $x(t)$

The signal written in unit step are:

$x(t) = (5-t)u(t-3) - (5-t)u(t-5)$

$h(t) = 2u(t-1) - 2u(t-3)$

So to convolve I first change the function to:

$h(t-\tau$$) = 2u(t-\tau$$-1) - 2u(t-\tau$$-3)$

The problem however is the signal I get after integrating the product of $x(\tau$$)$ and $h(t-\tau$$)$ doesn't even closely resemble the signal I get via $conv$ function on MATLAB.

I'm trying to do this through integration in the time domain.

If anyone can work it out and explain why, would help me understand this.

EDIT: $y(t)$ is the convoluted signal. The graph below shows one using the $conv$ function the other I manually calculated.

yt2 = (t-4).*(14-t).*u4 - (t-4).*(14-t).*u6 + (t+4).*(15-t).*u6 - (t+4).*(15-t).*u7;


texts:
I have $x(n)$ = {$4,-1,-3,1$} and $h(n)$ = {$2,1,3,5$}.

I would like to know how can I check whether

$h(n)$ is stable filter or not and
$h(n)$ is causal filter or not?


texts: What Is the Difference Between PCA and Karhunen Loeve (KL Transform)?

texts: I have been reading about Karhunen-Loeve or also known as KL transform and I see that when it is used to reduce dimension the procedure is identical to PCA, that is, for both methods the covariance matrix of the data is constructed and then the eigenvectors are calculated . I would like to know if there is any difference in this aspect that is being overlooked?


texts:
I'm trying to soft decode DQPSK. According to this paper the bits are:

$b_1 = \mathrm{Re}\{s_n s^*_{n-1}\}, b_2 = \mathrm{Im}\{s_n s^*_{n-1}\}$

If I put 0.707 as a 100% probability of 0:

$\begin{eqnarray}
s_n=(0.707 + 0.707i) \nonumber \\
s_{n-1}=(0.707 + 0.707i) \nonumber
\end{eqnarray}$

This gives:

$\begin{eqnarray}
b_1=0.707 * 0.707 - (0.707 * -0.707)=0.999 \nonumber \\
b_2=0.707 * -0.707 + 0.707 * 0.707 = 0.0 \nonumber
\end{eqnarray}$

And the result doesn't make sense:

The real part (0.999) is fine. Original data corresponds to 100% probability of 0 and real part now from -0.999 to 0.999.
The imag part (0.000) is weird. It can span from -0.999 to 0.999. Being 0.000 means bit 0 has 50% probability. Which is wrong, because original $s_n$ and $s_{n-1}$ correspond to 100% probability.

So it looks like the result constellation is no longer the same as original DQPSK: (0,0) bits now mapped to (1,0) point.

Questions are:

Am I right in assumptions?
How can I get soft symbols from the result constellation?


texts:
If I have transfer function coefficients,
I can analyze the transfer function in the s-plane and/or the z-plane.

If I wanted to demonstrate that the z-plane and s-plane responses are equivalent:

Would it be possible to use the frequency warping relation
(provided in the otherwise unused variable w.z in in the code)
to depict the z-plane plots in the s-plane or vice versa
in a manner in which the two would completely match?

Matlab code provided:

[Actual filter coefficients not provided.]

%% Input

A.s.cf = [ 0 0 0 0 1 ];
B.s.cf = [ 1 2 3 4 5 ];

T.z   = 2.5e-05;
f_z   = 1/T.z;
f_nyq = f_z / 2;

%% Processing: Matlab SP toolbox
%{
H.s.tf = tf ( A.s.cf, B.s.cf       );
H.z.tf = c2d( H.s.tf, T.z,   'zoh' );
%}

%% Processing: From scratch:
f.s = logspace(-2, +6, 5e3).';
w.s = 2*pi*f.s;
s.s = 1j * w.s;

z.z = exp( 1j * w.s * T.z);

w.z = 2/T.z * atan( w.s * T.z/2 );

switch 0
  case 0; s2z = 1/T.z * log(z.z);
  case 1; s2z = 2/T.z * (z.z-1)./(z.z+1);
end

A.s.jw = s.s.^4 * A.s.cf(1) ...
       + s.s.^3 * A.s.cf(2) ...
       + s.s.^2 * A.s.cf(3) ...
       + s.s.^1 * A.s.cf(4) ...
       + s.s.^0 * A.s.cf(5) ; %

B.s.jw = s.s.^4 * B.s.cf(1) ...
       + s.s.^3 * B.s.cf(2) ...
       + s.s.^2 * B.s.cf(3) ...
       + s.s.^1 * B.s.cf(4) ...
       + s.s.^0 * B.s.cf(5) ; %

A.z.jw = s2z.^4 * A.s.cf(1) ...
       + s2z.^3 * A.s.cf(2) ...
       + s2z.^2 * A.s.cf(3) ...
       + s2z.^1 * A.s.cf(4) ...
       + s2z.^0 * A.s.cf(5) ; %

B.z.jw = s2z.^4 * B.s.cf(1) ...
       + s2z.^3 * B.s.cf(2) ...
       + s2z.^2 * B.s.cf(3) ...
       + s2z.^1 * B.s.cf(4) ...
       + s2z.^0 * B.s.cf(5) ; %

H.s.jw = A.s.jw ./ B.s.jw;
H.z.jw = A.z.jw ./ B.z.jw;

If I plot the two,
the response will be different in either plane,
(as expected).

%% [Plot]: Setup
p = 0;

%% [Plot]: S plane
p = p+1;
figure(p)

% S plane: Magnitude
subplot(2,1,1)
semilogx(w.s, 20 * log10( abs(H.s.jw) ), '.')

title('Bode [From Scratch]')
ylabel('Mag [dB]')
xlabel('w.s [rad/s]')

axis( [ 1e+3 1e+6 -150 +0] )
grid minor

% S Plane: Phase
subplot(2,1,2)
semilogx(w.s, angle( H.s.jw ) * 180/pi, '-')
hold on
semilogx( [ 2*pi*1e-2 2*pi*1e+6 ], [ -180 -180 ], 'k')
semilogx( [ 2*pi*1e-2 2*pi*1e+6 ], [ +180 +180 ], 'k')
hold off

ylabel('Ang [deg]')
xlabel('w.s [rad/s]')

axis( [ 1e+3 1e+6 -200 +200] )
grid minor

%% [Plot]: Z Plane
p = p+1;
figure(p)

% Z Plane: Magnitude
subplot(2,1,1)
semilogx(w.z, 20 * log10( abs(H.z.jw) ), '.')
hold on
semilogx([2*pi*f_nyq 2*pi*f_nyq], [-50 +50], 'k-')
hold off

title('Bode [From Scratch]')
ylabel('Mag [dB]')
xlabel('w.z [rad/s]')

axis( [1e+3 1e+6 -50 0] )
grid minor

% Z Plane: Phase
subplot(2,1,2)
semilogx(w.z, angle( H.z.jw ) * 180/pi, '-')
hold on
semilogx( [ 2*pi*1e-2  2*pi*1e+6 ], [ -180 -180 ], 'k')
semilogx( [ 2*pi*1e-2  2*pi*1e+6 ], [ +180 +180 ], 'k')
semilogx( [ 2*pi*f_nyq 2*pi*f_nyq], [ -180 +180 ], 'k-')
hold off

ylabel('Ang [deg]')
xlabel('w.z [rad/s]')

axis( [1e+3 1e+6 -200 +200] )
grid minor

%% [Plot]: S Plane: Bode

%{
p = p+1;
figure(p)

bode(H.s.tf)
grid minor
%}

%% [Plot]: Z Plane: Bode

%{
p = p+1;
figure(p)

bode(H.z.tf)
grid minor
%}


texts: Smart way to sample in "time domain" for a known "frequeny domain"

texts:
I have an experiment in which every point in the "time domain" is very expensive to take. Good news is I know the center frequency and the bandwidth of the signal.

How can I sample (which times should I probe) in the time domain to get the most information about the frequency range of interest with the least amount of points?

Is there some known method I can look into? If anyone has any resources/tutorial on this, I would be very grateful if you could point me towards them.

EDIT 1:

  Please provide more details on your experiment, the actual expected frequency bandwidth and why the time domain samples are "expensive". Depending on the context the best answer will be very different. – Dan Boschen

In the actual experiment the "time domain" is a voltage that we use to manipulate electron spins. The Fourier transform of the spin response to the voltage gives us information about the spins surroundings (its just a blob in this case, but could have sharper features). I am attaching some real data (that encompasses the entire "frequency" range) so you get a sense of the center "frequency" and bandwidth.

The sampling is "expensive" because each pulse voltage needs averaging, taking a few minutes, and we would like to make this as fast as possible. Intuitively, its clear that we are oversampling as there as many points with a value of 0, hence the question.


texts: normalizing STFT output and calculating energy of signal

texts:
I have two questions regarding power spectra and energy of a signal:

For power spectra calculation I used STFT, and I saw that the value of STFT should be normalized before applying logarithmic function. The normalization is done by multiplying them by $2/\mbox{window length}$. What is the reason?
I also want to calculate the energy of the signal according to Parseval theorem and using STFT. Should I use the average of STFTs or the sum of them?


texts: FFT calculation suitable, Wav File, Java

texts:
If I have a wav file with this characteristics:

16bits (2 Bytes per Sample)
32000 Sampling Rate
60 Seconds
2 Channels

then I have:

Samples Per Channel: 32000*60 = 1920000
Samples Per File: 32000*60*2 = 3840000
Bytes Per File: 32000*60*2*2 = 7680000

It's a huge file!

If I change my file Characteristics to:

8bits (1 Byte per Sample)
8000 Sampling Rate
60 Seconds
1 Channel

then I have:

Samples Per Channel: 8000*60 = 480000
Samples Per File: 8000*60*1 = 480000
Bytes Per File: 8000*60*1*1 = 480000

If my FFT 2^12 = 4096 or 2^16 = 65536

Do I need to calculate the FFT 480000 times?

Does each how many samples I need to calculate the FFT?  (480000 - 65536 = 414464?)

what is the suitable calculus?

How to do to show the spectrum?


texts: Guitar Tuner or Guitar effects Algorithms?

texts:
I need to explore the algorithm of the digital guitar tuner or the digital guitar effects algorithms.
Is there any Book or articles that you recommend of?

All I know is that relate to DSP.

thanks


texts: Detecting peaks in a noisy 2D correlation plot

texts:
I'm doing pattern matching by means of normalized grey-scale correlation. My input image is noisy and so the peaks in the correlation plot cannot be found just by thresholding.

The original correlation plot displayed as a surface, the peaks are visible, but are not the only local maxima in the picture, I need to find only the ones that are cone-shaped. Note that I know how many I'm looking for.

I noticed that if I calculated the local gradient of the correlation plot, the peaks appear clearly on visual inspection, especially on the gradient orientation plot.

Gradient orientation:

I'm looking for a way to automatically detect the exact position of the peaks (2 in this case, but it is sometimes more) with a one pixel accuracy.

What I guess would help me is a way to identify the center pixel towards which the gradient direction points from all around. The center of the radial gradient, otherwise said.

I have tried to generate an ideal radial gradient such as this one:

and fit it to the image by moving correlation but with no success.

Does someone have an idea? I'm certainly not the first one needing this kind of method, but I couldn't find anything in the literature, am I missing something?

Edit: I found someone who had a similar question here but couldn't find a satisfactory answer.


texts: Linear Predictive Coding for general signals

texts:
I have a signal that's monotonic and roughly linear and have been looking at using Linear Predictive Coding to encode information and compress my signal. I guess my first general question is if this is even possible for general signals? I also see LPC mentioned in the context of speech so I'm wondering if LPC is particularly effective for one type of signal and would be really bad for another type of signal. Here's how I understand the process so far:

We have an input signal $x(n)$ that we look to model as $$x(n) = \sum_{i=1}^{p} a_i x(n-i) + e(n)$$

where the $a_i$ are our prediction coefficients and $e(n)$ is the error. We can find the error $e(n)$ by filtering the inverse transfer function $\frac{1}{H(z)}$ with our signal $x(n)$. Thus, if we can closely model the error with a small number of a parameters, our signal has been reduced to our prediction coefficients and these small parameters.

I see why LPC is often used for speech compression - because the error $e(n)$ is roughly white and so we can roughly model it by its variance. However, my intuition tells me that it should be much easier to compress something like perfectly linear signal than a highly nonstationary signal like a vocal sample. So I guess my questions are (provided I'm understanding LPC correctly):

Do people use LPC for things other than speech, and if so, how do they go about modeling $e(n)$ to achieve compression?
If my intuition is correct about it being easier to compress a linear signal using LPC, how would one do it in this case? I'm confused because since a perfectly linear signal is exactly predictable, it should produce $0$ for all $e(n)$. However, filtering a linear signal with prediction coefficients $[2,-1]$ produces a non-zero $e(n)$.


texts: Detection of wobbling sound with rising notes

texts:
In an audio file I would like to find calls of a specific bird specie, one that makes wobbling sound with rising notes. Can you recommend what would be the best audio feature to calculate in time or frequency domain to capture its essence?

Among others I tried MFCC, but the results were not promising (or perhaps my implementation poor) .

[Edit]

In this post I explained briefly the problem and included exemplary audio file with spectrogram: Example.

Here is another example. All marked areas contain kiwi (famous flightless bird from New Zealand) calls I am trying to identify. Mark that the original audio is very noisy; presented spectrogram is a result of applying high-pass filter and spectral subtraction on areas identified as noise-only.


texts: Are there any applications using absolute Savitzky-Golay filter to smooth data and preserve non-negative properties?

texts: I want to smooth peaks for non-negative matrix factorization (NMF). Since Savitzky-Golay has polynomial smoothing, I think it is good for narrow peak shape in my signal profile. But NMF needs to be non-negative, I wonder I can simply take absolute on Savtizky-Golay filter (or replace negative value with eps (1e-6))? Are there any published applications using only positive range of Savitzky-Golay to smooth data? Or are there other appropriate filter for non-negative data?


texts: For 2D signals can it be said that the frequency response is the same as the Fourier transform?

texts:
Say I have a signal: [1, 1 , 1; 1 100 1; 1 1 1];,

I wanted its frequency response, then can I take its FFT?
Why does MATLAB have two functions to do the same thing, fft2 and freqz2?


texts: Adaptive algorithm for active noise control

texts:
I'm looking at building a cost-effective active noise control system.

Basically, the system comprises two microphones - one reference and one error. A DSP unit samples and phase flips the reference microphone and plays the flipped signal against the input for destructive interference.

The error microphone then feeds-forward the resulting signal, which should be 0 (but is never quite). How might one develop a proportional control algorithm for the DSP unit to adapt the reverse signal from the error?


texts: How to build a realtime simple audio spectogram?

texts:
I'm trying to build a spectogram view of history time of 4 seconds and I'm stuck at the part where you load/draw FFT vectors into the 2D bitmap. I just don't know how to do that, but i got the FFT in realtime.

I'm doing it in swift(iOS)

Thank you!


texts: Two-sided impulse response and its inverse in real life

texts:
Dear StackExchange gurus,

I would like to pose this question from a practical standpoint rather than a theoretical one, although perhaps some context in the latter might come useful.

Real-life devices (say Digital-to-Analog converters for example) have a roll off in their frequency response due to analog front-ends. To find said response (and their inverse) one can use spectral methods (T(f) = out(f)/in(f)), or time-domain methods that will extract the impulse response (Wiener-Hopf equations, LMS, etc).

Now, in time-domain methods we can chose our solution to have the dominant coefficient at the center tap of the filter, or we can choose it to be the first tap (i.e. single sided and with zero delay). Both of these give similar results in some occasions, but there are times when the single sided inverse does not invert the system as well as the two-sided does.

Can someone give me some insight on why this happens? Which of the two approaches is better for equalization (regardless of implementation complexity)?


texts: Adaptive filtering: Optimum filter length and delay

texts:
I'm trying to find the optimum filter length for an Adaptive Filtering, using RLS Algorithm.

I'm using this design:

So the "error" signal is the signal without noise (and that's the signal that I want).

If I have $e(n) = d(n)-y(n)$ but $d(n)$ is my desired signal, I need that $e(n) \rightarrow 0$ so I find the optimum filter length (and the delay) using the MSE criterion, but now I have the signal that I want as the error, so I don't know how to find the optimum filter length, because I have NO idea what MSE do I have to get at the output!

Can anyone tell me what should I do?

Thanks!


texts:
I need help with iris localization. I have applied the imfincircle function to locate it, but it only locates the inner boundary of the pupil and does not find outer boundary.

How do I locate the outer boundary?


texts:
I have an imaging problem.  I essentially have a metallic plate with a number of holes through it.  I am looking at the number of holes with a thermal camera.  While looking at the plate I have two phases of this analysis, 1.a source of hot fluid is flowing through the holes and behind the plate 2.a source of cold fluid flowing behind the plate and through the holes.

During each phase I take a number of images N.  So in effect I have a set of random variables for each pixel in the image.  $H[x,y] = (x_1, x_2 .... x_n)$ and similarly $C[x,y]=(x_1, x_2 .... x_n)$ where H is the hot frames, C is the cold frames , x-y are the focal plane position of that particular pixel and 1,2,3 are the discrete time of the reading x.

As it stands I'm just taking the mean value for H and C then subtracting them to generate an image ($I=\bar{H} - \bar{C}$ where the bar denotes the expected value) to gain an image that has significantly decreased noise and removes much of the background effects.

This is really just to decrease the noise but I think there is more information there that I'm not utilizing considering this is not a wide sense stationary process, but is evolving with time as the temperatures reach equilibrium for each phase.

I'm looking for a way to utilize this data and for any recommendations on processing this time evolution to form an image.  Any suggestions?

I have tried Finite Difference Operator to approximate first and second derivatives, but this and various other algorithms have failed me.


texts:
I have confused question regarding the possibility of combination between QAM modulation and PAM in order to improve the spectral efficiency.  Is that possible?

As you see in the below Figure, it's the QAM modulation, where the four status of 1+j, 1+j, -1+j and -1-j are existed. where J represents the imaginary part.

My question, what's about using the amplitude with that? For example, in case of combination with PAM modulation, we can have eight status instead of four, which are:
1+j, 1+j, -1+j and -1-j  at the same time the four added status will be 1/2(1+j), 1/2(1+j), 1/2(-1+j) and 1/2(-1-j), where the amplitude is used to represent additional bits to transmit too.

Is that reasonable and logical to be used ? Is there any articles or books which explain that ?


texts: Calculating FIR filter coefficients with the windowing method and Convolution algorithm

texts: Please advice how to calculate FIR filter coefficients with the windowing method using the Rectangular Window function and convolution method? What utility or FIR design application can perform this task?


texts:
My receiver contains: AGC,timing recovery,carrier recovey blocks.

I am using gardner timing recovery with loop filter and farrow parabolic filter.
The input to the farrow filter is the fractional delay. I am assuming the fractional delay should be constant after preamble time may be only small variations. The fractional delay does not settle in the preamble time.
Between two data burst only pure carrier is sent so there any no transition for
timing recovery to work. So immediately after the data burst there is lot of variations in the frac. delay. which causes timing errors. any suggestions???


texts: What should the time-shift be when implementing a continuous wavelet transform on a computer?

texts:
I'm currently researching implementation methods of the Continuous Wavelet Transform(CWT). On paper, the CWT produces infinitely many outputs on a finite signal since the scaling and shifting parameters are continuous.
When discretizing the CWT for a computer application, it is common to limit the scaling based on the bandwidth of the signal and the number of desired voices per octave.
However, I'm confused as to how the time-shifting should be limited. Since we're operating on a discrete input signal, we could technically just time-shift each wavelet by 1 sample to obtain the highest resolution. Yet, this seems wasteful for low-frequency convolutions.
It seems like the time-shift should be a function of scale, such that the time-shift is smaller for high-frequencies and larger for low-frequencies.
How is the CWT time shift typically determined when implementing the algorithm on a computer?


texts: $\tt fft$ vs. $\tt fftshift$ in MATLAB: Please explain

texts:
I have been looking around on the internet for a really clear explanation for why fft and fftshift give different results.

I have a very long time series for which I have applied a Hanning window. The sample rate is $8\textrm{ Hz}$ and the recording is $\approx 9$ hours long ($262144$ point time series). The time series is assumed to be stationary.

When I do an fft in MATLAB, I get the following result:

When I do an fftshift in MATLAB, I get the following result:

Four questions:

How do I decide which one to use?
Which is "correct"
Why are they giving such different answers?
Why is the fftshift only real-valued while the fft is complex?

Pseudo-code:

t=load('t.mat'); %load time series
N = length(t);
dt = 0.125; Fs = 1/dt; df = Fs/N;
f = -Fs/2:df:Fs/2-df;

t_hann = hanning(t);
F_fft = fft(t);
F_fftshift = fftshift(t);
figure(1); plot(f,F_fft);
figure(2); plot(f,F_fftshift);

Any help is greatly appreciated.

Cheers


texts: Follow up question regarding: "Complex sampling can break Nyquist?"

texts:
I'm having some trouble understanding the sample rate limitations when considering a complex baseband signal.

I understand (based on the linked SE questions below), the either (1) physically sampling a signal with frequencies between $0$ and $B$ using 2 ADCs at sample rate $B$, one 90 degrees out of phase, or (2) sampling with one at $2B$, can satisfy the Nyquist criteria.

I'm working with a system where we have a real digital signal of bandwidth $B$, sampled at $2B$. It has then been converted to complex baseband representation by calculating the I/Q components from the existing samples. After that, the center frequency (also carrier in this example) is moved to zero and the signal now occupies $-B/2$ to $B/2$.

I suspect that even though the complex representation is twice the data, I still cannot down sample without breaking the Nyquist criteria because my I/Q components are dependent. For example, each real sample in the original signal is used to compute an I and Q component for that time.

Is this actually the case? Or can I somehow down sample by 2 and still fully represent the signal in complex form with a sample rate of $B$? To clarify, I believe that I could if the I/Q components were physically sampled independently.

I've done a few simple tests of this in python using a chirped signal, I could share some of the code if it will help. So far, if down sampling does violate Nyquist, I'm still able to reconstruct the original signal nearly as well as in the case with no down-sampling.

This question: Complex sampling can break Nyquist?

This answer on another question: Answer


texts:
I'm trying to figure out what a good DSP unit would be to try and learn how to do some embedded programming on? I would like something relatively inexpensive but that still has a decent amount of "kick" to it performance wise. I was looking at the TI ez-DSP but I don't know what to think of it.

I would like it to be programmable with C or C/C++

I would Prefer 32 bit or higher.(if possible)

Fixed point math.

Something in the Ghz clock speed.(once more if possible)

Support for hardware modulo would be a plus.

I know there are many other considerations to take into account and some of the things I have specified may not be possible to get all in one unit especially at a low cost.

I will be using it for Audio Processing and Audio Signal Generation So if a unit exists that is either compatible with ADC and DAC units or has them built in would be good.

I am open for any recommendations or any advice at all.

Thanks In advance, Alex


texts: How to find the value of regularization parameter/s when having one or multiple priors

texts:
I am working on a image reconstruction framework, specifically interferometry. In simple words we have data that is pertubated by noise, say

$$
V^o_i = V_i + \sigma_i
$$

Where $V^o$ is a complex value and usually called the measured/observed visibility which is corrupted by Gaussian noise, $V$ is the Fourier Transform of an image $I$ and $\sigma$ is the noise of each measurement.

We can solve the problem solving a non-linear equation such as:

$$
 \Phi(I)= \arg \min_{I} ||\frac{V^o - V^m(I)}{\sigma}||^2
$$

where we call $V^m$ model visibilities. I don't want to extend much more here, but since we have an infinity number of solutions and given that a chi-squared will return a noisy image we want an image that represents the highest possible fidelity to the data.

In that sense we want to add prior to regularize the solution. Here, we can think of Maximum Entropy, L1 Norm, Total Squared Variation, etc. Consequently, the equation turns into:

$$
 \Phi(I, \lambda)= \arg \min_{I} ||\frac{V^o - V^m(I)}{\sigma}||^2 + \lambda P(I)
$$

where $\lambda$ is the regularization parameter and P is the prior function.

I've been reading papers about how to find the $\lambda$ parameter and most of them use an estimation calculating the variance of the noise. But in that case then $\lambda = \text{Var}(\sigma)$? This is also called the discrepancy principle (See Regularization in Image Restoration and Reconstruction by W. Clem Karl) In that case, $\lambda$ would be the same for each prior right? But it has been seen doing the L-curve plots for same dataset and different prior that the "best" $\lambda$s for each one of the cases are different.

Moreover, since I am constructing a framework, the idea is to plug different priors, and then the equation will turn to:

$$
 \Phi(I, \lambda)= \arg \min_{I} ||\frac{V^o - V^m(I)}{\sigma}||^2 + \lambda_1 P_1(I) + \lambda_2 P_2(I)+..+\lambda_n P_n(I)
$$

On that case the calculation of the different $\lambda$s is not that trivial as before, I think....

I would appreciate if you can give me a slight hand with both of these cases.

Cheers


texts: I want to do some DSP and Machine learning experiments on Electrical and Acoustic signals, but, due to some language setbacks, I didn't know how to call the type of my signal, what I use to google now, is "Time series". Is there a more precise "key word" "terms" that may put me in the right direction?


texts:
Consider the analog signal,
$$x(t) = 2\cos(3000\pi t)+3\sin(4000\pi t)+7\cos(6000\pi t)$$
If the sampling rate is 8000 samples per second and quantized at 8 bits, find:

Discrete values at any two points
Quantization error at those points

I have been trying to solve this with no success so far. How can I solve this?
So far what I have done is convert $x(t)$ which is continuous time signal to $x[n]$ which is discrete time signal.

$$x[n]= 2\cos\left(\frac{3\pi}{8}n\right)+3\sin\left(\frac{3\pi} {2}n\right)+7\cos\left(\frac{3\pi}{4} n\right)$$ where $n$ is integer.
By any two points what I mean is for any two values of $n$. For example, suppose I calculate $x[n]$ for $n=1$ and $n=2$. How do I find quantization error at those points?


texts: Detection of similar segments of two signal

texts:
I have two signals, the first is a reference signal from which I selected a segment.

I'm trying to find which segment of the other signal is the most similar to this selected segment.

I'd be glad to have some ideas.


texts: IQ Sampling and Oversampling Differences? *Confused*

texts:
Two different thought process that produce two different alias band calculations .... What concepts are these?

If I have two 1 GHz ADC's and they sample with 90 degree offset (A comment said this is Over Sampling)

Each ADC has an Fs/2 of 0.5 GHz
The combined system sample streams have an Fs/2 of 1 GHz because the second ADC captures additional sample points between the samples of the first - doubling the highest frequency possible by one alone.
The system can be thought of as having Fs of 2 GHz, whereas a single equipment Fs is 1 GHz.
After an FFT, the entire spectrum repeats every 2 GHz, the 'system' Fs not the equipment Fs of 1 GHz.
After an FFT, The positive frequencies are 0 to 1 GHz and the negative -1 to 0 GHz.

I am still a bit confused by the thought that maybe everything above is false because maybe IQ sampling means ...

The FFT is being fed a complex number via IQ sample streams.
The FFT fed a complex number does not have negative frequencies - it can span the entire single equipment Fs of 1 GHz.
So Fs/2 is still 0.5 and Fs is still 1, each one is not doubled like the first bullet points suggest because the complex notion of FFT allows frequencies up to the single equipment Fs of 1 GHz when complex sampled.
So also, Fs is 1 GHz and not 2 GHz (3rd bullet point from first set of points) and so the spectrum repeats every 1 GHz and not 2 GHz.
Also in this case, there is no 'spectral inversion' from negative frequencies.

It makes a big difference to where I calculated the alias band and thus my analog filter whether the first set of statements is true or second set is true

Please can you offer advice. Still require assistance if you can 02/03/2020


texts: How use FIR filter for simulating an integral implemented using the trapezoidal integration?

texts:
Suppose I have this equation
$$
\phi = \frac{60}{T^5} \int_0^T \left( T^2 - 6T \tau + 6\tau^2\right) y(t-\tau) d\tau - \frac{30 \alpha}{T^5} \int(T-\tau)^2\tau^2 u(t-\tau) d\tau
$$
and I want to use digital Filter for calculating this integral.

How I can do this? What is the relationship between low pass filter and digital filter in this case?


texts: Sine wave amplitude modulation - Analysis using Audacity

texts:
I am trying amplitude modulation in C++ on a sine wave which is a 1000 Hz sine wave sampled at 48000 Hz. The original sine wave is like this (In Audacity) -

When I basically double the amplitude for each sample, it works except for the flipping up of the extremities -

When I multiple each sample with 4, this is what I get -

Why is this the case, when the equations that I am using in the algorithm don't have anything of this sort? What if I want to remove this 'error'?

P. S. When I import the .raw data, I am using Signed 32-bit PCM encoding, Big-endian Byte order and single channel, with zero offset and 48kHz sample rate.


texts:
I have obtained a set of three poles for a third order Chebyshev filter as shown below:

$$p_k=\rm -0.2471+0.9660j,\quad -0.2471-0.9660j, \quad -0.4942.$$

However I am unsure of how to actually implement these into the Chebyshev prototype filter equation?

Any help would be greatly appreciated!


texts: How to estimate the period of a pulse train?

texts:
I am trying to find the period of a pulse train like signal.

I have reviewed the following questions but I did not find a clear-cut answer to the problem, 12892, 16502, 6260.

I do not know the period nor the pulse width. I work with windows that will contain between 2 and 16 pulses. Pulses may vary in length but will usually be shorter than half the period. The signal may be noisy, although I would be happy with a solution for the case of the (relatively) clean signal. I can work out a minimal pulse width if it helps.

Here are plots of the type of signal I am working with.

I am primarily interested in the period of the signal. The analysis needs not to be online.

For the first example, I'm looking to find a result of about 136.5.

My first attempt was to look at the signal in the frequency domain. But a pulse train in time domain is also a pulse train in frequency domain. I got this plot for the clean signal.

Next I tried autocorrelation. I got the following:

It is much better. I need to get the value automatically though, and I'm not too sure how to find that second local maxima without introducing too much ad-hoc thresholds that may fail when there are more or less peaks in the window. Additionally, I'm looking for a floating point result.
Also, the autocorrelation for the noisy signal was much less clear cut, even after smoothing the input.

There is probably another approach by detecting the rising edges in time domain?

Question:

Is there a standard way to estimate the period of a pulse train?
If autocorrelation is the right approach, how to automatically find the floating point coordinate of the peak from the autocorrelation result without too many assumptions?


texts: Relationships between $f(t)$ signal and $f(t) \times \sin(\omega t)$ signal

texts: Let's assume that there is a signal $f(t)$, $0\le t \le T$. Then, is there any special meaning of the signal  $g(t) =f(t) \times \sin(\omega t)$? Like AM transmission. What is shown/heard when we play $g(t)$ in an audio player? Is there/what are meaningful differences  between $g(t)$ and $f(t)$?


texts: generating spectrograms in python with less noise

texts:
I'm trying to do something similar to Shazam (paper here) in Python. Basically the idea is to use FFT with sliding window to transform both WAV files and recorded audio to the same spectral representations for recognition purposes.

The only problem is the spectrograms I extract from certain songs just seem to generate lots of useless peaks that become false positives when heard in other songs. They look a lot like noise. Examples below.

Spectrogram of SONG_1 [channel 1]:

Spectrogram of SONG_2 (dominates and generates WAY more peaks) [channel 1]:

Spectrogram of SONG_3 [channel 1]:

Spectrogram of SONG_1, recorded through my laptop mic (looks like noise) [channel 1]:

Here's the code I use to generate the spectrograms:

import pylab
from scipy.io import wavfile

fs, frames = wavfile.read("song1.wav")

channels = [
    np.array(frames[:, 0]),
    np.array(frames[:, 1])
]

# generate specgram
Pxx, freqs, t, plot = pylab.specgram(
    channels[0],
    NFFT=4096,
    Fs=44100,
    detrend=pylab.detrend_none,
    window=pylab.window_hanning,
    noverlap=int(4096 * 0.5))

I've implemented the hashing, storage, recognition, etc as well but they are rather useless unless the extracted peaks are distinct to the song heard.

I should mention my peak extraction algorithm is not just a threshold - I find points in Pxx that are the highest value in a neighborhood of around 15-20 adjacent points.

What am I doing wrong? How can I get more distinctive spectrograms out of audio sample channels?


texts:
In Bishop book, page 4, section 1.1, there's a notation I don't seem to understand what's meant by it. The whole paragraph, with which the section begins, is:

  We begin by introducing a simple regression problem, which we shall use as a running
  example throughout this chapter to motivate a number of key concepts. Suppose
  we observe a real-valued input variable x and we wish to use this observation to
  predict the value of a real-valued target variable t. For the present purposes, it is instructive
  to consider an artificial example using synthetically generated data because
  we then know the precise process that generated the data for comparison against any
  learned model. The data for this example is generated from the function sin(2πx)
  with random noise included in the target values, as described in detail in Appendix A.
  Now suppose that we are given a training set comprising N observations of x,
  written $x \equiv (x_1, . . . , x_N)^T$, together with corresponding observations of the values
  of t, denoted $t \equiv (t_1, . . . , t_N)^T$. Figure 1.2 shows a plot of a training set comprising
  N = 10 data points. The input data set x in Figure 1.2 was generated by choosing
  values of xn, for n = 1, . . . , N, spaced uniformly in range [0, 1], and the target
  data set t was obtained by first computing the corresponding values of the function
  1.1. Example: Polynomial Curve Fitting 5
  sin(2πx) and then adding a small level of random noise having a Gaussian distribution
  (the Gaussian distribution is discussed in Section 1.2.4) to each such point in
  order to obtain the corresponding value tn. By generating data in this way, we are
  capturing a property of many real data sets, namely that they possess an underlying
  regularity, which we wish to learn, but that individual observations are corrupted by
  random noise. This noise might arise from intrinsically stochastic (i.e. random) processes
  such as radioactive decay but more typically is due to there being sources of
  variability that are themselves unobserved.

I'm not clear what the notation  $x \equiv (x_1, . . . , x_N)^T$ might mean here, especially the so called variable T.
I guess it means it's the set of input vectors together with its corresponding set of target vectors (t). If I am right, so is this an example of  supervised learning? Still, what's T?


texts: Nyquist noiseless channel capacity; how can bit-rate be two times the bandwidth?

texts:
I'm confused by the Nyquist channel capacity formula.  How can channel maximum capacity approach double the bandwidth.

$C = 2\times BW \times log_{2}(L)$ bits/sec

The way it was explained to me way back when when I took a telecommunications course was as follows.

C = Maximum channel capacity in bits/second

BW = (highest frequency - lowest frequency) that the channel can accommodate.

L = signal levels So in the figure below, where BW = f; the frequency of the sine wave; and L = 2 for -1 and + 1; the channel capacity is C
= 2B because we can consider the positive part of the sine wave as = to a 1 bit for half the period, and as a 0 bit for the negative half of the period.

Now i'm not sure that this is correct, as that doesn't seem like "information", being that it doesn't seem to decrease entropy, or uncertainty?

so if that's not correct, what is an example of a signal with capacity twice the bandwidth?
Or am I just not understanding the entire concept?


texts: How do you measure "detail" of a signal?

texts:
I have an image and I would like to measure the amount of detail in it. Another way to look at it is to measure how blurry an image is. One way is to analyse the high frequency components in the Fourier transform of the image.

Are there any other/better methods?


texts: How many bps would be required to reproduce the words in speech recognizably?

texts:
In Robert Gallager's lecture notes for 6.450 Principles of Digital Communications I, Exercise 2.1 asks an interesting conceptual question:

  Voice waveforms could be converted to binary data by sampling at 8000 times per second and quantizing to 8 bits per sample, yielding 64kb/s. [...] Modern speech coders can yield telephone-quality speech at 6-16 kb/s. If your objective were simply to reproduce the words in speech recognizably without concern for speaker recognition, intonation, etc., make an estimate of how many kb/s would be required. [...] (Note: There is clearly no “correct answer” here; the question is too vague for that. The point of the question is to get used to questioning objectives and approaches.)

My first instinct was to look up number of words used today in English language (~170,000), calculate how many bits would be required to brutely encode each word (~18), and look up, on average, how fast we speak in terms of words per second (~2), and come up with 36 bits per second as my crude answer. I'm not sure if not concerning with the speech waveform itself, and all the signal processing that would entail, is a reasonable way of thinking about the problem. I'm here for feedback and other ways of thinking about this problem.

(This is not a homework assignment. I was studying the notes myself and I'm curious.)


texts: Conditions for analytical signal reconstruction with the Hilbert Transform (HT)

texts:
I'm new to the Hilbert Transform, I was reading the paper: "SISAR Imaging - Radio Holography Signal Reconstruction Based on Receiver-Transmitter Motion", and the next sentence was written in the paper:

  The HT states that an analytical complex signal can be
  reconstructed from the signal’s real part under the condition
  that the real amplitude and frequency modulated signal in the
  form:
  $$
S(t) = I(t)\cos(\theta(t))
$$
  has the spectrum of $\cos(\theta(t))$ outside the spectrum of its envelope $I(t) \in [−f_0, f_0]$.

I have checked some documentation on the Hilbert Transform and I have not been able to find that condition, so, why is this condition needed? What is the reason behind it?


texts:
I'm trying to implement a software audio modem with FSK to transmit and receive short digital data between two devices (e.g. a PC and an iPhone) so I'm simulating simple cases study using Scilab/MATLAB.

While I can reconstruct the original signal when no noise nor background sound is present, problems arise when playing the afsk in presence of normal ambient background sounds: trying to reconstruct the digital data reading from the PC microphone is extremely difficult.

To generate the AFSK signal I use this simple script below, then I save the wav file to disk for the next steps.

function s = fsk (bit, f1, f2, fs, T)
% bit : Array of bit
% f1: Frequency 1
% f2: Frequency 2
% fs: Sampling frequency
% T : Tramsision period for 1 bit
    t=0:1/fs:T-1/fs;
    s=[];
    for i=1:length(x)
        s=[s (bit(i)==0)*sin(2*pi*t*f1)+(bit(i)==1)*sin(2*pi*t*f2)];
    end
end

The frequency f1, f2 I choose can vary: $3000\textrm{ Hz}$ and $4000\textrm{ Hz}$ or $3000\textrm{ Hz}$ and $5000\textrm{ Hz}$.
The frequency sampling is $44100\textrm{ Hz}$ and the period of transmission is $0.01\textrm{ sec}$ (for the moment I don't bother with higher trasmission rate).
The bit array is a very simple $[\textrm{0 1 0 1 0...}]$ sequence.

When I read-back the wav file I decode it using a simply coherent demodulation in these 3 step scheme:

Pass-band the signal with a FIR filter around $5000\textrm{ Hz}$ created using the command fir1(91,[2*4900/fs 2*5100/fs])
Multiply the filtered signal by the f2 (mark) carrier sin(2*pi*f2*t) (or multiply the signal with a 1 sample-delayed version of itself)
Low-pass the multiplied signal obtained in step 2 to isolate the low frequency component of the original digital signal using LP FIR filter generated with fir1(61, 2*500/fs)

The original sequence of bit is now almost perfectly reconstructed. I've played a bit around with the filters parameters to obtain a good recontruction.

Then I decided to play the wav file on my stereo with normal background noise present (actually very low), I recorded it with the pc microphone and saved with Audacity, read it back with MATLAB/scilab and I tried to reconstruct the original bit sequence with the 3 steps above but the signal is a bit different from what I expected and I can't see the original pattern and I don't know I to recontruct the bit sequence using any treshold value.

I suppose the cause is the SNR is awesome but I don't know where I'm wrong and I don't understand why the results are so different than the previous case.

How could I reconstruct the correct bit stream sequence in presence of background noise ?

UPDATE:
Do I have to low-pass filter the modulated signal before send it to the receiver ?

Any suggestion will be very appreciated.

Well...it seems that my recording system was buggy ,as Jim pointed out !
I played and recorded the sound using the speakers and mic on the same PC so this introduced some indesiderable feedback or filtering or something else I suppose.
Now I'm using an external device to play the sound (a smartphone) and the recorded signal is much better. This is the signal transmitted

and this is the decoded signal

I'm using a bandpass filter centered around the mark and space frequency before transmit it. The decoding steps are the same as before but in step 2 I'm using the slicing method correlating the signal with the delayed version of itself.

The reconstructed FSK signal has now a much better SNR than before. Now should be easy to recover the exact bit sequence....


texts: use of complex conjugate of weight vectors in beamforming literature

texts:
in beamforming literature, people generally quote the following as a starting point:

$y(k) = \sum_{i} w_i^*\cdot x_i(k)$

Here:

i varies from 1 to N where N is the number of sensors

$w_i^*$ is the complex conjugate of the ith weight vector

$x_i(k)$ is the kth sample received at the ith sensor

y(k) is the output of the beamformer for kth sample

The question is why is the complex conjugate of the weight vector $w_i^*$ used(instead of the complex weight vector $w_i$ iteslf). The literature seems to suggest that it simplifies notation but I can't see how.


texts: Interpret periods from fourier transfrom graph of periodic impulse signals

texts:
I have obtained the fourier transform of a signal with following graph

My question (though might be too junior for the professionals in signal processing field) is: how can we extract the period of 30 from the fourier transform mag graph?

The source code is

Fs = 1;
Ts = 1/Fs;
dt = 0:Ts:30*4-Ts;

 y1 = [
  100 ; zeros(Fs*29,1);
  100 ; zeros(Fs*29,1);
  100 ; zeros(Fs*29,1);
  100 ; zeros(Fs*29,1);
  ];

subplot(4,1,1)
plot(dt, y1, 'r')
title('impulse happen in time series', 'fontsize', 18);

X = fft(y1);
X_mag = abs(X)
subplot(4,1,2)
plot(X_mag)
title('fourier transform mag', 'fontsize', 18);


texts: Remove noise introduced while converting a Delta-Sigma Modulated (DSM) signal into a pulse-code modulated (PCM) signal?

texts:
I'm writing a tool that converts a dsdiff audio file (dsd64) into a flac audio file (DSM -> PCM). The project is both a programming exercise and a signal processing exercise; I'm learning a lot from it.
Currently, I'm doing the following to compute a sample in the PCM stream

Select a window of the DSM sequence that corresponds to the size of the window of a low-pass (FIR) filter.
Replace all 0s in the DSM window with -1s
Convolve (multiply+sum) the elements in the DSM window with the coefficients of the (FIR) filter and rescale the result to fir into the chosen bitrate (24-bit).
Move the window by downsampling factor many elements (16 in my case) and apply repeat this process.

I tested the conversion using a few (arbitrarily chosen) dsdiff files (songs by the band deep purple).
The conversion reproduces the signal (good) but introduces noise (bad) with the following profile (measured during the beginning of the track where there is a bit of silence in the original):

I've been trying to figure out how this is being introduced, and what I need to do to remove it. From my reading of random sources online I understand that this is introduced by integrating the DSM signal, which amplifies the noise floor and maps it to a spectrum similar to the one I see in my plot.
Did I identify the problem correctly? What can I do to verify this and, ultimately, fix it?
Any help is highly appreciated.
Edit:
Added a figure showing the filter's frequency response

Added the Periodogram of the first second of the DSM stream


texts: How to calculate linear index of fuzziness of an image?

texts:
i am trying to calculate  linear index of fuzziness of an image using this formula:
i have to see a positive value according to the paper. But i see a negative value, so i think i am doing something wrong. here is my try.

q=(pi/2)*(1-output_image/255); %I is output image and 255 is maximum gray level
y= sum(sum(min(q,1-q)))/(width*height); %M*N is width*height


texts: Autocorrelation and Power Spectral Density (Discrete)

texts:
The Autocorrelation, $\phi_{aa}[\kappa]$, of a discrete time random process, $a[k]$, is defined as:

$$ \phi_{aa}[\kappa] = \mathrm{E}\left\{ a[k+\kappa]a^*[k] \right\} $$

Taking its fourier transform would give the Power Spectral Density, which is given in my book as:

$$ \mathscr{F}\big\{ \phi_{aa}[\kappa] \big\} = \Phi_{aa}\left(e^{j2\pi fT} \right) = \sum\limits_{\kappa = -\infty}^{\infty} \phi_{aa}[\kappa] e^{-j2\pi\kappa fT}  $$

Now, why don't they denote simply as: $\Phi_{aa}[\kappa]$?

Is the form above with the exponential signify something related to sampling interval? I'm quite confused with this form and any explanation to better understand this form would be helpful.


texts: calculating estimate of a state for a system with two observations of the state from different times

texts:
I'm struggling with a problem that I just can't seem to get a grasp of. I'm supposed to calculate an estimate for state $x(k)$ at times $k=1,2,3$ from the state-space

$
\begin{align}
x(k+1) &= A  x(k) + B  u(k) + w(k) \\
y_1(k) &= x(k) + v1_(k) \\
y_2(k) &= x(k-1) + v_2(k) \\
\end{align}
$

where the covariances of the noises are given as

$E[w(k)w(k)']$ , $E[v_1(k)v_1(k)']$, $E[v_2(k)v_2(k)']$

and

$u(k)$ is given for $k=1,2,3$

$y_1(k)$ is given for $k=1,2,3$

$y_2(k)$ is given for $k=2,3$

What kind of technique am I supposed to use? Information Kalman filter?


texts: Cascade Transfer Function; Odd or Even Transfer Function

texts:
suppose we have a cascade realization of a transfer function of a higher order.

What difference does it make if the order of the transfer function is an Even or an Odd number?


texts: How to get the Fourier series using Python's $\tt fft$

texts:
I Would like to be able to reconstruct every individual sinusoid that makes up a Discrete signal. I Have the following signal: (I am working in Python)

The signal is essentially an array with about 400 elements that varies with time. In Python after calling the fft function on the data

 dft= rfft(dat)/len(dat)  #real fft

I receive the figure below:

I am aware that I can use the result of the fft to obtain the individual Fourier series components, but I am unsure exactly how. The code I have yields me an array of magnitudes for each harmonic of the signal, how do I use that to get each individual component of the Fourier series?


texts: FFT analysis of binned data of different length

texts:
I am conducting experiments to collect wind speed data from wind anemometers placed on a moving platform. Closely, there is a fixed wind mast holding a wind vane. Prior to analysing the wind data from the moving platform, I binned the data with 1-minute average wind direction. Each bin contains a number of continuous data sets (Set 1, Set 2, Set 3 etc.) which may not necessarily have the same number of data points. For example, Set 1 has 20,000 data points, Set 2 has 6,000 data points, etc.

I am performing an FFT on each of the sets to detect any periodicities in the signal. I am plotting the real amplitude of the sinusoids present in the signal on the y-axis vs. the frequency of the sinusoid. The FFT of each data set is yielding different amplitude results, on the y-axis and I'm not sure whether this is a result of the different length of the arrays. I am interested in an average value of the y-axis of the spectrum. However, I don't think that it is suitable to just calculate the average due to the different array size. Can someone please suggest a method? I am not used to analysis of signals in this way so I'm a bit confused. Thanks a lot.


texts: Estimate "speed" of repeating signal pattern

texts:
I have a digitized signal containing a repeating pattern. One could call it periodic, but the time for one period is not constant. The "speed" with which the period is advanced can change over time. How can i estimate/calculate this "speed" function if i know how one period with constant speed looks like?

An example:
Take one period of a sine wave. I know how this period of a sine wave looks like when plotted with a constant frequency. Now i record a signal where this sine wave changes frequency constantly. How can i calculate the frequency for each point?

My concrete example:
The signal represents the IR light reflected by the metallic plate in an energy meter (Ferraris type). There is a red colored area on this plate which i can see as a dip in signal intensity. One rotation of the plate is 1/75 kWh consumed. The interesting thing: My IR detector is so sensitive that i can see more patterns on this metal plate than just the dip from the red area. The "speed" at which this pattern evolves is basically the current power draw, which i want to calculate. If i just count the passing of the red area, i get only a very course resolution (one minute average) of the power consumption, but if a manage to look at the whole period in total i think i can get down to second-resolution. I hope this makes some sense.


texts: Control theory: how do you initialize input for a model predictive controller?

texts:
Although I think this is a pretty popular algorithm, I posted this question on electronics stacksexchange but received no response. I am posting it here to see if there are many signal/control experts who can answer my very simple question

On these slides, there is a description of MPC algorithm. I am confused as to how the input $u$ is initialized and how its values are known in the future.

How is $u$ initialized at the beginning?

In the  very beginning, we need all the future values of $u$ in order to predict the value of $x$. I am fully aware of how $x$ is predicted $N$ steps into the future.

But where do you obtain $u_k, u_{k+1}, \ldots, u_{k+N-1}$ to predict the states in the first place?
For instance how do you know $u_{k+1}$?

How to get the next $u$ values after receding horizon step?

After obtaining the optimal $u^*$ from the optimization program, we keep the first element of $u^*$ and discards the rest. This allows us to produce our output $y$.

But how are the subsequent $u$ generated in order to feed the next iteration of the program?
Again, we need every single value $u_k, u_{k+1}, \ldots, u_{k+N-1}$, but for instance how do you know $u_{k+1}$? (This is the same question as question 1 I suppose)

All in all, I am not sure how we know the future of our input $u$ in order to predict $x$. Can someone please enlighten me as to how this is done in simulation or in practice?


texts:
I'm working on an application where I want to grab high-resolution pieces of the spectrum, but I don't need the whole spectrum-- just a few bands. I've looked into using the Goertzel algorithm (too slow; with the number of points we want the regular FFT is faster), pruned FFTs (doesn't look like much of an improvement in performance), and polyphase filter banks (PFBs) (the problem being that it requires an FFT of the same size that I would have otherwise already computed, defeating the point).

The books that I've looked at so far ("Digital Signal Processing, 4th ed" by Proakis/Manolakis and "Signal Analysis: Wavelets, Filter Banks, Time-Frequency Transforms and Applications" by Mertins) both focus on either "critically sampled", where the decimation rate $D$ equals the number of subbands $N$ or "oversampled" where $N > D$.

Is there a way to produce a computationally efficient PFB such that $D > N$? In my application the loss of information is unimportant because I only need a few subbands of small bandwidth.


texts: difference between Gabor and log-Gabor function

texts: I am reading a paper using log-gabor filters for feature detection. I was thinking about the difference between Gabor filters and log-gabor filters. Can anyone tell me the difference(s), and a way to implement them? Thanks in advance.


texts: Image Interpolation Using the Yule Walker Equations

texts: I have been studying about the Yule-Walker equations for prediction of a time series data from knowledge of past values of the series. Is there any way I can use the same in an image to exploit the correlation between image pixels in order to perform image interpolation,, or to increase the resolution of the image?


texts: Trying to run an implementation of lsqcurvefit from the Optimization Toolbox from Matlab in python

texts:
trying to implement lsqcurvefit from matlab in python with no success. Below is the matlab code I am trying to port to python:

myfun = @(x,xdata)(exp(x(1))./ xdata.^exp(x(2))) - x(3);
xstart = [4, -2, 54];
pX = [2, 3, 13, 12, 38, 39];
pY = [12.7595, 8.7857, -11.8802, -10.9528, -15.4390, -15.3083];

try
    fittedmodel = lsqcurvefit(myfun,xstart,double(pX),double(pY), [], [], optimset('Display', 'off'));
    disp("fitted model:");
    disp(fittedmodel);
catch
end

Below is my matlab output:
fitted model:
    4.8389    3.3577   -2.0000

Below is my python code:
from scipy.optimize import curve_fit
pX = [2, 3, 13, 12, 38, 39];
pY = [12.7595, 8.7857, -11.8802, -10.9528, -15.4390, -15.3083];

def myfun(x, xdata):
            temp_val_1 = np.exp(x[0])
            temp_val_2 = np.exp(x[1])
            temp_val_3 = x[2]
            temp_val_4 = np.power(xdata, temp_val_2)
            temp_val_5 = np.divide(temp_val_1, temp_val_4)
            temp_val_6 = temp_val_5 - temp_val_3
            return temp_val_6

        popt, pcov = curve_fit(myfun, pX, pY, p0=([4, -2, 54]))
        print(popt, "\n", pcov)

and below is my python output:
myfun() takes 2 positional arguments but 4 were given

Any help on this would be greatly appreciated, thank you.


texts: Doing resampling, how can we choose lower cutoff frequency for LPF filter?

texts:
I would like to resample an audio signal.

Sampling frequency of an original audio signal is 50k Hz.

What I want to down sampling frequency is 30k Hz.

In this case, I followed the procedures as below using LabVIEW.

  '1. IIR low pass filter

  Butterworth / 6th / Lowpass / up to 15k Hz

  '2. Resample to 30kHz

Why I used LPF is to prohibit the anti-aliasing.

 from now on, I have some questions.

1. Do we need to do LPF for prohibiting anti-aliasing?

2. If we have to do, how can we choose the lower cutoff frequency?
   In my case, I choose as 15k Hz because this is the new Nyquist frequency.

Thank you in advance.

OR is there the standard procedure for resampling?


texts: Welch's Power Spectral Density - Time-Averaging Explanation

texts:
I am just branching out to taking the power spectrum of short-term audio frames ($20$ ms) in order to extract useful audio features.  I have been reading about Welsh's method which states that after computing the squared magnitude of the result, the individual periodograms are time-averaged which reduced the variance of the individual power measurements.

I was hoping somebody could explain what and how time-averaging is achieved?  Is it done by looking at the individual squared bins within a given frame, or does it compare a previous frame to a present frame like a Spectral Flux?  Also, I take it that a periodogram refers to the periodicity of a given audio frame?  Thanks.


texts: How can I detect pipeline cracks using OpenCV and Python?

texts:
I have developed a robot that captures images of the pipeline interior as it moves. The requirement was to be able to detect cracks inside. So far i tried several OpenCV codes that find the crack contours but i was not successful.

Code I'm working on:

import cv2
import numpy as np
image = cv2.imread('pipe_photo1.jpg')
blurred = cv2.pyrMeanShiftFiltering(image,41,91)
gray = cv2.cvtColor(blurred,cv2.COLOR_BGR2GRAY)
ret, threshold = cv2.threshold(gray,0,255,cv2.THRESH_BINARY+cv2.THRESH_OTSU)
_, contours, _ = cv2.findContours(threshold, cv2.RETR_LIST, cv2.CHAIN_APPROX_NONE)
print (len(contours))

cv2.drawContours(image,contours, -1,(0,0,255),6)
cv2.namedWindow("Display",cv2.WINDOW_NORMAL)
cv2.imshow("Display",image)
cv2.waitKey()

Below is the image i obtained from the camera. I want to detect only the crack shown at the bottom of the pipe and be able to draw it using red lines. Your help will really save me in achieving my objectives before its due.


texts: Reconstructing audio from a given spectrogram and phase data

texts: I have 2 seperate files called: spectrogram and phases. Spectrogram contains magnitude information of audio and phases contain the phase information of each windowed area. I need to reconstruct the audio using these information yet i am completely lost. Can you put me into the right direction. How can i do this?


texts:
Relation of Entropy and SNR : Based on this question and answer, I had another question that struck me and I am curious to know, if somebody can shed some light, on the following situation: $y= desired_{signal} + noise$ is received by the receiver and the $error = y - \hat{y}$ where $\hat{y}$ is simulated at the receiver by using guessed parameters with the knowledge of the kind of model used at the transmitter.

For increasing SNR (attenuating noise) would Entropy of error decrease? In the link provided, the entropy of the transmitted signal decreases with increasing SNR since the uncertainty decreases when noise is getting attenuated. But what about the entropy of error? Should error entropy increase with increasing SNR? I am confused about this.

EDIT: The way I calculated entropy of error:
Let the system model be AR(2): $y(t)= a1y(t-1) + b1y(t-2) + noise$. At the receiver end, I have $z(t) = y(t) - (a2y(t-1) + b2y(t-2))$ where $(a2,b2)$ are close guesses to ($a1,b1$). I calculated 2 entropies H1 for y(t)  and H2 for z(t).

z(t) will be close to zero when a2,b2 will be equal to a1,b1 in which case entropy of z(t) is found to be minimum among all different z(t) for a particular SNR value of noise.
 So, for a particular noise value, I have 10 pairs of (a2,b2) and for each I get 10 values of H2(z(t)). I chose the minimum among all H2 for a particular SNR level. But, for different noise levels, as I increase the SNR at the transmitter end, the minimum entropy, H2(z(t)), increases with increasing SNR.
I found that H1 decreases with increasing SNR but the reverse happens for H2.
Is this trend correct?