texts:
I am attempting to implement a cascaded IIR filter in MATLAB.  I used fdatool to design the filter and exported it as an object, Hd.

In my script, I import a log file of data I captured and store my data in a vector called serial.  I then use MATLAB's filter(Hd,serial) to filter it and plot the results.  The filter works well.

In order to check and improve my understanding of how to implement this filter in real life (before I attempt to embed it using C) I want to filter my data using my own code and compare the results with those obtained using MATLAB's filter function

My Hd object consists of the following parameters:
SOSMatrix:

0.25, 0.5, 0.25, 1, -1.9513, 0.95189
0.25, 0.5, 0.25, 1, -1.9584, 0.95907
0.25, 0.5, 0.25, 1, -1.9718, 0.97246
0.25, 0.5, 0.25, 1, -1.9896, 0.99024

Scale Values:

1.4055e13
1
1
1
1

In order to implement my filter I have stored each line of my SOS matrix under the following titles:

    Section1 = [0.25 0.5 0.25 1 -1.9513 0.95189]
    Section2 = [0.25 0.5 0.25 1 -1.9584 0.95907]
    Section3 = [0.25 0.5 0.25 1 -1.9718 0.97246]
    Section4 = [0.25 0.5 0.25 1 -1.9896 0.99024]

My attempt at implementing the filter is given below (as I said earlier, my input data is called "serial"

First I initialize the outputs from each of the 4 stages to zeros...

    output_section_1=zeros(size(serial,1),1);
    output_section_2=zeros(size(serial,1),1);
    output_section_3=zeros(size(serial,1),1);
    output_section_4=zeros(size(serial,1),1);

I multiply my input (serial) by the first scale factor.  As the rest of the values are 1, I don't bother multiplying any of the other sections.

    serial=serial*1.4055e13;

Then I run a loop for my filter...

    for i=10:size(serial,1)

% Section 1
output_section_1(i)=Section1(1)*serial(i) + Section1(2)*serial(i-1) + Section1(3)*serial(i-2) - Section1(5)*serial(i-1) - Section1(6)*serial(i-2);

% Section 2
output_section_2(i)=Section2(1)*output_section_1(i) + Section2(2)*output_section_1(i-1) + Section2(3)*output_section_1(i-2) - Section2(5)*output_section_1(i-1) - Section2(6)*output_section_1(i-2);

% Section 3
output_section_3(i)=Section3(1)*output_section_2(i) + Section3(2)*output_section_2(i-1) + Section3(3)*output_section_2(i-2) - Section3(5)*output_section_2(i-1) - Section3(6)*output_section_2(i-2);

% Section 4
output_section_4(i)=Section4(1)*output_section_3(i) + Section4(2)*output_section_3(i-1) + Section4(3)*output_section_3(i-2) - Section4(5)*output_section_3(i-1) - Section4(6)*output_section_3(i-2);
    end

The output obtained using MATLAB's filter function looks good and is nice and smooth as would be expected.

The output obtained using my own implementation looks just like the input (noisy) except that the values have changed and there is an initialization at the beginning.

I'm sure I've made a very schoolboy mistake but I cannot see where I have gone wrong.

Any advice anyone could give me would be much appreciated.


texts: Impulse Response with non-delayed feedback

texts:
I've been trying to solve a Discrete LTI System problem where the overall impulse response of the following System, whose input and output sequences are, respectively, $x(n)$ and $y(n)$, has to be found:

This is what I've tried:

Calling the output sequence that comes out of the adder $e(n)$, I know that $y(n) = e(n) \star h(n)$. I also know that $e(n) = x(n) - y(n)$, which leads to $y(n) = x(n) \star h(n) - y(n) \star h(n)$.
Thats the point where I got stuck. I don't know how to proceed in order to isolate $y(n)$ in the last expression. I would like to isolate it in order to make $x(n) = \delta(n)$ and evaluate the overall Impulse Response.

Also, I am not supposed to use any transform techniques as the book I am following has not yet arrived at those subjects. If possible, don't show the final calculations, only what I should do in order to separate $y(n)$. I would like to figure the rest out by myself.

Thanks in advance.


texts: Why we have to multiply the images with $(-1)^{x+y}$ to center the transform in frequency domain filtering?

texts: To perform the filtering infrequency domain we perform multiplication of $(-1)^{x+y}$ why??


texts: PSNR of two images of different size in matlab

texts:
I performed bicubic interpolation on a 256*256 image(img)

 dest = interp2(img,'bicubic')

and  i got a 511 * 511 image.I want to compute PSNR of a 512 * 512 image(original) and the 'dest' image as follows

original = double(original);
dest = double(dest);
[M N] = size(original);
error = original - dest;
MSE = sum(sum(error .* error)) / (M * N);

if(MSE > 0)
PSNR = 10*log(255*255/MSE) / log(10);
disp(['PSNR = ', num2str(PSNR)])
else
PSNR = 99;
disp(['PSNR = ', num2str(PSNR)])
end

But I'm getting error due to different matrix dimensions.How to avoid this error.Is it possible to calculate PSNR of images with different size?Please help


texts:
I have some code that collects amplitudes from a complex sound like a song. I now want to plot a graph of amplitude vs time, for the graph plotting I use JavaFX. However I get about 9 million amplitude shorts for each channel, so how to I select a few samples to plot (as plotting 9 million points is obviously out of the question)? Do I need to use the sample rate or sample size? Is there a particular way I should go about this?

Is this post more suited for StackOverflow?

AudioInputStream in = AudioSystem.getAudioInputStream(file);
AudioInputStream din = null;
 AudioFormat baseFormat = in.getFormat();
  AudioFormat decodedFormat = new AudioFormat(AudioFormat.Encoding.PCM_SIGNED,
  baseFormat.getSampleRate(), //sample rate
  16,                         //sample size (bits)
  baseFormat.getChannels(),   //channels
  baseFormat.getChannels() * 2,//frame size
  baseFormat.getSampleRate(), //frame rate
  false);                     //little-endian
  din = AudioSystem.getAudioInputStream(decodedFormat, in);
  byte[] temp = new byte[4];

  //....
  //leftDos and rightDos are streams of shorts
  //....

  int i = 0;
  while (din.read(temp, 0, 4) != -1) {
      if (decodedFormat.getChannels() == 2) {
          leftDos.writeShort(temp[1] * 256 + temp[0]);
          rightDos.writeShort(temp[3] * 256 + temp[2]);
          i++;
      }
  }

 //...retrieve shorts...//
 short[] left = new short[(bytes.length / 2)];

 int len = 0;
 try {
     while (1 > 0) {
         left[len] = leftDIS.readShort();
         len++;
     }
 } catch (EOFException e) {
      System.err.println("End of stream");
 }


texts: Determine minimum sample rate for continuous sinusoid

texts:
Consider a signal $$ x(t) = \cos(175\pi t) $$ which is sampled to produce discrete time signal $$ x[n] = x(nT_s) $$ The fundamental period of $x[n]$ is $$ N_0 = 7 $$

Given this, what is the smallest possible sample rate $T_s$? (Ans: 1.6327 ms).

I would assume that this is related to finding the Nyquist frequency. I was thinking:

Since, $$N_0 = 7 \implies f_0 = \frac17 \implies f_{\mathrm{Nyquist}} = 2 \frac17\implies T_s = \frac72 $$ However, this is obviously incorrect. I am not even using any information of the original signal. Any suggestions on what I could be doing wrong here?


texts: derivative filter for computing speed from robot rangefinder

texts:
I have a laser rangefinder on my robot and I need to compute the speed of the detected object. I have seen many different filters for that, but all of them apply only to offline data, or introduce delay as they rely on past and future points, like (Xn+1 - Xn-1)/2Tn.

I am using this signal for controlling the position of a drone so no delay is allowed. When I researched, the type of filter that I need is described as:

Online
Backward-looking
Causal

But I never found an actual implementation of such a filter.

I guess a Kalman filter would fit, but that's a big beast and I am not familiar with its math, so if something simpler existed I would be happy.

How do I implement such a filter ? Is there one such filter already implemented in a library, in any language (but preferably Python) ?


texts: DSP interview question: use of the identity in development of a significant transform

texts:
I'm preparing interview and found this question. But I don't really understand what is the question. Does it ask about Fourier transform or Z transform?

How the simple identity

$$xy=\frac{1}{2}x^2 + \frac{1}{2}y^2 - \frac{1}{2}(x - y)^2$$

was used in developing a significant transform in digital signal processing?


texts:
I am generating sine wave in MATLAB with following code

fs=2000e6;%sampling frequency
t=0:1/fs:(1e-6)-1/fs; %1us seconds duration
fc=501e6;%frequency of the cosine wave
x=cos(2*pi*fc*t);
plot(x)

The wave it generates is somehing like this

My question is: Why do we see amplitude modulation in the wave. Shouldn't amplitudes be 1 at all peaks of the wave? like this wave


texts: Finding the Best Gaussian Smoothing Kernel to Minimize the Discrepancy Between Two Images

texts:
Suppose we have two grayscale images, $A$ and $B$. $A$ and $B$ very strongly resemble each other, such that the mean of the absolute difference $\lvert A - B\rvert$ is fairly small. Suppose further that $B$ actually appears to be a blurred version of $A$, although it was not necessarily obtained through a blurring process. I'd like to find the ideal Gaussian kernel $g$ to smooth $A$ with to maximize its correspondence to $B$, in other words, select $g$ s.t.

$$\underset{g}{\arg\min}\sum{\lvert(g*A) - B\rvert}$$

where $g$, $A$, and $B$ are 2D arrays and $*$ is the convolution operator.


texts: What does double hat mean in Spectral theory?

texts:
I was confused about the double hat notation $\hat{\hat .}$ and my teacher and I thought it probably means an estimate of an estimate most of the time, but we could not be 100% sure of its meaning.
One hat $\hat .$ is an estimate but what does this double hat usually mean? Does it have many different meanings in DSP? Can it mean an estimate of unit vector sometimes? Does it have other signs whatever it means?
Example

You can see such double hat on page such as 213 in Prentice Hall - Modeling Of Dynamic Systems book. In the chapter on Periodograms: Blackman-Tukey's Method, they use this kind of spectral estimate for N length function:
$$\hat\phi_N(\omega)=\int_{-\pi}^{\pi}W_\gamma (\omega-\xi)\hat{\hat{\phi}}_N(\xi)d\xi$$


texts:
What is the fourier transform for finite aperiodic discrete signals and for infinite aperiodic discrete signals?

Are they different?


texts: Definition of Perfect Code (and example of Hamming Code)

texts:
Let's consider the following sentences about coding a message before transmitting it in a communication system (here you find the complete pdf)

First question: can you explain me these two definitions of perfect code? I do not understand what "bit patterns" mean (I think they are the received bits but I am not sure). And I do not understand what do they mean practically.

Now let's consider the Hamming Code (7,4): here the generator matrix G (such that C = B*G, by referring on the upper blocks scheme) is shown:

Second question: I do not understand these two sentences. The rows of G do not appear to be one bit distant from a codeword and the second sentence is completely obscure for me. Can you explain me them?


texts: I read that pole zero cancellation is a method to get rid of poles outside the unit circle. But I believe exact cancellation could be much harder because of finite precision effects. Are there any other reliable methods to cancel out poles? Also are there any methods to shift poles (and do some required compensation)such that the frequency response is unaffected? Thanks in advance.


texts: Why does the standard waveform from audio out get distorted when I reduce the amplitude?

texts:
In my attempt to generate a standard sine wave out of a sound card, I started using matlab which has a function soundsc for scaling the output

 t = 0:1/48000:10;
 y = sin(1000*2*pi*t);
 soundsc(y,48000,[a b]);

It's ok when my a and b are -1 and 1 respectively. When I try to make them even -0.9 and 0.9 respectively, the standard sine wave of 1Khz gets distorted, why?

Are our sound cards not able to handle standard sounds?


texts:
I try to detect rectangles in images. The background of the images is
one color (most of the time). I tried two methods to get an binary image (1 = background, 0 = edges), to do an Hough Transformation later on...

Sobel or Canny Filter
Smooth image A, Create difference image A - gauss,
Create binary image with threshold (Create Histogram, highest bin should be    background...)

The result is a binary image with edges. I don't really now which method works better
for a variety of different images. Any ideas?


texts: How can I get the function of a curve from a dataset without using a curve fitting tool?

texts:
previous relevant question here. I have a dataset of a signal captured from real hardware that produces a curve. Depending on the curve characteristics I can determine whether metal is present or not within a metal detector. The curve looks like this:

Each black line represent the average of a number of curves for each metal sample that goes through the system. The bottom curve being no metal present and the top, biggest curve being the largest metal piece present.
I want to remove the most significant portions of the curve and use definite integrals to help me decide whether metal is present. I can do this in matlab by putting the dataset of each average through the curve fitting tool, which generates the function (double exponential, 99.9% fit) of the curve and then from this I can put in my limits and integrate to get the area under the curve. I can then compare these areas to get good results.
My question: Is it possible to get the function of the curve without using a curve fitting technique?  Or is it possible to integrate the curve without curve fitting? I want to be able to analyse the signals with a (powerful) microcontroller in the end design and I do not want to write a curve fitting algorithm (time limited).


texts: Why is cepstrum used in MFCC instead of autocorrelation sequence?

texts:
I know that the inverse Fourier transform of the power spectral density $P(k)$ gives the autocorrelation sequence $A(n)$ (By Wiener-Khinchin theorem):

$${\cal F}^{-1}\left\{P(k)\right\}=A(n)$$

And if we replace $P(k)$ by its log value, we get the cepstrum sequence $C(n)$:

$${\cal F}^{-1}\left\{\log_{10}P(k)\right\} = C(n)$$

My questions are:

What's the reasoning behind using log value?
Why does that make a difference?
Why isn't the autocorrelation sequence used in MFCC, LPCC, PLPCC, etc.. instead of the cepstrum?


texts: What is the relationship between passivity of impulse response and passivity of the system?

texts:
What I mean by passivity of impulse response is that the energy of impulse response for this system is less than the input impulse signal energy.

For example, in a discrete-time system, considering the impulse response, the input is Kronecker delta function $δ[n]$, and the output is the impulse response $h[n]$. Suppose that
$$
h[n]=\sqrt2/2 , n=1,2
$$
and $h[n]$ is 0 for any other components. In discrete-time system, We calculate the energy of the signal with
$$
Energy = \sum_{n=-\infty}^{\infty} |h(n)|^2
$$
Therefore, the energy of the impulse response is 1, and the input signal energy($δ[n]$) is also 1, which means the system is passive for this input signal.

Then considering another input signal $δ[n] + δ[n-1]$, the energy of this input is 2. The output can be derived by convolution of impulse response and this input signal. We get:
$$
y[n]=\sqrt2/2 , n=1,3
$$
$$
y[n]=\sqrt2 , n=2
$$
and $y[n]$ is 0 for any other components. The energy of this output is 3, larger than the input energy. Therefore, It seems that the passive impulse response doesn't necessarily  make a system passive. But if a system is passive, the impulse response of the system should always be passive. Am I right with this conclusion? or can anyone give a rigorous proof for this?


texts:
Is there a transformation that will enable one to calculate the FFT in an arbitrary coordinate system? What I am interested is the following two cases:

The space is Euclidean and infinitely differentiable everywhere
The space forms a curved Riemannian manifold with metric:
$$ds^2 = \sum_{i, j} g_{ij} dx^i dx^j$$


texts: How do you normalize a sample-by-sample processed signal?

texts:
I am writing an audio effect based on a second-order filter that adds two samples. The resulting output is clipped. So I multiplied the signal by 0.5, which results in an output that is too low. Since I am process sample-by-sample, I don’t see how you can apply normalization because the peak-values in the normalization window will vary too much.
How do you do this? Calculate the right audible half volumes in dB and add those?


texts: Compute Hann window without cos function

texts:
In an environment with limited memory and computing power it is interesting to be able to generate a Hann window without using a cache or repetitive calling of expensive functions such as sine and cosine functions.
The following document demonstrate the use of an alternative approach in order to compute a pure tone by storing only a coefficient and the two last results:
http://ww1.microchip.com/downloads/en/appnotes/00543c.pdf
Rewrite using python code:
import math
import matplotlib.pyplot as plt

tone = [0.0] * 50
tone_freq = 1720.0
sample_rate = 44100.0
ffs = tone_freq / sample_rate;
k1 = math.cos(2 * math.pi * ffs)
k2 = math.sin(2 * math.pi * ffs)
tone[1] = k2
for i in range(2, len(tone)):
    tone[i] = 2 * k1 * tone[i-1] - tone[i-2]

plt.plot(tone)
plt.show()

Using the same approach is it possible to generate the Hann window ?
Note:
You can generate a Hann window with sin function like that:
tone = [(0.5 - 0.5 * math.sin(0.5 * math.pi +( 2 * math.pi * i) / (len(tone) - 1))) for i in range(len(tone))]


texts: The difference between MSK as CPFSK and OQPSK

texts:
I am working on MSK simulation. I have read the Chapter 5 of book entitled Digital Modulation Techniques by Fuqin Xiong. It is said that MSK can be viewed as a Sinusoidal Weighted OQPSK or CPFSK. If we write the MSK signal as

$$
s(t)=A\cos\left(2\pi f_c t + d_k \frac{\pi t}{2T} + \phi_k\right), \quad kT\leq t\leq (k+1)T
$$

From the first perspective (MSK as a Sinusoidal Weighted OQPSK), we have
$$
d_k=I_kQ_k , \quad \phi_k=\frac{\pi}{2}(1-I_k)
$$

From the second perspective (MSK as CPFSK), we have
$$
d_k=d_{k-1}\newcommand*\xor{\mathbin{\oplus}}a_k,  \quad \phi_k=\frac{\pi}{2}\left(\sum_{i=0}^{k-1}d_i-kd_k\right)
$$
where $\{a_k\}$ represents the original data stream and $\newcommand*\xor{\mathbin{\oplus}}$ is XOR. I have the following question:

Are the two methods for computing MSK and specifically $\phi_k$ equivalent?
Do they result in the same MSK signal?

In my simulations, the two MSK signals that I get from these two perspectives don't have the same phase.


texts: How to understand relationships between ellipse and second moment matrix of Harris corner detector?

texts: guys. I am really lacking in knowledge of linear algebra. I am reading slides of Harris corner. But I am really confused about one of them. I know that I can find corners by two large eigenvalues but I can not consider it as an ellipse....can anybody give me any suggestions?


texts: Estimate the Transfer Function of an Unknown System

texts:
Suppose you have a system, H, that you want to estimate its transfer function.

You have a finite number of complex input samples, x, and noisy complex (magnitude and phase) output samples, y:

In practice, how do you estimate it's transfer function given a mean-squared error criteria?
What algorithms are used to estimate transfer functions when you have limited and noisy data?
How would you estimate the transfer function using a Wiener filter?


texts:
My understanding of Nyquist's theorem is that you need to sample at a rate that's twice the bandwidth of a signal to fully recover it.

What exactly does the bandwidth mean in this context and what happens if you sample below this frequency?

I am trying to measure a physical signal that is filtered at 100 kHz, but which only has content I care about below 1 Hz. Since I only care about the low frequency components, can my bandwidth be considered 1 Hz for this purpose? Is there any reason to sample the signal at 200 kHz?


texts: Length of CDMA code and message length vs. number of participants

texts:
What is the minimum code size for CDMA codes for synchronuous CDMA (not asynchronous CDMA!) with N number of participants?

Furthermore, suppose N parties want to uniquely transmit a message of length K0 each, where K0 is small (e.g. 16 bit) and N is very large (e.g. 1000). What is the minimum size of the total message (i.e., the message including all N different messages of size K0).

Is the minimum size of the transmission K0*N = 16000 bit or can it be made lower?


texts: What's the purpose of augmenting an image with a random "gradient"?

texts:
I have the following Python code

import matplotlib.pyplot as plt
import numpy as np

def make_random_gradient(image_shape):
    x, y = np.meshgrid(np.linspace(0, 1, image_shape[1]),
                       np.linspace(0, 1, image_shape[0]))

    grad = x * np.random.uniform(-1, 1) + y * np.random.uniform(-1, 1)

    return grad

if __name__ == '__main__':
    grad = make_random_gradient((64, 80, 3))
    plt.figure()
    plt.imshow(grad)
    plt.show()

In particular, you should note the name of the function make_random_gradient. This function produces images similar to the following one.

In the past, I followed a computer vision course. I've also looked at the images in the Wikipedia article related to the image gradient. In a way, the function make_random_gradient seems to generate "random image gradients", even though they are not black and white, like the images in the linked Wikipedia article. I am looking for a rigorous and mathematical explanation of why the function above generates a "random image gradient" and precisely in which sense it generates a "gradient".

Furthermore, the "gradients" generated by this function are used to augment images. In particular, if x is an RGB image with shape (64, 80, 3), each channel i of x is augmented, in my case, in the following way

x[:, :, i] = (1 - alpha) * x[:, :, i] + alpha * grad

where alpha = np.random.uniform(0.05, 0.15) and grad is an object returned by make_random_gradient.

So, each channel is augmented with the grad using a convex combination. Why would one do this?

I am aware of certain image augmentation techniques, such as cropping, flipping, the addition of Gaussian noise, but I've never seen this technique. Has anyone ever seen it? If yes, why would one augment images in this way? I am looking for papers or, in general, works that I've used a similar technique and that provide a rationale behind it.

I've also quickly skimmed through the relevant pages of the paper A survey on Image Data Augmentation for Deep Learning, but the paper is a little bit too long.


texts: How to turn off calculating pyramids during image processing?

texts: I am using Erdas 2014 to run a 5x5 focal median filter on a batch of images (Raster Tab > Spatial > Focal Analysis (batch)). These images will eventually be resampled to a different pixel size, and therefore I do not want to spend the processing time calculating pyramids. How can I turn off pyramids while performing batch image processing in Erdas 2014?


texts: How do you measure the accuracy of a stereo system?

texts:
I have a philosophical question of sorts for users familiar with a multi-view stereo setup. Imagine that I need to measure the dimensions of a large building and I take several photos of the building from different viewpoints. After calibrating the camera (intrinsic and extrinsic), I compute matches between points in these images and map them to 3D points. The output of this entire process is a cloud of 3D points. If I specify a scale, origin and co-ordinate axes, then I can perform measurements between points in this cloud which relate to the real-world scene (the building) that I captured. My question is this: how do I determine how accurate was my reconstruction.

Basically, this camera apparatus forms a measurement device. To compare a measurement, one typically needs to use another measurement device that is much more accurate that it. How do I measure the accuracy of a such a system?


texts: packet encoder RF packets by using LPF for determining if it's zero or one

texts:
I'm transmitting RF packets with frequency 868MHZ ,
symbol rate (it's the same as bit rate) is 50KBaud , and Im sampling them by a dongle (fs=20480KHZ) to capture those RF packets, the captured samples are IQ RAW DATA and I save them into a file , by matlab Im using this file and plot all the captured data(packets), I see the RF packets as peaks where there's packet found, otherwise noise.
I have built demodulator which the output of its values ranges between 1 and -1 , so I have arrived to a file that has values that ranges between 1 and -1 , this demodulated IQ samples I want to convert them to binary 0 or 1 , so Im trying to build encoder in Matlab code by using LPF , how could I do that? for more clarification I want to input to my LPF the demodulated IQ samples and then the output would be something like:

I want this because if I have the output of LPF then I can walk through that output and decide where I have "zero-Crossing" this means when there's a step (zero-crossing) then the value would be changed from zero to 1 or 1 to zero (depend on the value what was before) ,
how can I do that? the cutoff frequency of my LPF  must be my bit rate frequency right? how can I also do LPF in matlab? thanks alot !!!!
second question, lets assume I have the output of my LPF , then how can I walk through the output in order to decide zero/one (deciding zero crossing-high to low or low to high-), this part called FLOW ESTIMATION , how can I do that?


texts: Why there is no paper optimize the consumption power of the circuit in energy harvest?

texts:
In the optimization of  energy harvest ,some paper want to maximize the harvested power ,rate or minimize the transmitted power,so when they write the formula of optimization,they may assume the beamforming,transmitted power  power splitter ,power allocation as variable.

However,i am confused that why there is no one assume the "consumption of circuit" as variable ? All of paper are either assume it is a constant value or just ignore them.

Does anyone know or have a good explanation about we don't assume  "consumption of circuit" as variable?

In my opinion,i think the purpose of optimizing something in the energy harvesting is that we want to let the life of battery become longer when we "transmit" the signal correctly,so we will focus on something which has a relation with the transmitted power,like SNR  or power splitter,beamforming direction,etc.

Does anyone know or have a better explanation about it?Because we do need a circuit to transmit the signal,and the signal does consumpt the power


texts: Processes/Transforms involved to get brainwave data from raw EEG? (Autocorrelation confusion)

texts:
Not clear on what the autocorrelation function of raw EEG means physically
why can't you take the FT of a the EEG itself and get frequency data?
With BCI & basic electrode setups you can understand your brain activity:
neuron-firing --> EEG voltage reading --> brainwave readings
How I understand that this is done:
Fourier Transform: Operation bringing you from TIME domain to FREQUENCY domain

Time Domain: EEG raw ... x= time , y= volts
Frequency Domain: brainwave data... x= frequency, y=intensity

Fourier Theorem ... basically that every line/wave can be broken down into pure components (single frequency sinusoids)
First you need the autocorrelation of the raw EEG
- Correlate EEG signal with itself - compare EEG signal at one time pt -- to another time pt... shifted tau
- Get the nature of the EEG signal through time ... how voltages are shifting
FT of Autocorrelation --> PSD
We treat the signal as stationary, even though in reality it is moving in time and apply a windowing method like the Welch's Method .. and slide this window along the entire EEG signal.
Outcome is a PSD plot <-- which gives intensity of each frequency band
Thanks!


texts:
I am developing a user group algorithm in NOMA, but I need channel gain for that.
I know that the channel gain consists of some components, but I have not found a way to calculate the gain.
Does anybody know?


texts: Confusion regarding calculation of fundamental period?

texts:
I am reading signal processing first and in chapter 3 ex3.8 i came across an example of fundamental period as shown  in attached photo

It apparently shows that signal $$x(t)=\cos^2(4\pi t)$$ has period 0.5 but then it also writes that fundamental period is 0.25

How is he doing that?

Also will be fundamental period if $$x(t)=\cos^n(4\pi t)$$ where n can be 3 or 4 or 5


texts: Fast/efficient way to compute Laplacian edge enhancement filter

texts:
I would like to implement a somewhat smarter Laplacian edge enhancement convolution. Right now it is implemented as (generic 3x3 convolution):

sum = 0
for k = -1 to 1 do
   for j = -1 to 1 do
      sum = sum + h(j +1, k + 1)*f(x - j, y - k)
   end for
end for
g(x, y) = sum

Looking for way to improve to implementation all I could find was separable filter optimisations (as in). However Laplacian edge enhancement is not separable:

octave:1> m = [-1 -1 -1 ; -1 8 -1 ; -1 -1 -1];
octave:2> rank(m)
ans =  2

Is there any litterature I could be missing ? When looking at the coefficients it feels as if somewhat could be done for matrix of type (CENTER=8, EDGE=-1):

E  E  E
E  C  E
E  E  E


texts: Library to differentiate people by their voice timbre

texts: Say there are 3 people in the room. Each of them have a unique voice timbre. I want to recognize people by their timbre and know who's speaking.


texts: Initialize a Room Impulse Response using reverberation time(T60)

texts:
I am doing Speech dereverberation using Non-negative matrix factorization. I need to initialize an impulse response($H$) such that its reverberation time ($t_{60}$) is $300$ ms. Let the length of $H$ be $10$.

I thought of doing that by identically initializing each row of $H$ using a linearly decaying envelope. Somewhat like this

x=1:10;
h = exp(-x/2);
H = repmat(h,600,1);

This will give me a $H$ of dimension $600$ X $10$. But, I don't understand how to use $t_{60}$ for the initialization in MATLAB. Please tell me.


texts: Basic Fourier Series Not Understanding how to interpret n values

texts:
I have a question regarding $n$ in the Fourier series. If a question states "find the Fourier series for (any function) and find the values of an $b_n$ etc", say I'm finding it from the trig way when $n<=4$. What $n$ values do I use, is it 1,2,3,4 and if so why does it start at 1 and not include 0. Also why doesn't the $n$ values keep going past negative so e.g 4 3 2 1 0 -1 .... What makes it start at 1 for this value of n.

Is it since in the formula n=1,2,3.....?

However if it stated for $-3<=N<=3$, would the n values be -3 -2 -1 0 1 2 3? if so why? Also would the zero be in it or not? My $n$ values were less than 0 would you then have to use the complex way?


texts: Obtaining normalized matrix for the Haar Wavelet Transform

texts:
I've been reading this article: http://aix1.uottawa.ca/~jkhoury/haar.htm which explains the Haar Wavelet Transform.

At a certain point, the author says:

  ... Since the transformation matrix $W$ is the product of three other matrices, one can normalize $W$ by normalizing each of the three matrices. The normalized version of $W$ is ...

and provides final matrix.

But he doesn't provide the exact procedure for matrix normalization he used. Probably someone with a deep understanding of the Haar Wavelet Transform can explain me, how matrices were normalized so that final result was obtained.

I've been searching the Internet for matrix normalization and couldn't find anything suitable. Also contacted author of article but didn't get any response so far.


texts: Smooth 2D Data with Discontinuous and Artificial Jumps

texts:
I have some data $(X,Y,Z)$, which is a set of measurements $Z$ over a $2D$ space $X$,$Y$. The $Z$ data on this space is continuous except for some discontinuous jumps in certain
domains of $X$,$Y$ (due to calibration errors).

For example, let's say $X$ and $Y$ are both between $[0,1)$, and let's say my underlying true function is $Z=f(X,Y)$. For most of the $X$,$Y$ space $Z$ will equal my true function, but perhaps for $0.25 < X < 0.5$, $0.25 < Y < 0.5$, my data will show $Z=f(X,Y)+4$, an offset of $4$ from the true function.

I can brute-force look for discontinuous jumps and then correct for them in the data set (by setting the edge points in the error region equal to their nearest points in the non-error region), but this can be problematic for certain functions and certain domains that I'm working with. I'm wondering if there's a signal processing method that is designed to work for these types of problems?

Edit: To be clear, I only need to learn the function $f(X,Y)$ up to a constant, so it doesn't matter what part of the space I consider the "error-free" region or not. I simply can't have the function exhibiting different offsets in different regions of the 2D $X,Y$ space.


texts: Attenuating peaks of a signal in a specific direction

texts:
Consider a segment of a time-domain signal shown in the following graph and the two marked peaks, the first one facing downwards (in red) and the second one facing upwards (in black).

I would like to attenuate the peak facing downwards only. Unfortunately, both peaks seem to have the same frequency, so a linear filter that attenuates that frequency will attenuate both peaks.
Is there a kind of a filter (perhaps a non-linear filter) that can attenuate peaks that face downward and keep peaks that faces upward (or vice-versa)?


texts: MATLAB: get members of histogram populations EXCLUDING those inbetween

texts:
The Problem

I have a (noisy) square wave signal and need to extract the data points that make up the 'baseline' and 'plateau' ONLY. I do not want data from the rise and fall. So far I have accomplished the following:

data points designated as 'baseline' or 'plateau' superimposed onto original signal

I have defined 'baseline' data points as those from the minimum to the mean of the lowest 25% of signals + 1 standard deviation of the lowest 25% of signals.

Similarly I have defined 'plateau' data points as those from the maximum to the mean of the highest 25% of signals - 1 standard deviation of the highest 25% of signals.

This works ok, but as you can see: sometimes 1 standard deviation is not enough to encapsulate all of the baseline noise. 2 Standard deviations often encapsulate too much of the rise/fall of the signal.

Can anyone provide a more elegant solution to this?

My aim is to measure the signal to noise ratio of the square pulse, so I need to extract the populations that will give me baseline noise and the mean pulse height.


texts: How to compare multiple polynomials in MATLAB?

texts:
I am trying to fit a second order polynomial to various curves on a binary image using MATLAB polyfit function.There are many such curves on the image, and the coordinates of the pixels forming the curve are known.

How do I compare the co-efficients, so as to arrive at one approximation to all such curves?

Here are some co-efficients of such curves:

 0.0010126     -1.5981       1074.8
 0.0024871     -3.3327       1866.3
 0.002662      -3.7067       1631.4
 0.0025574     -4.5129       2148.4
 0.002114      -3.3457       2195.6
 0.0015383     -2.6023       1682.9


texts: Raised-Cosine Filter - phase response graphs for three kinds of filter

texts:
In the previous post it has been mentioned that the phase response of Raised-Cosine filter is linear and is given with ϕ(ω)=−ωT, where T is the delay. But, I have also read in literature that it is also a minimum phase filter and non-linear. Which one is correct, and if it is minimum phase how I can calculate it?

I have three kinds of Raised-Cosine filters, as per picture below. What is the phase response for each kind of filter? How the phase response graph will look like for each filter?


texts: Limiter, compressor, AGC, dynamics processor?

texts: What are the primary differences between the usage or definitions of the terms “AGC”, “compressor”, “limiter”, and dynamics processor? They all seem to be forms of playing with gain knobs.


texts: Difference between Contrast and Intensity of an Image

texts:
I wanted to generate Gray Scale wedge image of 10 Levels in MATLAB and then increase and decrease its Intensity. By high intensity, I mean mapping increasing the level of gray scale intensity in an image from lower value to higher values.

The examples I have tried were performed on a Gray Scale wedge image of 10 levels:

0   28   57   85  113  142  170  198  227  255

I was using the function imadjust. For increasing the intensity levels I used:

 imadjust(grayStepImage, [.35 .75], [0.7 0.8])

The output:

179  179  179  179  184  192  199  204  204  204

and the result in image form:

Similarly, for lowering the intensity I used

imadjust(grayStepImage, [.35 .75], [0.1 0.39])

Original and modified levels are:

0   28   57   85  113  142  170  198  227  255

26   26   26   26   43   64   84   99   99   99

and the result in image form:

Was I actually increasing/lowering intensity of the image? How do those operations relate to the case of increasing and decreasing the contrast of an image?

What I know about contrast is,

  Contrast is defined as the separation between the darkest and
  brightest areas of the image. Increase contrast and you increase the
  separation between dark and bright, making shadows darker and
  highlights brighter.

Does that mean, if I want to increase the contrast of an image I should increase the high intensity levels of an image to much higher and lower intensity level to much lower?

I am confused with contrast and the intensity as they seem opposite to each other.
Kindly help me out.


texts: Minimum number of points to specify a sine wave with a given accuracy

texts:
I am generating a sine wave from a DAC at a constant frequency. If each full cycle is specified by N data points what is the relationship between the accuracy of the rendered sine wave and N?

N=1 will obviously be pulses. N=4 (peaks and zero crossings) will yield a triangular wave of sorts. Is there some equation that shows how good a sine wave approximation is being created? [Obviously I can put it through a post DAC filter to smooth it etc]

The kind of answer I am looking for is: "Specifying a sine wave with N points yields 1:X accuracy"


texts: Why wavelets based transmitter/receiver digital signal processing aren't common?

texts:
I have seen this thread: Difference between Fourier transform and Wavelets

AFAIK there is no common usage of wavelets in the real-time DSP world (excluding image and video processing).

I am curious why is it so?

Where on the other hand Fourier analysis is widely common and all the transmitter/receiver block diagrams that I have seen contains an FFT module in it:

Is it related to runtime complexity issues? I guess not, since Wavelets transform can be formalized as convolution, right?

e.g aren't wavelets a good candidate to replace the Short-time Fourier transform method ?


texts: Is there any difference in Implementing sigma delta modulator using filters and state space model in FPGA?

texts:
Sigma delta modulation is extensively used in quantization to reduce quantization noise.

In the literature one can see different architecture for example python-deltasigma to implement a modulator. One can learn the output of the modulator by implementing the filters or one can convert the problem into a state space model and design the modulator based on A,B,C,D parameters of the state space model. The question: is there a difference or advantage/disadvantage in working with state space model instead of implementing the filters as per the  architecture, specifically for FPGA implementation.


texts:
I would like to use time averaging technique for denoising vibration signal, using the below function, how do we choose the appropriate parameters D and N for optimal denoising  !
% sigav.m - signal averaging

%
% y = sigav(D, N, x)
%
% D = length of each period
% N = number of periods
% x = row vector of length at least ND (doesn't check it)
% y = length-D row vector containing the averaged period
% It averages the first N blocks in x
function y = sigav(D, N, x)
  y = 0;
  for i=0:N-1,
      y = y + x((i*D+1) : (i+1)*D);       % accumulate i-th period
  end
  y = y / N;


texts: How to measure the short-time change in an audio signal?

texts:
The choice of frame length for calculating the Mel Frequency Cepstral Coefficients of a signal, depend on the how much the signal change over a small time scale.
I guess that human speech change rapidity while some animal sound change more slowly, but what is the measure of small time scale of an audio signal?


texts:
In LTE systems the in the scrambling section, a gold sequence is used. But I've studied that gold sequences are used in spread spectrum communication. This generally increases the bandwidth. So doesn't this mean that scrambling results in an increase of bandwidth.?

regards,
phani tej


texts: Controlling one signal with another one (Rocket engine effect - Subtractive synthesis of sound)

texts:
I am programming some sound effects in Java and exporting them into .wav files. Currently, I am trying to program a rocket engine sound effect. I want to do it in the following way:

The sound of a rocket engine may be synthesized with a red noise generator controlled by a second red noise generator. The parameter of the first generator modified by the second one is the number of interpolated samples, influencing the spectral content of the generated noise. In order to simulate changes of sound intensity (e.g. during launch) the envelope generator should be used.

I am wondering how can it be done, e.g. what does it mean that one signal controls another one. Probably this part is explaining it, but I am not sure what to do now:

  The parameter of the first generator modified by the second one is the
  number of interpolated samples, influencing the spectral content of
  the generated noise.

Is it about this parameter describing how many values from white noise are taken and linearly interpolated while creating the red noise? (see my simple drawing explaining this process below)

I have a red noise generator, which returns an array of doubles with values between -1 and 1 (it is generated from the white noise as described). What am I supposed to do now? How can I control the second red noise? I guess that it does not mean that I should control the amplitude of the second signal. Does it?  Schema of steps required to obtain the rocket engine sound effect is attached below.


texts: realtime software 90 degree phase shift a sinusoidal signal

texts:
I would post this to stack overflow but it seems to me more like a signals and systems problem, which EE's thrive on.

I need a realtime system implemented in software that can generate a sinusoid that matches the input signal's frequency but with a phase lag of 90 degrees. The input signal will be a pure tone but can drift a bit in frequency (0.4 to 0.8hz) over time. Eventually I'll want output amplitude to scale with input amplitude as well but that's for another day. I'm aiming to be able to respond to input signal changes such that I'm under 5% phase error within about 10-15 seconds (a few input signal cycles).

My first attempt at this has used a hilbert transformer FIR filter to give me the "instantanous" relative phase between my input signal (100hz sample rate) and output, then a PI controller with set point =-90 to drive frequency f_out of output_signal = math.sin(2*math.pi*f_out*t). My hilbert transofmer FIR has a latency of 2.5 seconds and the dead time seems to prevent me from having any kind of a P constant without instability.

I'm currently doing a deep dive into phase locked loops (I'm new to them) but it's not clear to me yet if this is the right tree to bark up.

I'd appreciate any manner of direction/advice.


texts: What's the noise in this signal ? (Beginner question)

texts:
There is this raw ECG signal that I have obtained its frequency content in Matlab. The signal was sampled at 600Hz. The following is signal and the frequency content.  I want to know what the interferences are. I'm very new to DSP, and I feel like I'm lacking some common sense to understand a signal by looking at it.
It looks like to me the isolated peaks are the noise. How can I tell which one the noise is?


texts:
I'm trying to solve the following problem. I got a "low frequency" input 2D signal over a square region. I'll collect a few samples, somewhere around 10-30 maybe, the exact sample count will be dynamic. Then I want to compress this to much less data and later reconstruct the square representation. The compression will only occur once but the decompression multiple times. Performance is essential. So I need a method that compress fast, but decompress faster and to a compact format so the data transfer won't become the bottleneck.

The actual application: It's for real-time computer graphics processed on a GPU. I compute samples of lighting in a vertex shader. Then the data transfers to a pixel shader where it gets reconstructed. I will accept a certain amount of error as long as the result is "blurry", i.e. no sharp artifacts. Also transcendental functions are somewhat slow on a GPU (even though a slim transcendental representation will beat any verbose solution ofc)

I've been thinking about several methods. Maybe something fourier based. I could use a standard fourier series in its power series form, e^i(mx+ny) but that would require n^2-1 params to be handled right? Does it exist some basis function that's less verbose in 2D? Or in general more suited for this kind of thing?

Maybe I could do something related to curve fitting (plane fitting). Or maybe something completely different altogether. Anyone with a suggestion?

Thanks!


texts:
I'm doing work on SDR for GNSS and I'm tracking the signals I/Q using PLL.
The calculations for the PLL is done like the following for error estimation:

double PLLAtan::compute(const CorrelationResult &correlation_result)
{
    sampleType I_P = correlation_result.getPromptCorrelator().real();
    sampleType Q_P = correlation_result.getPromptCorrelator().imag();

    double discriminator = static_cast<double>(std::atan(Q_P / I_P)) / (2 * M_PI);
    return discriminator;
}

in each epoch time frame:

pll_fll_filter_.setDiscriminator(pll_disc_->compute(correlation_result));

My problem is the PLL filter is oscillating between 0.98 to 0.96 and it never locks to 1

Should I do something like accumalting the error in the correllation result from the I/Q calculate the error between last I/Qs and current and use it as input for PLL phase error correction ?

    double LoopFilterPLLFLL::updateLoop()
    {
        if(filter_order_ == 2)
        {
            double out_velocity = velocity_integrator_.setInput(w0f_ * fll_discriminator_ + w0_ * w0_ * discriminator_);
            loop_output_        = a2_ * w0_ * discriminator_ + out_velocity;
        }
        if(filter_order_ == 3)
        {
            double out_acceleration
                = acceleration_integrator_.setInput(w0f_ * w0f_ * fll_discriminator_ + w0_ * w0_ * w0_ * discriminator_);
            double out_velocity = velocity_integrator_.setInput(a2_ * w0f_ * fll_discriminator_ + a3_ * w0_ * w0_ * discriminator_
                                                                + out_acceleration);
            loop_output_        = b3_ * w0_ * discriminator_ + out_velocity;
        }
        return (loop_output_);
    }

void LoopFilterPLLFLL::setLoopParameters()
{
    LoopFilter::setLoopParameters(fp_.filter_order_pllfll, fp_.update_rate, fp_.narrow_noise_bandwidth_pll);

    calcLoopParameter(0.0, fp_.wide_noise_bandwidth_fll);
}

void LoopFilterPLLFLL::calcLoopParameter(double pll_noise_bandwidth, double fll_noise_bandwidth)
{
    if(filter_order_ == 2)
    {
        w0f_ = fll_noise_bandwidth / 0.25;
        w0_  = pll_noise_bandwidth / 0.53;
    }

    if(filter_order_ == 3)
    {
        w0f_ = fll_noise_bandwidth / 0.53;
        w0_  = pll_noise_bandwidth / 0.7845;
    }
}

Gain Parameters:

 "pll": {
                "filter-order": 2,
                "discriminator": "PLL-ATAN",
                "wide-bandwidth": 20,
                "narrow-bandwidth": 5
            },
            "fll": {
                "discriminator": "FLL-DP",
                "wide-bandwidth": 10,
                "narrow-bandwidth": 2
            }


texts:
I have to detect my elbow. I'm able to produce an image in which I'm facing the camera, standing before it and my arms are in my body plane(body plane is plane containing my head,chest,and legs)(i.e. to say that my arms don't stick out of this plane though they may be placed on my chest parallel to the body plane). Then using skin detection I'm getting my both arm, hand and palm.(I'm wearing half shirt so little portion of my arm is detected between elbow and shoulder). Then using convex hull defect method the point giving the largest defect should be my elbow.

This thing is fine theoretically but does not work practically. Before using the convex hull and defects method my hand and arm contour is a bit wavy (and if the waves become a bit more high then that point is detected as my elbow than the real one). Or if i use polynomial approximation method then also it is not detecting.

Skeletonization of image takes considerable time and I'm doing this real time.(My image has two contours of my both hands, that's it).

I'm doing this in opencv c++. Any Help is appreciated.


texts: Which method should be selected to enhance the compressed image?

texts: I have a raw infrared image with 16384(2^14) intensity levels.The raw image must be compressed into one with 256(2^8) intensity levels. In order to retain the original details( small blocks) , the compressed image must be enhanced. Which the best method I should select?Someone advice that "histogram matching" method should be taken.But the both images are different in intensity levels,how can do it?Any advice will be appreciated!


texts: Deconvolution (Linear Convolution) with an Under Determined System of Equations?

texts:
If I have a measured signal $\mathbf{y}$ which is the result of the true signal $\mathbf{x}$ convolved with another function (linear and not circular convolution), I always seem to get an underdetermined system. Convolution expressed as a matrix equation:

\begin{equation}
\mathbf{y} = \mathbf{A} \mathbf{x}
\end{equation}

I am interested in finding $\mathbf{x}$, given $M=5$ observations and knowing that the convolution matrix is defined by a boxcar function of window length $n=3$. The corresponding $M\times N$ convolution matrix is defined as the following, where $N=M+n-1$.
\begin{equation}
\mathbf{A} = \begin{pmatrix}
1 &   1 &   1 &   0 &   0 &   0 &   0 \\
0 &   1 &   1 &   1 &   0 &   0 &   0 \\
0 &   0 &   1 &   1 &   1 &   0 &   0 \\
0 &   0 &   0 &   1 &   1 &   1 &   0 \\
0 &   0 &   0 &   0 &   1 &   1 &   1 \\
\end{pmatrix}
\end{equation}

In the literature, the convolution matrix is often $N\times N$ (or $M\times{}M$, however you want to define it). But it seems that with a finite number of observations, the convolution matrix corresponding to the available observations would have to be $M \times (M+n-1)$, which is always underdetermined (except when $n \equiv 1$, in which you don't have convolution of your signal).

In most convolution/deconvolution examples, the matrix is often represented as being "square," in which there exist rows in the matrix where the number of non-zero elements are less than $n$. This does not correspond to any possible measurement. Padding is common for the input signal ($\mathbf{x}$) but not for the output ($\mathbf{y}$), so how can the system not be underdetermined?

Edit for clarification:

I understand that convolution matrix is traditionally written as $N+m+1$:

\begin{equation}
  \mathbf{A} = \begin{pmatrix}
    1 & 0 & 0 & 0 & 0 & 0 & 0\\
    1 & 1 & 0 & 0 & 0 & 0 & 0\\
    1 & 1 & 1 & 0 & 0 & 0 & 0\\
    0 & 1 & 1 & 1 & 0 & 0 & 0\\
    0 & 0 & 1 & 1 & 1 & 0 & 0\\
    0 & 0 & 0 & 1 & 1 & 1 & 0\\
    0 & 0 & 0 & 0 & 1 & 1 & 1\\
    0 & 0 & 0 & 0 & 0 & 1 & 1\\
    0 & 0 & 0 & 0 & 0 & 0 & 1
  \end{pmatrix}
\end{equation}

In practice, you never measure the signals corresponding to the first two and last two rows, so it would seem that you cannot use this matrix for creating the deconvolution filter.


texts: what is the kernel ker in the ImageCorrelate of Mathematica Program?

texts:
Thank you for consideration. In image Processing with Mathematica Program, it is an order to ImageCorrelate for an image. In this order there is a "kernel ker", I don't know what is that exactly?
Thank you for help


texts: How to Normalize the ideal Tx and Equalized Waveforms in LMS algorithm

texts:
Below mentioned steps has been done but I don't know if I am doing the right way of Normalizing the 2 signals.
Ref_1 = load('Tx_16QAM_25_New.txt'); % transmitted data
Tx(:,1)=Ref_1(:,1)+1i*Ref_1(:,2);
Tx = Tx(:,1);
Tx1=rms(Tx)
normTx=Tx./Tx1

Eq=c'*X;   % received symbol estimation
Eq1=rms(Eq)
NormEq=Eq./Eq1

Normalized Tx,Rx plots are attched while no plot is obtained for the normalized equalized constellation.
Rx1 =3.1859, Eq1=Inf, Tx1=3.16 (RMS Values)


texts: How to select the ARMA model parameters?

texts:
I have a series of data containing 120,000 points. The mean of each N(=60) point is zero. I want to forecast the next 60  points using the ARMA model. My question is, specificaly, how to choose the appropriate model parameters (p and q) and the size of the training set.
Particularly, I run the ADF Test while providing the 100,000 points and getting the following results:
ADF Statistic: -52.49
p-value: 0.0
Critial Values:
   1%, -3.43
   5%, -2.86
   10%,-2.56

You can see the ACF and PACF and the absolute value of the FFT of the empirical ACF as follows. It seems that almost all lags are critical.

Then, I used a small subset of the data set (400 points) and applied the ADF test, and plotted ACF and PACF, the results are as follow:
The last figure shows the absolute value of the FFT of the empirical ACF as follows.
ADF Statistic: -6.21
p-value: 5.38e-08
Critial Values:
   1%, -3.44
   5%, -2.86
   10%,-2.57


texts: How to design a narrow band pass filter without DSP toolbox

texts:
I've two signals X and Y. These are proximity sensors reading measuring vibrations of a rotating shaft. Now I want to plot orbit plots for these signals. I can plot unfiltered orbits with just X and Y. But since there is lot of noise in the signal I filtered it. Now there are dominant 1X frequency component and may be other sub frequency components. I want to plot the filtered orbit of 1X say 656 Hz. The sampling frequency I used is 1000000 Hz. I know band pass filter could give me only the 1X component. But I don't have DSP toolbox. Could someone help me how to design and implement narrow band pass filter without toolbox.
Any suggestion or hints or idea would be helpful.
Thanks in advance.


texts: How to Extract a Period of a Periodic Pulsed Signal?

texts:
I rercorded a sequence of zeros and ones which is repeated many times. However, the first and the last repetitions may not be complete because I have to start recording at a random time, so the sequence can be started when I hit record and the sequence may have not finished when I stop recording.

I'd like to know how to extract, for instance, the 3rd period of the complete recording, because this is something I know it would have the complete sequence.

here is a example of the signal I'm working with. I computed the FFT and I think the frequency is 96000Hz (half of the sampling rate), but still I'm not sure how to extract a complete period from it.

How could I do this?


texts:
I have access to the step response of a system and I want to find its poles and zeros without knowing the order of the system.

Consider an example of a step response of the system shown in the figure below.  Here is the link to the raw data if someone wants to try it (https://pastebin.com/BpT4d9RW).

I want to build an IIR filter (with arbitrary number of poles and zeros) that has the same step response (or pretty close to it).

I have tried using the prny() function in MATLAB for that purpose, and the closest result I could get was the following. I am actually not sure what I am doing wrong here and if there are any considerations for using the prony() function that I am not aware of.

I would appreciate if someone can help me on correct use of the function OR can come up with some other approach that works.

Please also note that the approach should be general enough to accommodate any type of step responses (over-damped or under-damped).


texts:
Here I have a MATLAB code for Ideal Low pass filtering an Image. In the below code they have used D=sqrt(U.^2+V.^2); which finds some distances. Then this distance matrix D is compared with the cut-off frequency to create a filter function H ( H=double(D<=P); ) that has be multiplied with the fourier transformed image F. i.e, G=H.*F; to get the output image.

My questions: What is this distance ? Why we need this distance matrix D ? What information this gives to us? How it helps in Low Pass Filter an Image? (I understood that it helps to create the filter function H) but Why we need the distances to create the filter function? I need to it understand clearly. This may be a simple question but could anyone please explain me this with little patience.

Thank you for reading my question.

Reference: MATHWORKS:IDEAL LOW PASS FILTER

%IDEAL LOW-PASS FILTER

    function idealfilter(X,P) % X is the input image and P is the cut-off freq
    f=imread(X);  % reading an image X
    [M,N]=size(f); % Saving the the rows of X in M and columns in N
    F=fft2(double(f)); % Taking Fourier transform to the input image

    u=0:(M-1);
    v=0:(N-1);
    idx=find(u>M/2);
    u(idx)=u(idx)-M;
    idy=find(v>N/2);
    v(idy)=v(idy)-N;
    [V,U]=meshgrid(v,u);
    D=sqrt(U.^2+V.^2);   % finding out the distance

   H=double(D<=P);       % Comparing with the cut-off frequency
    G=H.*F;               % Multiplying the Fourier transformed image with H
    g=real(ifft2(double(G))); % Inverse Fourier transform
    imshow(f),figure,imshow(g,[ ]); % Displaying input and output image
    end


texts:
Actually I am doing an RF measurement and this includes phase and power. i generate a list of powers ranging for example from (-50 to 40) dBm, while the phase rotate to make it for example spiral shape as shown in the attached graph. Then I measure that spiral.

I wanted to calculate the error between the measured and targeted signals.

The normalised mean square error (NMSE) has been suggested to be used. But This one shows me high error even when i see virtually a good agreement through the graph.

The Equations that I have used is shown below in the attached graph. is there better method to calculate such kind of error?


texts:
How can I calculate 8 point DFT of signals whose length is less then 8 ( say 2,4) then what will I assume other members in formulae "0" or rotation or keep repeating the same two numbers

Consider I want to find 4 point DFT of [1,2]


texts:
So, I've created a simple sound analysis application and one of the features I've implemented is the spectral centroid (as explained here https://en.wikipedia.org/wiki/Spectral_centroid). In order to get reliable results, we need to set a threshold on our spectral centroid data because the FFT always returns some energy in bins that really don't have any energy.

Now, I would like to be able to manipulate the spectral centroid (e.g push the values up or down, which to some extent resembles pitch shifting) and then use an inverse FFT operation to convert the values to waveform data. This way, we'd be able to hear the changes but I'm not sure how to do this.

The spectral centroid is calculated by using magnitudes, not raw FFT data. The magnitudes are calculated as given below

mag[i] = 10*log10(sqrt((pFFTLReal[i] * pFFTLReal[i] + pFFTLImg[i] * pFFTLImg[i])/fftNorm));//fftNorm is basically a correction factor that depends on a type of window (blackman, hamming etc)

So, given all this, how would I change the spectral centroid by, say, 10 percent and see that change in the generated (manipulated) waveform?

Is this even possible?

Please, if this is not clear enough, I'll try to provide some more information.

Thank you


texts:
I have a stream of audio blocks at input, each having 64 samples.

Sending a previous block of samples to output, i.e. delaying one block, means applying a unity gain FIR filter?

out[n] = 0 * in[n -  0] +
         0 * in[n -  1] +
         0 * in[n -  2] +
         ...
         1 * in[n - 64].

What's the frequency response of this filter?
Is there a use case for such design, or does this make no sense?


texts: First order hold from zero order hold filter

texts:
I'm stuck in this question for more than 4 hours. (h is the sampling period)

I simplified the block diagram but it just seems impossible to get the FOH transfer function from the reduced block diagram... Am I missing something?
My try:


texts:
I am trying to understand the functionality of fft.m in matlab .

For example: if I input a signal of length $7$, and use the FFT function, the output is of the same length.

How does that occur if the radix-$2$ takes the nearest power of $2$, which is $8$?

So how does it return an output of length $7$, assuming it operates the way I described?


texts: What do high and low order have a meaning in FIR filter?

texts:
I'm new DSP fields, and slowly study of FIR Filer and with Matlab which is fir1 function.

My question is that

What do high and low order have a meaning in fir1 filter? I can't understand what does it have a meaning when N=2 or N=48.
In cut-off frequency, I found some explanation from the MATLAB as the below,

fir1   FIR filter design using the window method.
    B = fir1(N,Wn) designs an N'th order lowpass FIR digital filter
    and returns the filter coefficients in length N+1 vector B.
    The cut-off frequency Wn must be between 0 < Wn < 1.0, with 1.0
    corresponding to half the sample rate.  The filter B is real and
    has linear phase.  The normalized gain of the filter at Wn is
    -6 dB.

B = fir1(N,Wn,'high') designs an N'th order highpass filter.
You can also use B = fir1(N,Wn,'low') to design a lowpass filter.

If Wn is a two-element vector, Wn = [W1 W2], fir1 returns an
order N bandpass filter with passband  W1 < W < W2. You can
also specify B = fir1(N,Wn,'bandpass').  If Wn = [W1 W2],
B = fir1(N,Wn,'stop') will design a bandstop filter.

But  I want to know should we have to know the sample rate to set cut-off frequency? usually in practically, How does it handle to set the cut-off frequency?

I've seen from the book, in the Liner system, if we put the impulse to the Liner system then we can totally analyze the system.
In practically, how to implement the impulse signal to the liner system to get the filter information?


texts:
I have an audio buffer with 1024 samples per second.
The sample rate is 44100.
I am trying to identify the frequency tones.
Using the Gorzel algorithm, I can get the frequency but it's not that accurate.
I was trying some other scope program and I got better results with the same hardware.
How can I decrease the id-bandwidth (to be sharper) so the id will be accurate for the target frequency?
(I can't use FFT because it's real time and takes some memory.)

float goertzel_mag(int16_t* data ,int SAMPLING_RATE ,double TARGET_FREQUENCY,int numSamples )
{
    int     k,i;
    float   floatnumSamples;
    float   omega,sine,cosine,coeff,q0,q1,q2,magnitude,real,imag;
    float   scalingFactor = numSamples / 2; // -2

    floatnumSamples = (float) numSamples;
    k = (int) (0.5 + ((floatnumSamples * TARGET_FREQUENCY) / SAMPLING_RATE));
    omega = (2.0 * M_PI * k) / floatnumSamples;
    sine = sin(omega);
    cosine = cos(omega);
    coeff = 2.0 * cosine;
    q0=0;
    q1=0;
    q2=0;

    for(i=0; i<numSamples; i++)
    {
        q0 = coeff * q1 - q2 + data[i];
        q2 = q1;
        q1 = q0;
    }

    real = (q1 - q2 * cosine) / scalingFactor;
    imag = (q2 * sine) / scalingFactor;

    //double theta = atan2 ( imag, real); //PHASE
    magnitude = sqrtf(real*real + imag*imag);
    return magnitude;
}


texts: Need of $\tt abs()$ method when plotting a power spectral density for a given dataset

texts:
I am a newbie in signal processing and would like to know the significance of using the abs() function and squaring the values received as an output of fft() function of Python's scipy.fftpack library, used when trying to plot a power spectral density for a dataset. I have found that many of code examples to plot a power spectral density do use an abs() and then square the values obtained thereafter.

Can anyone please provide me a reason for doing so?
Can't we just directly plot the values obtained from fft() function in Python's scipy.fftpack library?

Here is the code I have written till now to plot a power spectral density by referring some of the code examples,

import scipy.io as sio
import numpy as np
Import matplotlib.pyplot as plt
import pandas as pd

df = pd.read_csv("denoised.csv")
data = df.values

from scipy.fftpack import fft,fftfreq
dft= fft(data)
PSD = np.abs(dft) ** 2


texts: How do I estimate the derivative of the current position, when I have only values from past to present?

texts: If I have a discrete real-time signal $x[n]$, with its latest value $x[i]$ and all its past values $x[i-t]$, how can I estimate the derivative at $x[i]$?


texts: Smoothing a Signal which contain White Noise

texts:
I'm trying to smooth a signal in order to reduce white noise as a result of changes in amplitude across sequential frames of an audio signal.

I've been investigating using a moving average or boxcar filter to help minimise the issue.  In my investigations it would appear normal implementations of a boxcar filter only smooth the contents of a single frame/buffer of audio.

Can a boxcar filter be implemented in order to smooth the transition across frames as well? or is there a better approach?


texts: How can I know if these two signal are orthogonal to each other

texts:
So the question is to determine if this statement is true or false. The question can be seen below.

I know that two signals are orthogonal if their inner product is equal to 0. I know how this calculation is done with vectors but not with signals. Would be really helpful i someone could explain this to me. Is it also possible to just directly see it in the two graphs if they are orthogonal to each other?


texts: Finding the minimum phase h[n] and its Z transform

texts:
Hello, this is one of my homework questions and i have already solved the first question but im having trouble gettin a relation that helps me solve the second one. From the question i understand that R(z)=H(z).H(1/z) but i dont know how to develop this.
This is my response to the 1st question

Thanks


texts: I'm currently trying to understand decimation and possible ways to retrieve back the original signal. If the signal of supposed frequency $X$ Hz was decimated by a factor of four, would there be $4$ distinct peaks in the decimated signal spectrum at frequencies $X/4$, $X/2$, $3X/4$, and $X$? Furthermore, would I need to upsample (interpolate between the $2$ data points) the signal by $4$ times to get back the original signal $X$ (Would the upsampled signal of the decimated signal spectrum look similar to the spectrum of $x$)?


texts:
A little background: I am simulating the response of an antenna to a pulse. The pulse has a wide frequency range, but the antenna only responds to a fairly narrow frequency. I can extract voltage and current from the antenna: they look like sinusoidal signals quickly decaying in time.

I now want to calculate the power spectral density (PSD) of that received signal. If I calculate the FFT of my voltage and current signals I see a nice isolated peak at the resonance frequency, so I would expect the power to also be concentrated around that frequency.

My first attempt was simply to multiply the voltage and current in time domain and then calculating the PSD of that. It did not work, the peak of the PSD happens at a frequency roughly double the expected one and at the desired frequency there is actually a "notch". I suspect that it is because by simply multiplying current and voltage I have effectively modulated one onto the other, therefore mixing their frequencies.

I can get a PSD that seem to make sense if instead I multiply the voltage and current in frequency domain and then calculate the PSD from the resulting spectrum.

So my questions are:

Why can't I just calculate power as v*i in time domain?
Does multiplying voltage and current spectra to then calculate the PSD make sense?


texts: Combining channels into single FFT left/right vs mid/side

texts:
I am creating a spectrum analyser to visualize some stereo audio data in real time. I want to display a single spectrum for the input and a single spectrum for the output. The stereo data is the usual two channels: left (L) and right (R).

I am combining them into a single channel by taking the FFT of L + R. Is this conventional or would it be more appropriate to combine FFT(L) + FFT(R) rather than use FFT(L + R).

Also if I encode the left and right channels as mid (M) and side (S) do some processing on each of M and S and convert back to L and R then I observe the following:

If the side channel is modified, then taking FFT(L + R) shows no visible difference between the input and output (because L + R is pretty much the mid channel without side) so the question is: is there a way to combine the mid and side channels together so that taking the FFT of them gives a good representation of the sound. Or should I be doing FFT(L) + FFT(R) for both L/R and M/S processing modes? I would rather not have to do two FFTs if one would suffice.


texts:
I want to do one project and I just want to make sure that I've understood correctly.

Human's voice frequency is around 85-255Hz.

When we talk on phone, carrier is modulated with voice signal in order to shift frequency spectrum to higher frequency for transmission. If we take for example GSM900 which is 876 - 915 MHz paired with 921 - 960 MHz. That means that our speech will be in this GSM900 range after modulation?

Moreover, how much bandwidth user gets in 876 - 915 MHz paired with 921 - 960 MHz range? Any articles about that?


texts: What are the eigenvalues of the 8 point DFT matrix?

texts: I  know that the eigenvalues for 4 point DFT matrix can be found from $F_4^4=I$. Is this also valid for 8, 16 and higher orders? For example with  8 points, will it be $F_8^8=I$ ? If not, how  can I can compute them?


texts: Dynamics processor for speech - apply gain on delayed samples?

texts:
A Dynamics processor (limiter, compressor, expander, and gate) is calculating the output gain at some point in time, and it is applying this gain at some delayed point in time.

This allows for Dynamics processor to "look ahead" and act on some audio events a few mili-seconds before they actually arrive. This is usefull in cases where sharp, transient signals occur. However, my use case is speech, in a bit noisy environment, but nothing "sharp", no transients, no drum kicks and similar.

Would adding lookahead in such a "speech" use case be beneficial?


texts: Is my Kalman Filter model reasonable? Using for a 3 wheels rover platform

texts:
Background: I'm building a 3 omniwheel rover platform that looks something like this:

It has 1 IMU sensors on each of the wheels (3 in total). So in theory, I can get gyroscope and accelerometer data from each of the wheels and use that information to calculate the current speed of the wheels.

So let's say this is my model:

$x_k = \begin{bmatrix} v_k\\a_k \end{bmatrix} = Ax_{k-1} + Bu_{k-1}$

where: $A = \begin{bmatrix} 0&\Delta t \\ 0&0 \end{bmatrix}$ and $B = \begin{bmatrix} 1&0\\0&0 \end{bmatrix}$ and $ u_{k-1} = \begin{bmatrix} u_{velocity}\\0 \end{bmatrix}$

I can "control" the speed of each wheels, but of course that's not accurate. But I cannot control the acceleration of each wheels, so I let $u_{acceleration}$ in matrix $u_{k-1}$ to be 0.

$z_k = \begin{bmatrix} v_{measure}\\a_{measure} \end{bmatrix} = H\begin{bmatrix} \omega_{measure}\\a_{measure} \end{bmatrix}$
with $H = \begin{bmatrix} R&0\\0&1 \end{bmatrix}$

$\omega_{measure}$ is the data from gyroscope and $a_{measure}$ is the data from accelerometer (of each wheels).

Predict:

$x_k = Ax_{k-1} + Bu_{k-1}$

$P_k = AP_{k-1}A^T + Q$

I haven't really make up my mind about what Q should be.

Update:

$K_k = P_kH^t(HP_kH^T + R)^{-1}$

$x_k = x_k + K_k(z_k - Hx_k)$

$P_k = (I - K_kH)P_k$

I haven't made up my mind about choosing R.

I'm also aware that I dropped all the negative sign and the hat sign but hopefully it won't bother the readers too much.

What do you think about this module? Should I change something? Is something is not correct? This is my first time trying out Kalman filter (or any filter at all) after spending a full day learning about it.

In the future I would also like to add an encoder on each wheels to get a better reading.


texts:
I'm trying to find an object in a picture, my solution is to take the picture and a photo of the object and find the maximum of the mutual covariance, this is my project on GitHub: ImageFinder

It works but if the object is rotated or transformed it doesn't work well. Which other method should I use?

I thought that a good way can be to look for the object through the colors and then try to rotate the image or to resize it while does not match, but this idea need an improvement.


texts: Seamless looping of a signal without pops

texts:
I'm trying to seamlessly loop an audio signal. When I loop the current audio signal I have, I hear audible clicks or pops right as it loops back to the beginning. I want to be able to loop the audio signal without having any pops.
I have tried changing the signal near the end by writing an interpolation which follows this formula
$$
B_n=\frac{(A_n d_1) + (Cd_2)}{d_1+d_2}
$$
Where $B$ is the output, $A$ is the audio signal itself, and $C$ is the first signal value within the audio signal. $d_1$ equals to how many signal values until the first sample. $d_2$ equals how many signal values since the specified start. And I choose when the specified start is based on how many milliseconds I wish for the interpolation to be. In this case I gave 15ms a try. So, the specified start will be 15ms before the end of the signal.
And in simple terms; it's a crossfader that is linearly transitioning from the signal source $A$ onto the single value $C$.
The end result is this.

The transition points are highlighted via the cursor. Also, do not mind that both signals are not exactly the same. For the upper signal, I cut a portion of the start so that it would start at a non-0 value.
Since the crossfade is fairly quick for human ears (15ms), I did not manage to hear any disintegration of audio. But is this good practice? I am new to signal processing and cannot tell whether or not this will work for all types of signals I might have. Is there any other way this could be done in a more "proper" manner?


texts: How do I implement the Structure tensor (in Matlab)?

texts:
I am trying to figure out the details on how to implement the structure tensor in Matlab and need some advice!

For an image $\ I(x,y) $ the structure tensor S is given by:

$$ S=\begin{pmatrix}
W \ast I_x^2 & W \ast (I_xI_y)\\
W \ast (I_xI_y) & W \ast I_y^2
\end{pmatrix}$$
where W is a smoothing kernel, $\ \ast $ denotes convolution and subscript denotes partial derivative with respect to.

I want to use a Gaussian as smoothing kernel. This is separable:
$$ W = G_{\sigma}*(G_{\sigma})^T$$
The constant $\ \sigma $ controls the outer/integration scale.
E.g. if I truncate it at $\ +/- 2\sigma$ and use  $\ \sigma = 1$ I get:

$$ G_{\sigma} = [0.0103,0.2076,0.5642,0.2076,0.0103] $$

However the input to the Gaussian might have a mean different from 0.

How do I compensate for this?

I was thinking substracting the mean first.

I could use central differences to estimate the partial derivatives:
$$ I_x \approx [-0.5 \hspace{2 pt} 0 \hspace{2 pt} 0.5]*I $$
However since the image is noisy it is probably better to use a more robust differentian operator.
The Sobel operator is such an operator. However according to wikipedia it is not perfectly rotational symmetric. The Scharr operators are supposed to be perfectly rotational symmetric and robust approximations of the partial derivatives:

$$ I_x \approx \frac{1}{32}\begin{pmatrix}
-3 & 0 & 3\\
-10 & 0 & 10\\
-3 & 0 & 3
\end{pmatrix} * I$$

$$ I_y \approx \frac{1}{32}\begin{pmatrix}
-3 & -10 & -3\\
0 & 0 & 0\\
3 & 10 & 3
\end{pmatrix} * I$$

The Scharr operators are separable:
$$ \begin{pmatrix}
-3 & 0 & 3\\
-10 & 0 & 10\\
-3 & 0 & 3
\end{pmatrix} = [-1,0,1]*[3,10,3]^T $$

$$ \begin{pmatrix}
-3 & -10 & -3\\
0 & 0 & 0\\
3 & 10 & 3
\end{pmatrix} = [3,10,3]*[-1,0,1]^T$$

The elements of the structure tensor now becomes series of consecutive horizontal and vertical convolutions with the original image.

EDIT: Another option (thanks to nikie!) is the Gaussian derivative.
Introducing
$$ GD_{\sigma}=\frac{-x}{\sigma}G{\sigma}$$

Is this then correct?
$$ Ix = G{\sigma}^T*(GD_{\sigma}*I) $$
$$ Iy = G{\sigma}*(GD_{\sigma}^T*I) $$

How do I Implement this easily in Matlab?

How about in C++?

Any comments regarding choice of smoothing kernel and differentiation operator?

EDIT 2: My Matlab code so far:

% I: image
% si: inner (differetiation) scale,
% so: outer (integration) scale
function [Sxx, Sxy, Syy] = structureTensor(I,si,so)
[m n] = size(I);

Sxx = NaN(m,n);
Sxy = NaN(m,n);
Syy = NaN(m,n);

% Robust differentiation by convolution with derivative of Gaussian:
x  = -2*si:2*si;
g  = exp(-0.5*(x/si).^2);
g  = g/sum(g);
gd = -x.*g/si; % is this normalized?

Ix = conv2( g',conv2(gd,I) );
Iy = conv2( g,conv2(gd',I) );

Ixx = Ix.^2;
Ixy = Ix.*Iy;
Iyy = Iy.^2;

% Smoothing:
x  = -2*so:2*so;
g  = exp(-0.5*(x/so).^2);
Sxx = conv2( g',conv2(g,Ixx) );
Sxy = conv2( g',conv2(g,Ixy) );
Syy = conv2( g',conv2(g,Iyy) );

How do I normalize the Gaussian derivative?


texts: Magnification $M$ is defined as $\frac{-d_i}{d_o}$, where $d_i$ and $d_o$ denote the image distance and object distance, respectively. For a given camera, is the value $M$ typically the same for all values of $d_o$? That is, is it true that $\frac{d_{i_1}}{d_{o_1}} = \frac{d_{i_2}}{d_{o_2}}$ for all $d_{o_1}, d_{o_2}$ where $d_{i_1}$ and $d_{i_2}$ are the corresponding image distances?


texts: Origin of brightness artifacts in intensity resolution compressed images

texts:
when compressing this image

(It corresponds to Fig2.20 from the book Digital Image Processing by Gonzalez and Woods and the original .tif document can be downloaded through this
link

By reducing the intensity resolution with Matlab I obtain a similar image that exaggerates the following fluctuations brightness fluctuations in the background.

Just out of curiosity, does anybody know what the origin for those is? I am puzzled as to whether they originate

at a physical level, like small stains on the lens that took the photograph, (although they suspiciously stop near the watch)

at the transition of the signal into digital form when the image is captured by the pixels of the camera for example

at some step of compressing/ converting the image from one format to the other to obtain the .tif image. Thanks for your help!


texts: What is the maximum allowed transition band of the low-pass filter used in the reconstruction kernel for a CD with maximum speech frequency of 4 kHz?

texts: My question is if the transition band in the question below would be from 4Khz to the Nyquist frequency of 22.05Khz or from 4Khz until 44.1Khz-4Khz=40.1Khz?


texts: Why are vector quantizers not used in video codecs?

texts: Using a vector quantizer (VQ), groups of image samples (or motion compensated and/or transformed data) can be quantized as vectors. I'm curious about why VQs are not preferred in the modern, highly popular video codecs such as H.264/AVC and H.265/HEVC. A quick Wikipedia search showed that there are some codecs (such as Daala) using VQs. However, the best codecs, in terms of compression efficiency, do not use them. Is that because of the difficulty of designing a good codebook and efficient searching of the codebook?


texts:
This questions seems to be a little strange but bear with me.

I have an application that needs a digital filter to be implemented in a microcontroller.

What I want to do

I want to use the DMA to get the data from the ADC every 300us (approximately). and store it in a buffer and I'm using FreeRTOS. I want to have a task check the DMA for half transfer completion and whenever it does reach that point, then the microcontroller would execute the LP filter on half the array, while the DMA keeps on filling the rest of the buffer. Due to my microcontroller's specs I'm sure I can run this filter in those points way before it has filled the other end of the buffer.

After this it would go to a waiting mode until it finishes the later half  and when it's done it will execute the filter on the second half.

Does this make sense?

What I have done in the past

Program an interruption , when the interruption occurs put the data in a FIFO and keep filling the buffer until you reach its limit while doijg the calculations

Question

Is "what I want to do" ? Of course it can be done, but will it actually be a proper filter? It sounds strange in my head when I think about it, but it somehow makes sense (ever get that feeling?). Or do I have to do it the way I've done it in the past?

Another option would be to do it in a mixture kind of way? Fill the buffer with DMA and execute FreeRTOS every milissecond and do the filter the same way I used to do it (fill FIFO, do calculations, etc)?

EDIT:

Thank you @Ben for your response. I'm adding the same info I put as a comment here in the edit so it can be better visualized.

Desired filter - IIR filter initially

Sampling rate- 10kHz. An ADC acquisition every 100us (waaay above the Nyquist frequency, but it's necessary for what I want .. a power meter)

Tap - between 2 and 10, it doesn't need to be a high order filter

Acceptable latency- Well I'd need to use the output to do al the calculations related to power and then store it and send it to storage somewhere, so the latency can be high


texts:
Is the following system invertible or not?

$$y[n]=x[n] * x[n^2]$$ where $*$ stands for the aperiodic convolution operator. I have not been able to find a mathematically sufficient argument for it...


texts: Linear Difference Equation and Method of Least Squares

texts:
I'm reading the book "Fault-Diagnosis Systems" by Isermann in the par. 9.2.1a.
The author explains how to estimate the parameter of a linear difference equation using Least Squares.
We start with a standard ordinary linear equation with an error term:

the equation is called 9.15.
Then,

The first thing i don't understand is the syntax used by the author. What's the meaning of the | inside $\hat{y}(k|k-1)$ ?

ps: I'm sorry if what i ask is trivial, but when i have a kind of doubt like this, i don't know how to solve it if not asking it in this site...


texts:
if i have 2 sensors, sensor1 and sensor2, that sample a signal on complementary points so that sensor2 samples always between sensor1s sampling points. can I achieve the doubled sampling rate with the combination of both ?
so with merging both samples, although they are aliased individually.


texts: Pulse Doppler Radar: time delay in down conversion

texts:
My understanding of Pulse Doppler Radar is the system measures Time Of Flight for the pulse.

My confusion is that during the down conversion process prior to ADC the lowpass filter causes some time delay in addition to the time of flight. In principle the delay should be known beforehand but I'm sure it's not exact and in general the delay is not completely constant even among the band of interest.

How is this delay dealt with?


texts: Detecting objects in real world video feed?

texts:
I have an image which consists on an object in a plain black background. This is the reference image of the object. Now, My task is to take a webcam and obtain the video feed and I have to detect if the object is there in the video feed.

For instance consider I have the image of a pen with a plain black background. Now when I take move my webcam around the office space I should be able to find the object.

What is the best approach for this?


texts:
I'm trying to learn about Karhunen-Loeve transform. My source is this YouTube video : https://www.youtube.com/watch?v=Rrc2wEoCcJo&t=97s
It says the first step is to form a vector of the object in the image. They show that it is done by finding the places of 1's in the image.
But what if the image is greyscale, not binary? How are the vectors then formed?
And what if there is no clearly distinguishable object in the image?


texts: Turn acausal Butterworth output into output of a causal Butterworth without input signal

texts:
I have a dataset which contains a large number of time-series which have been filtered with acausal Butterworth filters. For a real-time application, I can only use causally filtered data.
Is it possible to turn the output of an acausal Butterworth filter into what would have been obtained if the original time-series had been run through a causal Butterworth filter?
(I do not have the original, unfiltered, time-series)

I do know the order of the filter that was applied, the corner freuquencies, and some parameter 'nroll' of which I don't know what it is ... (nroll is usually equal 2.5)


texts: finding the duration of a sequence of data

texts:
I need to find the duration of a sequence of data values that represents head yaw

movements .where the positive values are head looking to right and negative values

are head looking to left .the below signal has different duration of looking right

and left. How can I find these duration ?


texts: Finding the output for an input sequence

texts:
So I think I have the answer I just want to make sure I'm doing things correctly.

Let a be in the open unit disk, $F_s > 0$, and $0<f_0 < \frac{F_s}{2}$.  Let $h[n]=0$ for $n < 0$ and $h[n] = a^n\cos(2\pi\frac{n}{F_s}f_0)$ for $n \geq 0$ be the impulse response of an LTI system. What is the output for the input sequence $x[n]=0$ for $n < 0$ and $x[n]=\cos(\frac{2\pi f_1 n}{F_s})$ for $n \geq 0$.  Assume that $0<f_1 < \frac{F_s}{2}$.

So basically I form the convolution sum:

$y[k] = \sum_{n=0}^{\infty}a^n \cos(2\pi\frac{n}{F_s}f_0) \cos(\frac{2\pi f_1 (n-k)}{F_s})$

It would take far too long to type it all up and latex but basically I split the cosines out into their exponential forms, FOIL them, pull the exponential functions which only depend on k outside of the infinite series, and then split up the series and perform the z-transform for each one to obtain:

$$y[k] = \frac{1}{4}e^{\tfrac{2\pi ikf_1}{F_s}}\left[ H \left(\frac{1}{a}e^{\tfrac{-2\pi ik(f_1-f_0)}{F_s}}\right) + H\left(\frac{1}{a}e^{\tfrac{-2\pi ik(f_0+f_1)}{F_s}}\right)\right] + \\ \frac{1}{4}e^{\tfrac{-2\pi ikf_1}{F_s}}\left[ H\left(\frac{1}{a}e^{\tfrac{-2\pi ik(f_0+f_1)}{F_s}}\right) + H\left(\frac{1}{a}e^{\tfrac{-2\pi ik(f_0-f_1)}{F_s}}\right)\right]$$

Is this correct?  Thanks.


texts:
I'm using a cell-averaging constant false-alarm rate (CA-CFAR) detector to detect peaks like the first one in the following plot:

The magenta line is the adaptive threshold, which is the local standard deviation multiplied by a contrast factor. What other peak-detection techniques could I try in this context?


texts: Non-negative matrix factorization for audio separation - why does it work?

texts:
Non-Negative Matrix Factorization aims to factorize a matrix $\mathbf V$ into the product of two matrices, $\mathbf V = \mathbf W\mathbf H$, where $\mathbf W$ represents a set of basis vectors and $\mathbf H$ their activations or weights. Applied to audio source separation, $\mathbf V$ is the spectrogram (which is known to us), and $\mathbf W$ is the basis non-negative spectra, $\mathbf H$ is their temporal activations.

This seems to be me to be able to express a spectrogram as a sum of different frequency bins and their magnitudes - which would help to filter or separate different pitches (seemingly doing the job of a normal band-pass filter?). I do not understand how this would actually separate, say, a piano and a singer singing the same or similar pitch.


texts:
I'm working on an algorithm that counts blinks of a LED. I have series of non-neagtive integers. They represent brightness of a LED in time. I give you two examples of these integer series below. Visually it's very obvious that LED blinked 7 times in example 1 and 8 times in example 2. But how can this be calculated by a computer. Any ideas? I don't even know if I'm asking this question on the right site. I'm new to this stuff but any help would be appreaciated.

Example 1:
Data of Example 1

Example 2:
Data of Example 2


texts: Downsampling and Gaussian Filtering in the Context of Scale Space Pyramids

texts:
In the context of scale space image pyramids using Gaussian filters I noticed that it's common to downsample the image after blurring with $\sigma = 2*\sigma_{init}$ . My question is: What is the relation between variance and downsampling? What is the proof that I'm not losing information when the image is downsampled after blurring with this factor?

I think this can easily be shown using Nyquist and the cut-off frequency of the Gaussian filter, but I'm not sure how.


texts: Is the value of a frequency bin of a DFT-output the average of the 'real' frequency values within that bin's range?

texts: I am wondering whether the value of a frequency bin with a certain resolution is the average of the fourier transform values of the 'real' frequencies within that bin's range.


texts: How to remove the noise without destroying the main edge?

texts:
I am working on a project. I want to segment a logo from a car and the picture's size is 3648*2432. The logo is selected by the red rectangle.

I extracted this area and turn into gray level. The area's size is 249*173. Then, I use histogram equalize to enhance the edge. Finally, I use adaptive canny to find edge.

However, the result is not perfect.
So, I used median filter to smooth the left picture. The result is as following.

I remove the small parts.

The result has been improved, but it is also not good. Besides, I have tried  bilateral filter, mean filter and Gaussian filter, but the result is not good, too. I know that after passing histogram equalize, the noise can be enhance. But, in some case, low contrast, I need to adapt this step to enhance my picture. Is any one who can provide me any different comment? Thank you very much.


texts: plotting phase transformations of a signal

texts:
I would like to find the phase variations of the signal with respect to time axes:

S = sin(2*pi*100*t + pi/4) + cos(pi*500*t) + sin(2*pi*100*t + 5*pi/2);


texts:
In my group, we have developed an algorithm which shows abstract information from quantum mechanical systems as images. This way, given a quantum system, we obtain an associated image which has the same information and makes some features visible.

One of the important features is obtained using a "cross-correlation matrix": We divide the image into $L\times L$ sub-images and find the "overlap" between all pairs. So, the entry for sub-images $i$ and $j$, $A_{i,j}$ is a number stating how similar they are. The matrix dimension is $L^2\times L^2$.

The question is: is this matrix, or a close relative, used in image-processing? If it does, does it have a name? Does it have any interesting properties, or does it help for any useful algorithms?


texts:
As we know, if we want to get the time domain response of the channel $h[n]$, we can first get the frequency domain response of the channel $H[k]=Y[k]/X[k]$. And then we do IFFT.

But now one OFDM symbol includes some null subcarriers. So we don't need to get the $H[K]$ on these null subcarriers. But under this circumstance, how do we get the time response of the channel $h[n]$?


texts: Normalize data for transfer function determination?

texts:
I have the input and the output of a system and I want to determine the transfer function of it. I am using Ident Tool from Matlab.
How is it more ok: to normalize the output and the input between 0 and 1 or to leave it as they are?


texts: Want to do shot boundary detection via SVM, what are some good features?

texts:
I want to do shot boundary detection via SVM's.

I'm dividing the frame into nxn blocks. Per block I'm finding these features:

shannons entropy
edges (H,V,Diag)
standard deviation

for consecutive frames and corresponding blocks, the above features are compared via differences and L2 norms, I threshold over that, if above threshold then given a 1 else 0. That is, the feature value belonging to that feature for that particular block becomes 1.

This assignment happens per feature for every block. For a frame, the feature vector is the summation of all the feature vectors of its blocks.

Then I take +-N/2 frames around it and concatenate their feature vectors to create a single feature vector (for that frame).

I trained the SVM using this. There are so many normal frames and only a few shot boundary frames. I'm dividing amongst cut,dissolve, fade, other and normal frames.

I'm unable to properly detect the shot boundaries.

Can someone suggest which features I should use and how I should select them for shot boundary detection? I think I may also have not done quite well on the threshold selection, so any ideas how I could do that better?

I'm using TRECVID's 2001/2002 database.


texts: Efficient implementation of FIR comb/decimating filter

texts:
I'm working on trying to find an efficient way to implement a FIR comb/decimating filter. For a given "prototype" low-pass filter $
H(z)$, a comb that "repeats" $L$ times in the frequency domain can be constructed by upsampling (stuffing $L-1$ zeros in between each sample), yielding $H(z^L)$ (Proakis and Manolakis 342). We will denote the impulse response of the comb filter by $c[k]$.

I also wish to decimate the comb filter output by a factor of $D$. I do not wish to cascade an integrator to make a CIC filter, and I am aware that aliasing will be introduced. I am okay with this. To implement the decimating FIR filter, I have chosen to use a polyphase implementation. This means that I will have $D$ branches, each of which will contain a downsampled and delayed version of the base FIR filter (Proakis and Manolakis 767), i.e.

\begin{equation}
    p_i[n] = c[nD + i]
\end{equation}

However, since $c[k]$ is $L$-sparse, i.e. only one out of every $L$ samples is nonzero, each polyphase branch, even though they are all downsampled, will still have a lot of zeros in them. In general, $L > D$ for my application, so effectively each branch filter $p_i[k]$ is a resampled version of $h[k]$ by a factor of $L/D$ (with no interpolation/decimation filter in between the upsampling and downsampling steps).

Is there a computationally efficient way to take advantage of the fact that most of the taps are zero? Does multiplying by zero in floating point arithmetic take less computational effort or the same amount as multiplying two nonzero numbers? If so, then perhaps a naive implementation would be sufficient.

Works Cited:

Proakis, John and Dimitris Manolakis. Digital Signal Processing: Principles, Algorithms, and Applications. 4th ed. Upper Saddle River, NJ: Pearson 2006. Print.


texts:
Could anybody tell me a method to calculate the "4 Hz modulation energy" for a speech signal?

Thanks.

EDIT

I get some details : 4 Hz modulation energy: Speech has a characteristic energy modulation peak around the 4 Hz syllabic rate [3]. We use a portion of the MFCC algorithm [4] to convert the audio signal into 40 perceptual channels. We extract the energy in each band, bandpass filter each channel with a second order filter with a center frequency of 4 Hz, then calculate the short-term energy by squaring and smoothing the result. We normalize each channel's 4 Hz energy by the overall channel energy in the frame, and sum the result from all channels. Speech tends to have more modulation energy at 4Hz than music does.


texts: How to implement Minimum mean-squared error (MMSE) equalization?

texts:
I have the channel estimation and noise power estimation and I need a matlab implementation of time domain siso channel Equalization.

Theory:

Consider the FIR filter of length $N+1:$
$$
\hat{s}_{k}=\sum_{m=k-N}^{k} h_{k-m} y_{m}=\sum_{i=0}^{N} h_{i} y_{k-i}
$$
We need to find the coefficients $\left(h_{i}\right)$ that minimize the MSE:
$$
E\left(s_{k}-\hat{s}_{k}\right)^{2} \longrightarrow \min
$$
To find $\left(h_{i}\right),$ we can differentiate the error. More conveniently, start with the orthogonality principle:
$$
E\left[\left(s_{k}-\hat{s}_{k}\right) y_{k-j}\right]=0, \quad j=0,1, \ldots, N
$$
This gives
$$
\begin{array}{l}
\sum_{i=0}^{N} h_{i} E\left[y_{k-i} y_{k-j}\right]=E\left(s_{i} y_{k-j}\right) \\
\sum_{i=0}^{N} h_{i} R_{y}(i-j)=R_{x y}(j)
\end{array}
$$
In matrix form:
$$
\left[\begin{array}{cccc}
R_{y}(0) & R_{y}(1) & \ldots & R_{y}(N) \\
R_{y}(1) & \ddots & \ddots & \ddots \\
\ddots & \ddots & \ddots & R_{y}(1) \\
R_{y}(N) & \ddots & R_{y}(1) & R_{y}(0)
\end{array}\right]\left[\begin{array}{c}
h_{0} \\
\vdots \\
\vdots \\
h_{X}
\end{array}\right]=\left[\begin{array}{c}
R_{x y}(0) \\
\vdots \\
\vdots \\
R_{x y}(N)
\end{array}\right]
$$
$$
R_{v} h=r_{x y} \quad \Rightarrow h=R_{y}^{-1} r_{s y}
$$

How can I introduce the noise power estimation in my computations, and fine the delay of the wiener solution ?

Any help would be really appreciated

Here is my own matlab implementation

clear all;
clc;
depth = 50; % solution depth
ntaps = 8;
MQAM = 2;
N = 1e2;
delay=floor(ntaps/2);
s = randi([0 1],1, N);
h = [1 0.5*exp(1i*pi/6)  0.1*exp(-1i*pi/8)]';
smod = qammod(s,MQAM,'gray','InputType','bit');
r = awgn(filter(h, 1,smod),0, 'measured');
A=convmtx(r(1:depth).',ntaps);
R=A'*A;
X=[zeros(1,delay) s(1:depth) zeros(1,ceil(ntaps/2)-1)].';
ro=A'*X;
coeff=(inv(R)*ro);
decision = filter( r(1:depth),1, coeff);
dem = qamdemod(decision,MQAM,'gray','OutputType','bit');
cmp = [dem(1:10)', s(1:10)']

Here cmp is the concatenation of the input signal and dem is the output of the MMSE Equalization :we can see that because of the delay, dem and s don't match !