texts: Why does time-domain convolution correspond to frequency-domain multiplication? (visual)

texts:
I seek a visual explanation of this. I've already seen the maths, and can derive the proofs - they amount to nill for an intuitive understanding. Any amount of math is welcome, as long as serving to ultimately explain it visually.

Examples of excellent explanations:

Fourier Transform: the "winding around the circle".
Convolution: the "input-side algorithm", and how it ties to output-side.

Continuous vs Discrete: both are welcome, but ultimately the discrete case must be explainable. I've long thought of the input-side algorithm myself for continuous convolution, but never completed the picture; real analysis is tricky.
Circular vs linear: a complete explanation ought to cover both, but I'm primarily interested in linear convolution (rather, how the circular of padded signals is equivalent to linear).
Duality: ideally, should be covered (conv in freq domain <=> mult in time domain).

Ideas:

Only right-padding works, and both inputs are padded; this looks like forcing inputs to correlate with lower and fractional frequencies (and more frequencies) relative to unpadded's frame.
Something about convolving with shifted deltas and observing modulations of complex sinusoids in other domain; awaiting clarification from @AndyWalls.


texts: sketching frequency domain representation

texts: I've been given this signal and asked to sketch a frequency domain representation how do i do that? how do you know which frequency components are present in the signal?


texts: How to find the noise power in my simulation

texts:
I am currently working a tricky communication chain that I want to simulate into Matlab.

The chain is like a classic BPSK modulator/demodulator chain, but where the BPSK modulated signal has been merged with an other carrier using PM (Phase modulation) modulation method.

So, after the PM modulation step, the signal can be expressed like this ($\mathrm{h}$ is the modulation index and $d(t)$ is the NRZ data signal):

$$\cos\left( \omega_c t + \mathrm{h} d(t) \cos(\omega_{sc} t + \varphi_{sc}) + \varphi_c \right)$$

My goal is to find the curve of the $BER = f\left(\frac{E_b}{N_0}\right)$ (Which should fit the theorical BPSK BER curve).

Well, I used to work with this kind of chain without the PM step. To simulate it, I usually just use a theorical noise power $P_n = \frac{\sigma_c^2}{\frac{E_b}{N_0} * \mathrm{BPS}}$ to create a gaussian noise of standard deviation $\sqrt{\frac{Pn}{2}}$. But I have some issue in my current case: I made this chain in Matlab, but I just cannot simulate it because I don't know how to find the standard deviation I have to apply to my gaussian noise when my signal is between the two PM steps.

Actually, I have no idea about how to find it... I just have a big mental breakdown because of the non-linearity of the PM modulation/demodulation step, and I don't know what to do to find $P_n$...

I already calculated the theorical signal autocorellation/PSD just after the PM modulation step, which use first kind Bessel functions (Something similar to the expression (9) of this paper). But honestly I don't know what to do know...


texts: Speech recognition using mfcc and neural network

texts: I am doing the project on speech recognition.I have implemented both the mfcc features and neural network.I am a bit confused in joining the two blocks.For example i have three different words with their '.wav' files.I extract 13 mfcc features from these three words respectively.There will be a matrix of mfcc e.g (13x100) for each word.I convert this matrix to a vetor.Than i will have three vectors of 1300 length (assuming each word is of same length).Than i pass these vectors as input to my neural network.I am confused about my training output.Whether it will be the same wav files from which i extracted the mfcc features or will it be something else.


texts:
I was looking for a straightforward proof for Bedrosian's theorem which says the Hilbert transform of the baseband signal  times the passband signal is the original baseband signal times the Hilbert transform of the passband signal in which by assuming $x(t)$ is a passband signal and and $m(t)$ is a baseband signal, then we have :
$$c(t)=\widehat{x(t)m(t)}=\widehat {x(t)}m(t)$$
Now for the sake of proof since we know Fourier transform of Hilbert  equals to  :
$$\hat X(f)=-j\operatorname{sgn}(f)X(f)$$
Then we say :
$$\hat c(f)=-j\operatorname{sgn}(f) [x(f) \star m(f)]$$
Now can we say( because $x(t)$ is a passband signal) $X(f)$ is only in positive f part and $\operatorname{sgn}(f)$ value is constant and bring it in? and then we have:
$$[-j\operatorname{sgn}(f)x(f)] \star m(f)=\widehat {x(f)}\star m(f)$$
And so we can say
$$\hat c(f)=\widehat {x(f)}\star m(f)$$
$$\hat c(t)=\widehat{x(t)m(t)}=\widehat {x(t)}m(t)$$


texts: MIDI out from freeware sequencer to PD or Max/MSP

texts:
I'm looking for a way to send MIDI from a cheap or free MIDI sequencer to Pure Data or Max MSP so I can do realtime MIDI processing in PD/MAX. So I guess I have a couple questions:

What is a good cheap/free midi sequencer (that runs in Windows 7) that supports MIDI out to another program
How do I configure MIDI out?

Is this something normal midi sequencers can do? It seems like it should be pretty basic, but I was having trouble finding information on it.


texts:
This question is regarding applying minimum mean square operation at the receiver side of communication system.

I have read in paper that

Given an output signal from communication system $$R(k)= H(k)S(k) + N(k)$$

an MMSE filter $$G(k)= \frac{H^*(k)}{|H(k)|^2+\frac{\sigma_n^2}{\sigma_s^2}}$$

My incomplete solution

I think we need to find G(K) such that

$$\mathbb{E}[G(k)H(k)S(k)+G(k)N(k)- S(k)]^2$$

Does anyone know how should I continue, thanks alot.


texts:
I have been given the task to implement the 5/3 CDF transform for image compression.

Given that the impulse response for the low and high pass are:

$h_1 = [-0.5, 1 ,-0.5]$ (High Pass)
$h_2 = [-0.125, 0.25, 0.75, 0.25, -0.125]$ (Low Pass)

Say I want to start by computing the horizontal 1D DWT for an image. Now given that the pixels of a row obtained are $x[n]$, we can compute the convolution $y[n] = x[n]*h[n]$. At the beginning of the convolution process we have partial overlap, then full overlap and then partial overlap again. But due to the last part of the partial overlap we obtain extra values.

For the high pass filter I end up with 2 more and for the low pass 4 more. This means that the values obtained for the row pixels are now slightly larger. What am I do about this with respect to programming? Do I eliminate the last few values so $y[n]$ will match $x[n]$ in size?

Another very important question, I have implemented this regardless but I face an even greater problem. While I am fairly new to python, the algorithm I implemented is very very slow. Is there are faster way to implement this rather than using convolution? I will refrain from uploading the code here, because I am not sure if my problem is applicable to this forum or stack overflow. I shall upload it upon request.


texts: The benefits of a fractionally spaced equalizer

texts: I often read that one of the benefits of FSE vs a symbol rate equalizer, besides avoiding any aliasing in the roll off region, is that it is  less sensitive to the timing phase. What does that really mean? Does it mean that if at the output of the matched filter instead of taking one sample (assuming it is not exactly at the center of the eye diagram) , two samples are chosen, then the performance in terms of SNR will improve? How? by interpolation? I don't really understand the relation between the equalizer and what is called the timing phase ( not sure what that expression means either)


texts: How to equalize frequency contents using filter?

texts:
Suppose there is a bounded continuous signal $f(t)$ that can be fourier-transformed into $F(\omega)$ of frequency contents. $|F(\omega)|$ is non-zero for only finite number of frequencies: $f_1,f_2,..f_n$. Using a filter, we want to retain zeroes of $|F(\omega)|$ at frequencies other than $f_1,..f_n$, but we want to make $|F(\omega)| = \alpha$ for some constant $\alpha$ at all frequencies $f_1,..f_n$.

Is there such a filter? It is fine to assume that phase of $F(\omega)$ is zero or constant for all frequencies.


texts: What is appropriate windowing technique for sine burst excitation?

texts:
I had to give sine burst to get good SNR.

In above Blue color -> Actuation signal or input.
Orange -> Response or output.
I did use hanning window for both.

The frf real & imaginary plot got swapped. Which means this is not the real frequency spectrum looks like when I do sine sweep.

I used Rectangular window, I got below result

….. I tried searching, I couldn’t find any solution that solve this problem. Can you please guide me what am I missing?


texts: Compressive sensing: numerical generation of RIP matrices

texts:
The restricted isometry property (RIP) states that:
\begin{equation}
  (1-\delta_K)||x||_2^2 \le ||A x||_2^2 \le (1+\delta_K)||x||_2^2
\end{equation}
for any $K$-sparse vector $x$ of length $N$. The corresponding restricted isometry constant is $\delta_K$, $0 < \delta_K < 1$.

I would like to perform a numerical experiment based this expression to produce a RIP-compliant matrix $A$ (Matlab code is supplied below).
The matrix is constructed with Gaussian-distributed RV's (zero mean, variance $1/\sqrt{M}$). Having specified $N$ and $K$, the number of measurements $M$ is obtained with $M \geq K \log(N/K)$. This is a hit-or-miss procedure, I let the code run until the isometry constant falls between zero and one.
Is this the proper way to generate $A$? Am I missing anything? Further analysis on $A$ (eigenvalues, etc.) doesn't seem to behave properly...

function [delta] = RIP_numerical(n,s)

% delta: restricted isometry constant; n: # columns; s: sparsity

% choose m (# measurements/rows) based on measurement criterion

m = ceil(s*log(n/s))

flag = 0; incr = 0; % initialize computing parameters

while flag ~= 1  % run until condition met

    incr = incr + 1
    % generate random, m x n normalized matrix N(0,1/m)
    A = (1/sqrt(m))*randn(m,n);
    sumA = sum(A);
    for k = 1:n
        A(:,k) = A(:,k)/sum(A(:,k));
    end

    % generate vector of length n, sparsity k
    x = zeros(n,1);
    list = randperm(n);
    for l = 1:s
        x(list(l))= randn;
    end

    y = A*x;    % obtain measurements

    norm_x = norm(x,2);  % l_2 norms
    norm_y = norm(y,2);

    delta = 1-(norm_y^2/norm_x^2);  % look at lower bound

    if (delta >0 & delta < 1) % condition met; save vars and kick out of routine.
        flag = 1
        save('RIP_config.mat', 'A','x','n','m','s');
    end
end


texts: Audacity Spectogram view: how to see exact frequency and intensity under the cursor?

texts:
In the Spectogram view, simply clicking anywhere with the selection tool gives a vertical line to the top ruler, and shows at the bottom the exact sample, H:MM:SS, or other time coordinate you'd want.

I'm looking for the exact same functionality for frequency (Hz, or say musical pitch such as C7+14 cents) and intensity (-25dB).

Should Audacity not be able to do this, is there another tool that can read a WAV and do it?  I've looked at Spec and it can't but it should be easy enough to modify (I'm a software engineer).


texts: OFDM modulator: Transmission time calculation

texts:
I have a problem with a certain task.

  The frequency of the signal sampling at the OFDM modulator output is
  20 MHz. The length of the protection period is 2 microseconds and it is
  20% of the whole symbol, including the cyclic prefix. The size of the
  transmitted data is 6 MB. The 16-QAM constellation is used in the
  subchannels.

My questions are

how to calculate the data transmission time,
how many bits are transmitted in one OFDM symbol,
what is the FFT size,
what is the bandwidth of the transmitted signal and
what will be the effective transmission speed if 25% of subchannels will be
used for transmitting pilots.

Each data stream is divided into groups of "n" bits. Then for example QPSK constellation will have 2 bits at a time, 64-QAM constellation: 6 bits at a time and 16QAM: 4 bits at a time. And we are using 16-QAM constellation, so the answer to the second question will be 4 bits, right?


texts: what Graphical DSP programming tools are there for ASIO4ALL?

texts:
Many years ago I had a PCI sound card (I think it was some kind of SoundBlaster but not sure, this was maybe 15 years ago) with a 5 1/4 disk sized front panel with various ins and out.
I remember using ASIO4All but there was also a graphical component where I could program the DSP directly by creating a visual circuit and add effects and other features, and I used this to create a multitrack studio with talkback capabilities to the artists.
But I can't for the life of me remember what it was, or find anything similar on the webs.
is the 10k DSP chip still a thing? If not, what is the modern-day equivalent, and is the visual programming still available for it.
Thanks in Advance.


texts: Minimum possible sampling frequency for continuous time signal

texts:
A continuous-time signal $x(t)$ has the magnitude spectrum $\lvert X(F)\rvert$ shown in Figure 6. The signal is sampled to obtain $x(n)$ using the sampling frequency $F_s$. Determine the minimum sampling frequency, which ensures that the information in $[0; 1000 \ \rm{Hz}]$ is not distorted.

Immediately, I would think that $2000 \ \rm Hz$ is the correct answer.

Can someone explain to me why it is not and what formulas/intution is required to answer the above task?


texts: What is right and full Frobenius canonical form?

texts:
I'm having a trouble here. I'm supposed to learn Frobenius canonical representation form for finding statespace matrices, but I found many different forms.

Let's suppose we have a system with this transfer function: $$H(s) = \frac {b_ns^n + b_{n-1}s^{n-1}+ \cdots + b_1s+b_0}{a_ns^n+a_{n-1}s^{s-1}+ \dots + a_1 s +a_0}$$

For this system I'm able to find matrices as:

$$A = \left| \matrix {0 & 1 & 0& \dots & 0\\0 & 0 &1&\dots&0\\ \vdots &&&& \vdots\\0&0&0&\cdots&1\\-\frac{a_0}{a_n}&-\frac{a_1}{a_n}&-\frac{a_2}{a_n}&\cdots&-\frac{a_{n-1}}{a_n}} \right|$$
$$B = \left| \matrix{0\\0\\\vdots\\\frac{1}{a_n}} \right|$$
$$C = \left|\matrix{b_0-\frac{a_0b_n}{a_n}&b_1-\frac{a_1b_n}{a_n}&\cdots&b_{n-1}-\frac{a_{n-1}b_n}{a_n}}\right|$$
$$D = \left|\matrix{\frac{b_n}{a_n}}\right|$$
Then matrices that look the same as this but $a_n$ is probably meant to be $1$, 'cause it's missing in every matrix.
But also:
$$A = \left| \matrix {-a_{n-1} & -a_{n-2} & \cdots & -a_1 & -a_0\\1 & 0 &\dots & 0&0\\ 0&1&\cdots&0&0\\\vdots &&&& \vdots\\0&0&\cdots&1&0} \right|$$
$$B = \left| \matrix{1\\0\\\vdots\\0} \right|$$
$$C = \left|\matrix{1&0&\cdots &0}\right|$$
$$D = 0$$
or:
$$A = \left| \matrix {-a_{n-1} & -a_{n-2} & \cdots & -a_1 & -a_0\\1 & 0 &\dots & 0&0\\ 0&1&\cdots&0&0\\\vdots &&&& \vdots\\0&0&\cdots&1&0} \right|$$
$$B = \left| \matrix{1\\0\\\vdots\\0} \right|$$
$$C = \left|\matrix{(b_{n-1}-b_na_{n-1})&(b_{n-2}-b_na_{n-2})&\cdots &(b_0-b_na_0)}\right|$$
$$D = \left|\matrix{b_n}\right|$$

Somewhere it's shown as for
$$H(s) = \frac {b_ms^m + b_{m-1}s^{m-1}+ \cdots + b_1s+b_0}{a_ns^n+a_{n-1}s^{s-1}+ \dots + a_1 s +a_0}$$
and  is wrote that $$D = \left|\matrix{\frac{b_m}{a_n}}\right|$$
without denote what exactly these indexes $m,n$ mean.

Please, what is the right Frobenius canonical form?


texts: Image processing: How to find center of biggest blob in an image?

texts: Let's say I have an image like this one, and have masked it so all yellow pixels are 1 and all black pixels are 0. I would like to (1) fit the largest possible circles containing solely yellow pixels inside the yellow blobs; multiple circles may fit in irregular blobs, (2) from all such circles, identify which has the biggest radius, and (3) find the location of its center. What image processing techniques can I use for this? I am familiar with some morphological operations (like these), but can't get to where I want to go yet.


texts: Analyze the transient response of a system

texts:
Suppose I have a sytem with the following transfer function :
$$H(s) = \frac{N_H(s)}{D_H(s)}$$
I would like a general method which is not dependant on the order of the system to analyze what would be the step response.
I know what is my transfer function and the order of the system is high (superior to 5). In this transfer function I would like to modify certain parameters to see how evolve the response time. Nevertheless I do not want to make my eyes work and say by looking the step response "Well that is what I wanted to have " . I want something that could optimize the different parameters via a programm. I was thinking of using the quality factor but I do not how to evaluate it on a hight order transfer function. Do you have an ider or other idea that could help me ? (I would prefer to use python rather than Matlab)
Thank you very much and have a nice day :)


texts: Comfort Noise Insertion in VAD algorithms

texts:
I am in a intention of Porting Voice Activity Detection algorithm in real time hardware,so i was successful in classifying speech and silence/noise but now the problem is replacement of silence frames with comfort noise(preferebly any stationary noise)

so can any one suggest any opensource CNG algorithms or how a comfort noise can be of any type(pink,white etc)

to do this usually we use a rand() functions but where as these type of functions are not available with real time DSP so how to port such type of Algorithms?


texts: How to compute dBFS from rms-values if the samples exceed 32bits

texts:
For a first approach on audio-sampling I'm currently working with the raspberry pi's pcm pins and MEMS-microphones, using the ALSA library.
A problem arises when I want to calculate dBFS from a buffer's rms-value because i have to assume that a sample can be as big as (2^32)/2. In order to obtain rms, I would have to square every single sample, some of which might not be computable due to their size.
What would be a good way of obataining the dBFS of a finite array of discrete values without causing overflows?


texts: How can $a_n$ for the trigonometric Fourier series be $0$ for a non-odd square signal?

texts:
Can $a_n$ for the trigonometric Fourier series be $0$ for a non-odd wave? Since at the moment the square wave is not even nor odd. Square wave is given by:
$$
s_2(t)=\begin{cases}
 -4, & 0\leq t\leq 2.5\\
6,& 2.5\leq t<5
\end{cases}
$$

The signal looks like this:

When I find $a_n$ it gives me $0$ but when you get a value of zero, doesn't this mean that the function is odd?
OR DOES THE ZERO ALSO IMPLY THAT IT IS NOT AN EVEN and clearly the wave isnt even?

EDIT:

If this wave is odd, how is it odd because I thought it had to be $s(t)=-s(-t)$ so for e.g. $6$ doesn't = $-4$?


texts: How do I apply a function window to a signal?

texts:
How can I apply a window function like Hamming or Lanczos to a signal, using its coefficients?

I mean, which method can I use to do this? FFT? Convolution? Which method has the better performance?


texts:
We used to classify signals as 1D and 2D etc ie one dimensional and two dimensional. For example a periodic square wave signal is 1D and an image is a 2D signal etc (reference - Signals and systems by Simon Haykin and Barry Van Veen, 2nd edition , page 2).

But the same periodic square wave signal can be decomposed using fourier series to infinite sinusoids with different frequencies. In the linear algebra terms these infinite orthogonal sinusoids forms the basis and the the dimension of a periodic square wave is infinite.

So actually which is the actual dimension or what is dimension?


texts: How to detect missing bit position in a bit-stream?

texts:
I'm trying to establish a one-way communication and facing problem in the receiver end because of some missing bits. For example, if the sender sends 100 bits the receiver received 95/96 bit. Now the position of these missing bits in the received bitstream is non-deterministic. Sometimes bit get missing at the beginning (for example after 50 bits), sometimes at the middle or sometimes at the end. Basically,

My received bitstream length is less than the sending bitstream length.
Is there any way I can detect, at which position these bits get missing?


texts: I know very basic image processing. Please suggest me how do I and where do I start learning Gabor Filters ?


texts: Determing inverse Z-transform using impulse response?

texts: In Matlab there is a command iztrans for finding inverse Z-transform. But how can we find inverse Z-transform using impulse response? The Matlab command impz gives Impulse response of digital filter but how can I then proceed towards finding inverse Z-transform?  Actually i want to determine inverse ztransform of H(z) shown in attached photo. But i want to determine inverse z transform using impulse response since i read somewhere that inverse ztransform can be determined by impulse response using matlab but i am unable to find any google resource that guides much in this regard. I have also attached a table of common z tranforms. I have also highlighted the two most relevant cases as far answer of GKH is concerned


texts:
Because of interesting spectral phenomena, I would like to shift the >20kHz content of my audio recordings to audible frequencies, so that I can listen to it.
I tried the Audacity pitch correction feature, but didn't succeed.
Any help much appreciated.


texts: Uncertainty Relation between time and frequency

texts: Typically, humans can hear sound waves in the frequency range 20 Hz to 20 kHz. If one wants to make digital record of sound such that no audible information is lost, what is the longest time interval between samples that could be used?


texts: books on audio processing open source projects and software

texts:
Because I would like to

understand how open source projects work and give my contribution
improve my programming skills and
understand better and implement algorithms in a real project

I'm looking for books.

The programming language is not the most important aspect: C/C++ would be perfect but any of the most common language is good.

A good example of open source project, and a very famous one in audio, is Audacity: the source code is easily available. However I could not find any book guiding or with implementation examples on existing source code for Audacity. Something on Audacity would be perfect.

Does anybody have a suggestion?


texts: Question about poles and zeros in AR filter

texts:
For AR HP filter

  zeros, to the right of the imaginary axis
  poles outside the unit-circle
  zeros on the real axis
  poles, to the left of the imaginary axis

Apparently, the right answer is (4). Why?

We did learn in class some properties of the zeros and poles, and I probably need to relate it to the fact that it's a AR HP filter.

I'd be glad for your help.

Thanks!

EDIT:
I know that AR is an all-poles filter so it must be (1) or (4). What is the meaning of the position of the pole with relation to the imaginary axis?

EDIT 2:
Is it because the large $\omega$'s (frequencies as angles of the unit circle) are at the left side of the imaginary axis and so we get an HP filter?


texts: Computing shifted signal without first reconstructing

texts:
Looking for a solution to the following problem:

  A signal $x(t)$ is band limited to $B$ Hz, and sampled above the
  Nyquist rate, with corresponding $f_s = 1/T$. If the sampled signal is
  given by

  $$ y(n) = x(nT), $$

  describe how to compute the phase-shifted sampled signal

  $$ z(n) = x(nT + T/2), $$

  from the samples $y(n)$ without first reconstructing $x(t)$ from
  $y(n)$ and resampling.

Have not been able to find one on the web.


texts: Determining the 'amount' of speaking in a video

texts:
I'm working on a project to transcribe lecture videos. We are currently just using humans to do the transcriptions as we believe it is easier to transcribe than editing ASR, especially for technical subjects (not the point of my question, though I'd love any input on this). From our experiences we've found that after about 10 minutes of transcribing we get anxious or lose focus. Thus we have been splitting videos into ~5-7 minute chunks based on logical breaks in the lecture content. However, we've found that the start of a lecture (at least for the class we are piloting) often has more talking than later on, which often has time where the students are talking among themselves about a question. I was thinking that we could do signal processing to determine the rough amount of speaking throughout the video. The idea is to break the video into segments based on the amount of typing they require, as opposed to splitting based on time.

I've done a little research into this, but everything seems to be a bit overkill for what I'm trying to do. The videos for this course, though we'd like to generalize, contain basically just the lecturer with some occasional feedback and distant student voices. So can I just simply look at the waveform and roughly use the spots containing audio over some threshold to determine when the lecturer is speaking? Or is an ML approach really necessary to quantify the lecturer's speaking?

Hope that made sense, I can clarify anything if necessary.

Appreciate the help as I have no experience with signal processing.


texts: How do you change pitch with wavelet transforms in MATLAB?

texts:
I'm trying to increase the pitch of a wav file using MATLAB. (I don't expect it to sound very good, but I'm trying to figure out how wavelets in MATLAB work.) Here's what I'm doing:

scales = {2.34375,0.5,21};
result = cwtft(b5,'scales',scales);
result.scales = 2*scales;
b5high = icwtft(result);

(b5 is just the opening of Beethoven's Fifth, but it could be any sound file. It happens to have a sample rate of 48 kHz.)

To my thinking, I'm doing a wavelet transform of b5 using a Morlet wavelet. The scales that I'm using capture the human-audible frequencies for a 48-kHz signal. I then double the scales vector--resulting in a shift upward of 1 octave, and then convert back into a signal. However, when I play b5high, it sounds exactly like b5 - no pitch shifting at all.

What's wrong with my logic?

Thanks.


texts: How does one know how many times to downsample in Harmonic Product Spectrum?

texts:
I have been working on a Harmonic Product Spectrum Algorithm. Now, all the literature I've read about the subject tells me to downsample an N amount of times. How does one determine what this N value should be? Here is my implementation of Harmonic Product Spectrum so far. Feel free to tell me if I have gone wrong somewhere.

private int HarmonicProductSpectrum(Complex[] fftData, int n){
    Complex[][] data = new Complex[n][fftData.length/n];
    for(int i = 0; i<n; i++){
        for(int j = 0; j<data[0].length; j++){
            data[i][j] = fftData[j*(i+1)];
        }
    }
    Complex[] result = new Complex[fftData.length/n];//Combines the arrays
    for(int i = 0; i<result.length; i++){
        Complex tmp = new Complex(1,0);
        for(int j = 0; j<n; j++){//multiplies arrays together
            tmp = tmp.times(data[j][i]);
        }
        result[i] = tmp;
    }
    //Calculates Maximum Magnitude of the array
    double max = Double.MIN_VALUE;
    int index = -1;
    for(int i = 0; i<result.length; i++){
        Complex c = result[i];
        double tmp = c.getMagnitude();
        if(tmp>max){
            max = tmp;;
            index = i;
        }
    }
    return index*getFFTBinSize(fftData.length);
}


texts:
I'm attending a very practical course on signals and i have some doubts, i hope to receive answers in layman terms.

1) My prof said i can use the autocorrelation of the output of a process to diagnose faults. I'm trying to understand how. Maybe by calculating the autocorrelation in different $\Delta t$ and notice changes ?

2) Is it possible to calculate the autocorrelation in subsequent $\Delta t$ and then compare the Power Spectral Densities of each  $\Delta t$ to notice faults ?

3) Why autocorrelation is considered a statistical approach? To me, it's simply multiplications (and a division) operations applied to signal samples...

4) In which cases i'm forced to use autocorrelation instead of FFT to diagnose faults ?


texts:
Is there a rule of thumb for calculating the approximate filter order of a FIR/IIR based on passband vs sample rate?  I seem to remember reading about it once in Undergrad. I'd think you could make a rough estimate like scaler * passband/stopband.

I have been working on implementing different FIR/IIR filters and would like some rough formula that says if one variable is changed, this is approximately how many more taps you would need to get similar performance.


texts:
When capturing microphone audio, I noticed that my static buffer eventually overflows.
So I timed the samples, and low and behold I'm getting an actual sample rate of 44102-ish Hz (my timer is very accurate).
Is this the actual sample rate recorded by my microphone, or is windows sending me samples twice?
Doesn't this cause all sorts of synchronization problems. Say if you're listening to a live recording, wont a noticeable lag start to accumulate if sample rates aren't exact? How is this problem solved traditionally? Is dropping a sequence of samples, or padding zeros normal?


texts: Why are linear phase filters called so, if they provide Constant delay instead of linear

texts:
Linear phase filters delay all frequencies by the same amount. Why aren't they called Constant phase filters instead of Linear phase?

As I understand, if there is an input signal with two components $f_1$ and $f_2$, and it is passed through a linear phase system which introduces a delay $t$, $f_1$ is delayed by $t$, and $f_2$ is also delayed by $t$.

To what property does the phase have a linear relationship with? Could someone provide some insight into this? Some Math would be helpful.


texts: Get position of different frequencies in sound signal

texts:
Assumed I have a signal where different frequencies appear depending on the time position, i.e. a signal of $440\,\textrm{Hz}$ appears at $t=0$, a signal of $441\,\textrm{Hz}$ appears at $t=1$, a signal of $442\,\textrm{Hz}$ appears at $t=3$. Is it now possible to find out when the different signals appeared, and when they leave again, and relate that to the time position? The normal fft in one dimension just tells me that there are different frequencies, but not which.
My application is, that I have a sound file taken from an interferometer. The sound is generated by the interfering light beams. Different wavelengths generate different interference frequencies, but at different times. Thus I want to get the different frequencies in the resulting signal, and the time when they begin/end. The interference pattern generates a sine wave while interfering.


texts: What could be the optimal filter given this situation?

texts:
I am using a camera and a Deep Neural Network to predict one angle. This network received as input the frame, and calculates the mean and the variance associated to the prediction (which is basically its uncertainty), so that the output of the systen consists of a normal distribution.
I want to use this extra information (the uncertainty of the prediction) to better aproxmate the angle I want to predict. Since I have no experience with this, I don't know which filters could be better:
I have implemented a Kalman Filter that works really well with linear environments, as expected. However, if the angle starts to change non-linearly, which filter could use the information I have obtained (mean and variance for each sample) to improve the accuracy of the system?


texts: Video: Who crossed the finishline first?

texts:
Okey, so there was a boat race during the weekend. It was a very tight race, but a few minutes after the race, the green boat was declared winner.

Given only this video images, how can they determine this?
http://kvf.fo/kapprodur?sid=70276
Crossing finish line around 5:17
(Download best quality press "Tak niður" -> "Hágóðska")

Iv tried with ffmpeg to fiddle around with colors but i have no luck to determine if in fact green or yellow won. In fact the colors seem to very distracting.

Any suggestion to how I can render a more determinable result?


texts: How do I check for linearity for the following piecewise-defined system?

texts:
The problem at hand:

Where I'm currently stuck:

I'm not entirely sure about how to move on from this point, I'm trying to find the superposition of the responses of the two individual signals so I can later check if that's equal to the response of their (the two signals') superpositions. The intervals in each response involve different signals and thus I'm unable to find the superposition of responses.


texts: Generate impulse from random frequency response

texts:
This code supposed to generate random frequency response which could be used as impulse...But it fails to pack correctly when going from fft domain to real signal time domain, the other half is left inside the real signal...

clear all,
close all,

mn = 37.2700     %minimum amplitude value in dB
mx = 140.3640    %maximum amplitude value in dB
N = 44100        %length of the impulse

ampr = mn + (mx+mx)*rand(N,1) %random values for amplitude

a=db2mag(ampr)        %convert dB to magnitude
a1=a(2:end)
a2 = flipud(a(2:end))          %create a symetric mirror for negative frequencies
a3=[a1;  a2]            %combine positive and negative frequency amplitude arrays

pmn =-180        %minimum phase
pmx = 180        %maximum phase

phaser = pmn + (pmx+pmx)*rand(N,1)   %random phase values

ph=wrapToPi(deg2rad(phaser)) %convert phase degrees to radians values between pi and -pi
ph1=ph(2:end)
ph2 = flipud(-ph(2:end))    %flip and create mirror phase values for negative frequencies
ph3=[ph1;  ph2]                %combine two mirror symetric halves of positive frequency and negative frequency phase arrays

mk2=real(ifft(pol2cart(ph3,a3)))  %convert polar values to cartesian and finally to real time domain signal
fs=44100                          %assign sampling rate
% mk2 = mk2(1:floor(length(mk2)/2)) %remove the other half caused by negative frequencies
mk2=mk2/max(max(abs(mk2)))       %normalize the signal

audiowrite('RandomGeneratedImpulse.wav',mk2,fs,'BitsPerSample',24)


texts: Which transformation in frequency domain equals a x-axis shift of a signal in time domain?

texts:
I have discrete Fourier transformation results from measurements. Looking at the signal from the time domain perspective, I want to shift the signal on the $x$-axis to the left or right.

Which transformation in the frequency-domain results in a left/right shift in the time-domain?


texts: Question about the derive of filter roll-off

texts:
I don't understand the follow step in my notes:

It is about the filter roll-off.

$H(w)={1\over 1+j(w/w_c)}$

$\Rightarrow$
$|H(w)|= {1\over \sqrt(1+(w/w_c)^2)}$

where $H(w) = {v_o \over v_i}$

$w_c$ = cut-off frequency


texts:
Say, I create a Hermitian complex signal using,

import numpy as np
t = np.arange(-4, 4)
z = np.exp(1j * t)

Here z should be a complex signal with Hermitian Symmetry, as you can see below.

In [2]: t
Out[2]: array([-4, -3, -2, -1,  0,  1,  2,  3])

In [3]: z
Out[3]:
array([-0.65364362+0.7568025j , -0.98999250-0.14112001j,
       -0.41614684-0.90929743j,  0.54030231-0.84147098j,
        1.00000000+0.j        ,  0.54030231+0.84147098j,
       -0.41614684+0.90929743j, -0.98999250+0.14112001j])

But when I take the Fourier Transform of this signal, I don't get a real spectrum,

In [8]: np.fft.fft(z)
Out[8]:
array([-1.38531768+0.7568025j, -7.02599565+0.7568025j,
        2.57935201+0.7568025j,  0.93952139+0.7568025j,
        0.41344309+0.7568025j,  0.08151870+0.7568025j,
       -0.22205190+0.7568025j, -0.60961892+0.7568025j])

However, when I take the Hermitian FFT using the hfft function, I get

In [10]: np.fft.hfft(z[:5])
Out[10]:
array([-1.38531768, -7.02599565,  2.57935201,  0.93952139,  0.41344309,
        0.0815187 , -0.2220519 , -0.60961892])

This is the real component of the result I got using the regular fft function.
I don't understand what I am doing wrong here. Shouldn't fft be giving a result with a zero or almost-zero imaginary component when the input is Hermitian?

I feel like I am doing something wrong here regarding how the signal is presented to fft.


texts: do i have an image processor to compare an image with a number of images and choose the identical image??

texts: I require an image processor which compares an image with a set of images and finds the equivalent image out of them.


texts: High frequencies disappear when applying discrete wavelet transform

texts:
Trying to decompose and reconstruct a signal using a to some extent self-made implementation of DWT for some reason fails. The result looks highpass filtered and/or shifted. I wanted to write the code in Python not using the PyWavelets package (except for getting the wavelet object).

import numpy as np
from scipy import interpolate, signal
import pywt
import matplotlib.pyplot as plt

# Test signal and wavelet
sig = np.zeros((10*512))
for i in range(256):
    sig += i * np.cos(2*np.pi*i*np.linspace(0, 10, 10*512)+i**2*2*np.pi/((i+1)/1024))
w = pywt.DiscreteContinuousWavelet('db6')

I convolved the test signal with the decomposition filters from the wavelet object and used scipy's signal.decimate(...). As far as I understood, the decomposition filters do not lowpass filter the signal (or do they and is this where I went wrong?). However, a simple downsampling instead of decimation did not yield better results.

# Downsampling by decimation, decomposition
coefs_own = []
cA = np.copy(sig)
for i in range(1):
    cA_last = cA
    cA = np.convolve(cA_last, w.dec_lo, 'same')
    cD = np.convolve(cA_last, w.dec_hi, 'same')
    cA = signal.decimate(cA, 2)
    cD = signal.decimate(cD, 2)
    coefs_own.insert(0, cD)
coefs_own.insert(0, cA)

And this is my reconstruction using linear interpolation and a convolution with the recomposition filters provided by the wavelet object.

# Helper function
def upsample(data):
    x = np.array(range(data.shape[0]))
    x_double = np.linspace(x[0], x[-1], x.shape[0]*2)
    f = interpolate.interp1d(x, data, kind='linear')
    data_upsampled = f(x_double)
    return data_upsampled

# Upsampling by interpolation, reconstruction
cA = coefs_own[0]
for i in range(1):
    cA_up = upsample(cA)
    cD_up = upsample(coefs_own[i+1])
    cA = np.convolve(cA_up, w.rec_lo, 'same') + np.convolve(cD_up, w.rec_hi, 'same')
sig_rec_own = cA

For comparison, see the multilevel decomposition by PyWavelets, which does a perfect reconstruction.

# PyWavelets
coefs = pywt.wavedec(sig, w, mode='zero', level=1)
sig_rec_pywt = pywt.waverec(coefs, w, mode='zero')

# Plot both decomposition versions
fig, axes = plt.subplots(2, 1)
axes[0].plot(coefs_own[0])
axes[0].plot(coefs_own[1])
axes[0].set_title('Own')
axes[0].legend(['Approx.', 'Detail'])
axes[1].plot(coefs[0])
axes[1].plot(coefs[1])
axes[1].set_title('PyWavelets')
axes[1].legend(['Approx.', 'Detail'])

# Plot both reconstructions
fig, axes = plt.subplots(3, 1)
axes[0].plot(sig)
axes[0].set_title('Original')
axes[1].plot(sig_rec_own)
axes[1].set_title('Self-made reconstruction')
axes[2].plot(sig_rec_pywt)
axes[2].set_title('PyWavelets reconstruction')

plt.show()

Does anybody have an idea, why the decomposition and the reconstruction both obviously are erroneous? What did I miss in my implementation?


texts: Decomposition of $H(z)$ as maximum-phase, minimum-phase

texts:
The frequency response is:

$$H(z) = 2-7z^{-1}+7z^{-2}-2z^{-3}$$

I see that it has $3$ zeros: $z_{01} = \frac 12$, $z_{02} = 2$, and $z_{03} = 1$; and $3$ poles in: $$z_x = 0$$

Now, I have to write it like:

$$H(z) = H_{\rm min}(z) H_{\rm max}(z) H_{\rm uc}(z)$$

where $H_{\rm min}(z)$ is the minimum-phase frequency response, $H_{\rm max}(z)$ is the maximum-phase frequency response and $H_{\rm uc}(z)$ only has zeros on $\lvert z\rvert=1$.

For $H_{\rm uc}(z)$, I have one zero on $\lvert z\rvert=1$, so:
$$H_{\rm uc}(z) = z-1$$
Is that OK?
For the minimum-phase frequency, I have all the poles/zeros that are inside the unit circle, so:
$$H_{\rm min}(z) = \frac{z-1/2}{z^3}$$
For the maximum-phase frequency, I have all the poles/zeros that are outside the unit circle, so:
$$H_{\rm max}(z) = z-2$$

But I also know that I can find $H_{\rm max}(z)$ as:
$$H_{\rm max}(z) = H_{\rm min}\left(z^{-1}\right) z^{-M_i}$$
where $M_i$ is the quantity of zeros of $H_{\rm min}(z)$. So:
$$H_{\rm max}(z) = \frac{z^{-1} - 1/2}{z^{-3}} = z^3 \left(z^{-1} - 1/2\right) = z^2 \left(z-1/2\right)$$

So I have $2$ different expressiones for the same $H_{\rm max}(z)$. What am I doing wrong?


texts:
Is it accurate to say in digital audio that, when a fader is down, then its value is "-$\infty$"?

What has confused me is, whether it's really $-\infty$ or some lower bound of the amplitude. it could be that the $-\infty$ is just a "joke" by the developers and it doesn't have mathematical validity.


texts:
Given a length-n vector (the elements of vector describe samples of any time changes).

There are $m$ elements, and other values $n-m$ (coming right after the first n elements)

I need to predict what their value of other $n-m$ elements .

I would like to know which prediction method is recommended to use for solving the above problem?

for example:

%create vector of sampling
Fs= 1/200;
n=200;
m=50;
t=0:Fs:2*pi;
x=sin(t);

%call prediction function
fprediction(x, Fs ,m); % x:vector, Fs: frequency sampling, m:number of element to predict
...

I need to predict what their value of other n−m (50) elements, relying only on the first 150 elements .

Comment:

The use in sampling rate to resolve the problem  is optional

Thanks


texts: What is proof of Basic Hilbert Transform?

texts: Attached is the question. Need to find the output for signal x(t).


texts:
I'm very new to digital signal processing, I've been reading a text book to try to figure out this assignment at work (avionics related).  I'm supposed to implement a two-pole low-pass filter with a fixed cutoff frequency (of 0.02Hz)

I'm not sure if my current output is correct, given my understanding of how low pass filters are supposed to work.

The top is my input signal. 3Hz to the left, 1Hz to the right.

Given my current level of understanding of digital signal processing, my expected output of a low-pass filter with a frequency cutoff of 2Hz is in the middle. 0 to the left, and the original 1Hz signal on the right.

The bottom is the output of my current implementation. I don't know if my current implementation is wrong, my understanding, or my expectations of how close to ideal the filters actually are.

This actual implementation is a single pole recursive filter, just trying to get a working example (and understanding) I can test against, before I design everything the way they want it.

The sample rate is 100Hz, so for f = 2Hz, I tried poles of both x = e^(-2*pi*2/100) and x = e^(-2*pi*2/200) (the formula I saw wanted a value from 0 to .5, so assumed it was some fraction of sample rate?)

Any ideas?


texts: How to evaluate a close-loop control during transient?

texts:
I have implement a close-loop control in order to perform the tracking of the desired signal R.
The basic scenario shall be summarized as follow:

I have to evaluate the performance of my control G with a moving average value of the Error E.
When there is no transient the error E shall be less then a certain threshold.

But, how to relax the threshold during transient ?
Addition info:

Frequency of R is limit to 2 HZ
Rising time is 200ms
Control G has a 1 Hz band (but this is not a problem)
I can acquire R E U and Y
Real-time FFT of R is not possible during execution.

My idea is to find a threshold as function of r signal which is small when reference is stationary (constants from at least 200 ms) and high when the reference is fast changing.

Thank you all
Regards


texts:
I want to generate and play sound for tones with frequency 100 Hz and 200 Hz for 10 sec and 5 sec respectively using MATLAB with sampling frequency equal to 20500 Hz.

Also I want to plot their frequency spectrum. Is it possible to do this using MATLAB?


texts: Why are numerator coefficients named $b$ and denominator coefficients named $a$?

texts:
This is not a technical question, it is rather a question about the implicit notation used in MATLAB and some digital signal processing books to refer to the numerator and denominator coefficients of impulse responses or system functions.

When I think of a fraction the first thing that comes into my mind and the first thing I write is the numerator and then the denominator. When I spell the alphabet the first letter is a and the second is b.

When I started reading DSP books and working with Matlab I found:

$$ Y(z)=\frac{b(1)+b(2)z^{-1}+\cdots+b(n_b+1)z^{-n_b}}{1+a(2)z^{-1}+\cdots+a(n_a+1)z^{-n_a}}X(z)$$

and y = filter(b,a,x) respectively

In these cases the first things to write were labeled b (numerator coefficients) and the second things to write were label a (denominator coefficients). This caused me a little bit of confusion at first as it seemed a little bit counter-intuitive. Is there any reason why this was selected in this way?

I'm thinking that maybe poles are considered more important than zeros (due to stability, etc.) so they most be though of first and a comes first in the alphabet... but this is just an educated guess.


texts: Removing outlier data points from frequency-domain signal

texts:
Hello signal processing stack exchange,

I have a complex frequency-domain signal representing the impulse response of a physical system over some range of frequency space. As an artifact of the simulation, some of the points are clearly off from what would otherwise be a smooth and Kramers-Kronig-obeying curve. Any recommendations on how to interpolate this curve while entirely ignoring these outliers?


texts: How to model/estimate channel in discrete time with different sampling time?

texts:
I want to estimate a channel based on LTE 3GPP EVA with given power delay profile (set of average power and delay of channel taps).

tau = [0, 30e-9, 150e-9, 310e-9,  370e-9,  710e-9, 1090e-9,  1730e-9, 2510e-9]; % relative delay (s)
pdb = [0,  -1.5,   -1.4,   -3.6,    -0.6,    -9.1,    -7.0,    -12.0,   -16.9]; % avg. power (dB)

Unfortunately, the relative delays of channel taps are not multiples of my sampling time Ts=3.3333e-07 s.
Here, I assumed the number of taps to be estimated is Ds/Ts, where Ds is the delay spread.

First of all, please let me know if there is anything wrong in my approach.
Also, How should I measure my estimation error in this case? Should I compare my estimates with channel sampled at my sampling rate, or should I subtract the reconstructed frequency response of estimate and the original one?
Anyhow, due to the difference in the sampling time and delays, there is an intrinsic error in my estimation no matter how I measure it.

Thanks in advance :)


texts: Why do we follow Shannon's Theorem in audio sampling?

texts:
My professor told me that in sampling a signal, is not usual to respect the Shannon theorem because the interpolation formula contains sines and fractions which are difficult to calculate with a microprocessor.

Why is audio sampled at 44.1-48 KHz (Nyquist frequency of audio signals)? Is there a dedicated  hardware for accomplishing this? I can sample sound with my PC, so does this mean that I have this hardware in my PC?


texts:
I am trying to understand the LSM algorithm applied to grayscale image segmentation.

There are essentially 2 things that are blocking me:

1) From my point of view - the level sets - i.e. ''moving'' the 3D surface (that is, conceptually visualizing the grayscale image (say 8 bits/pixel) as a surface whose height is given at each point (pixel $x,y$) by the intensity value of that pixel, through a fixed plane, or the other way around, is exactly the same as performing a thresholding from 0 to 255 (i.e. each "level set" would be all the pixels at threshold value 0,1,2 ...255)?

2) Let $f(x,y)$ be the image function (associating an intensity between 0 and 255 to each point $(x,y)$. The LSM
is defined as $$\Gamma = \{(x,y) | \phi(x,y)=0\}$$

However, it does not make sense to me, (from my point of view which is probably wrong, hence my question), that $$\phi(x,y)=0$$ as as the surface has its own height at each point and it would make more sense to me to consider that the plane intersecting the surface is moving from bottom to up, i.e. I can ''accept'' better the formulation given for the level sets:

$$\Gamma = \{(x,y) | \phi(x,y)=c\}$$

and for the example I use, $$c \in [0,255] $$

Question 2 is more like a ''conceptual'' problem, but question 1 is more fundamental, since I cannot seem to see any difference with a simple thresholding.

Rephrased:

Question 1) Can we consider (loosely speaking, as indeed as the li is a set of points not a matrix) that each level curve $$l_i$$ of a level set of an image I is (in programmatic terms, Matlab syntax): li= (I==i) ? (where i is the height, in grayscale 8 bits, i is a number between 1 and 255, and as Matlab uses matrices everywhere, I is a $n \times m$ matrix and the operation == checks for every pixel in I if it is equal to i and the result is a binary image with 0s or 1s, i.e. this is called thresholding, hence my question is LSM the same as thresholding with operator == (and the single '=' is the assignment operation of course))

Question 2) Is it essentially (conceptually) the same to imagine moving the shape through a fixed plane, or moving a plane through a fixed 3D shape ?

What I was (in an incorrect way) trying to say is that we can visualize the level curves $$l_i$$ as the white circles on the thresholded image (although not perfect since it is discretized an maybe rounding errors/noise etc but suppose they are continuous) and ideed i do understand that the result of thresholding is a matrix so strictly speaking not equal to a lvl curve. But "loosely" speaking, visualizing the lvl curves do you agree?


texts: Why does my sinusoid look "AM" in shape?

texts:
My code is :

Fs=200e6;
Ts=1/Fs;
NFFT=2^14;
Runtime=(NFFT-1)*Ts;
t=0:Ts:Runtime;
f_in=90*1e6;
y_in=sin(2*pi *f_in *t);
plot(t,y_in)
ylim([-1.5 1.5])

Then why does my plot look like amplitude modulated when you zoom into it?


texts: Is it common to find systems having signal to noise ratio (SNR) less than 1 i.e. noise is more than signal? I am currently working with SNR of -5db to -20db. And how useful is it to work with sensor readings having negative SNR? I am asking it because the algorithm I am using to solve a problem works well with negative SNR but at positive and higher SNR, it is numerically unstable. This instability arises because of computation of hypergeometric function $_pF_q\left(;;\frac{z}{σ^2}\right)$. At high SNR, i.e. low $σ^2$, the argument $\frac{z}{σ^2}$ becomes very high, making the convergence of hypergeometric function very slow and the algorithm becomes numerically unstable.


texts: Why are Gaussian filters used as low pass filters in image processing?

texts:
In 1d signal processing, many types of low pass filters are used. Gaussian filters are almost never used, though.

Why are they so popular in image processing applications? Are these filters a result of optimizing any criterion or are just ad hoc solution since image 'bandwidth' is usually not well defined.


texts:
I am using complex morlet wavelet to estimate frequencies and damping ratios in a 1-D signal. I have done it so far, and the results are accurate enough. Now my question is how can I modify or optimize my method to get better result. like mixing it with other method such as Neural networks or anything else. I know this question is a bit weird But I need to modify my method somehow.
Kind regards.


texts: Constructing a Gaussian kernel in the frequency domain

texts:
I'm currently learning about Fourier transform, but find the differences between spatial domain and frequency domain a bit confusing at times.
Let's say I would like to perform convolution of an image with a Gaussian kernel. As far as i understand the Fourier transform of a Gaussian is also a Gaussian i.e the Fourier transform of
$$
g(x,y;\sigma) = \frac{1}{2\pi\sigma^2} \cdot e^{- \frac{x^2+y^2}{2\sigma^2}}
$$
is
$$
G(x,y;\sigma) = e^{- 2\pi\sigma^2(x^2+y^2)}
$$
However if I directly construct a kernel using $G(x,y;\sigma)$ the output kernel looks completely different from first computing $g(x,y;\sigma)$ and then doing the FFT. What is the correct approach here? I hope the question is clear, but it is essentially how I can construct a gaussian kernel directly in frequency space so i do not have to FFT the kernel.


texts: Multiresolution discrete wavelet transforms: First steps

texts:
I have been looking at trying to use wavelets for multichannel audio dynamic range processor. I am in the very early stages of trying to get decomposition and reconstruction working. In this thesis, it is suggested that a multiresolution dwt can be calculated using a tree of complementary filters. Being a hands-on person, I decided to jump straight in with the following Matlab code . . .

n = 31;
flipVec = [repmat([1 -1],1,floor(n/2)) 1];

% Get analysis and reconstruction filter weights using order and sign flips
bH0 = fir1(n-1,0.5);
bH1 = fliplr(bH0).*flipVec;
bF0 = fliplr(bH0);
bF1 = fliplr(bH1);

x = [1; zeros(1023,1)];

y0 = fftfilt(bH0, x);
y1 = fftfilt(bH1, x);

%skip stright to adding zeros every other sample
y0(1:2:end) = 0;
y1(1:2:end) = 0;

y2 = fftfilt(bH0, y0);
y3 = fftfilt(bH1, y0);

%skip stright to adding zeros every other sample
y2(1:2:end) = 0;
y3(1:2:end) = 0;

% =============
% RECONSTRUCT
% y0 = 2*fftfilt(bF0, y2) + 2*fftfilt(bF1, y3);
z = 2*fftfilt(bF0, y0) + 2*fftfilt(bF1, y1);

figure; freqz(z,1)
figure; plot(z)

If I try to go deeper than the first level of the tree (by uncommenting the line towards the end of the snippet), I no longer get perfect reconstruction. This makes sense to me as the signal down the left-side of the tree has gone through more processing stages, so it has greater delay relative to the signals processed by the right-side of the tree. I'm obviously not using the correct type of filter, or am fundamentally misunderstanding the examples available online. A nudge in the right direction, or basic code example would be great help.


texts: When Two Sampled Sinusoidal Are Orthogonal?

texts:
If two analog sinusoidal are orthogonal with duration T, their minimum frequency difference should be 1/2T in case of no phase offset between them.

And if there is phase offset, the difference is 1/T.

And in N sampled domain the orthogonality depends on sampling frequency.

The condition is $f_{s}=Nf_{k}/k$ obtained from $f_{k}=kf_{s}/N$.

($f_{s}$:sampling frequency, $f_{k}$: frequency difference, $k$: integer bigger than 0, $N$:# of sample)

Ex.)$N=2048, k=1, f_{k}=15kHz$ then $f_{s}=30.72MHz$

I like to know how $f_{k}=kf_{s}/N$ is derived

Thanks


texts:
I read from material that coherent integration can be provided by DFT. But it doesn't give further explanation on this remark. As far as I know, coherent integration is the summation of discrete samples with phase considered and it is usually used in time domain. How to understand that DFT's capability of providing coherent integration?

Thanks.


texts: I have just started learning some DSP, I keep seeing things like memory in regard to filters.  I am unsure exactly what this means and why it is useful.  I understand that increasing the number of taps increases the memory etc.. and a tap is essentially just a coefficient for changing the signal.


texts:
I am a beginner again to signal processing, I did my degree over ten years ago, and I am relearning all my signal processing stuff.

I am a software developer by trade but now I am looking at some LTE PHY signal processing.

I am trying to understand how QPSK signals get recovered in LTE. As I understand things, the I/Q signal is passed through a low pass filter, which results in a sort of NRZ complex waveform.

However if I try to get a constellation diagram of this data it is a mess. So how do I recover this LPF baseband to see the constellation diagram?

All I get by Google is 'NRZ' or 'Zero-forcing' or 'DC offset'. I am awaiting my LTE matlab toolbox, but I should be able to recover the I/Q right? Do I need to do some sort of convolution?


texts: Granulometry by morphological opening and closing

texts:
I am trying to understand granulometry. I understand that it deals with determining the size distribution of particles in an image and uses the principles of opening and closing.

I have a smoothed image with marbles in different sizes and call it A and a structuring element B, which is a disk with variable radii. I get that when opening the image with the structuring element B, smaller marbles disappear and by increasing the radius of B, bigger marbles start to disappear.

What I do not understand is what will happen if I use closing instead of opening. Do the bigger marbles disappear with a smaller radius B? So, does the exact opposite happen with closing?


texts: Eigenvalues of correlation matrix which have the form of an harmonic function

texts:
As a continuation to this question, I took the matrix $C_{2 \times 2}$ which is:

$$C=\left[
    \begin{array}{}
      a& ace^{-\frac{|\phi_1-\phi_2|}{2}}\\
      ace^{-\frac{|\phi_1-\phi_2|}{2}} & a
    \end{array}
\right] $$

And I found its eigenvector and eigenvalues which are:

$$v_1 = \left[
    \begin{array}{}
      1 \\
      1
    \end{array}
\right]
,v_2 = \left[
    \begin{array}{}
      -1 \\
      1
    \end{array}
\right] \\ \lambda_1 = a+ace^{-\frac{|\phi_1-\phi_2|}{\rho}},\lambda_2 = a-ace^{-\frac{|\phi_1-\phi_2|}{\rho}}
 $$

So, it is still not clear to me why the author says that the eigenvectors have the form $e^{in\phi_i}$.

This form doesn't even have the shape of a vector, and I specifically found the eigenvectors, and there isn't any complex component in it.


texts: Converting cascaded biquads to parallel biquads


texts: How to prove a filter has linear phase response?

texts:
I have designed an FIR filter to have linear phase response using odd-symmetry design. The coefficients of this filter are {2,1,3,1,0,-1,-3,-1,-2}. I am now being asked to prove it has linear phase response.

I don't think saying odd symmetry design is deemed a sufficient answer.


texts:
I was reading this document and it shows the computation of the magnitude and phase of $h[n]=-\delta[n]$.

We can get the DTFT as:
$$H(e^{j\omega})= -1$$

So the magnitude will be $1$, and according to my understanding, the phase is calculated as $$\tan^{-1}\left(\frac{\Im\{H\}}{\Re\{H\}}\right) = \tan^{-1}\left(\frac{-0}{1}\right)=0$$

But in the document, the phase is $\pi$. Where is my mistake here?


texts: Receiver function, frequency domain deconvolution not giving logic results

texts:
I'm working on some code for receiver function method in seismology. For anyone one not into the topic, it's just a deconvolution of two time series (seismograms). This can be done in the time domain as well as in the frequency domain, however, the method I'm working on is doing it in the latter. The function I'm applying is this one, as described in Langston, 1979:

The phi term is what we call water level, and what it does is that it "fills" the low frequencies in the DFT so it can ignore the noise and smooth the result. The Gaussian filter is also intended for smoothing of the result.
Both R and V seismograms looks like this, an usual seismogram:

However, after applying the deconvolution I get something like this:

I'm very confused by those two high spikes at the beginning and the end of the result, it should be smooth and not spiky.
This is my code:

verfrec = sp.fft(stv.data, numpoints)[:numpoints//2]
radfrec = sp.fft(stradial.data, numpoints)[:numpoints//2]
frec = np.linspace(1, numpoints//2, numpoints//2)

gauss = np.exp((-frec**2)/((4*60)**2))

waterlevel = np.abs(verfrec * np.conjugate(verfrec))
waterlevel = np.maximum(waterlevel, 0.1*np.max(waterlevel))

erefrec = radfrec * gauss * np.conjugate(verfrec) / waterlevel

inverse = sp.ifft(erefrec)

fig = plt.figure(figsize=(12, 4))
ax = fig.add_subplot(1, 1, 1)
ax.plot(frec, inverse, 'sandybrown')
plt.title('Sismograma')
fig.show()

I'm splitting the transforms range in the middle since the real input will produce symmetric outputs. The 60 in the gauss definition is the gauss parameter and it should be between 0.1 and 1.5 or so, with those values I get wrong results so I'm guessing its a magnitude thing.
Any help is very appreciated.


texts: select meaningful part out of a signal (for goertzel)

texts:
I have an on going input signal, most of the time with no activity.
Once in a while there is an activity (oscillations in different freq).
I need to take good 205 samples out of the signal.
what i mean by good?
really oscillating signal for 205 samples (not just crossing threshold).

here is an image to clarify:

bubble 1 - not good selection.
bubble 2 - good selection.

I need an efficient way with low computation complexity (for example STFT is heavy), to measure signal's sections oscillation level and chose the best one for further processing, someone knows?


texts: Calculating phase with noninteger periods of a sinusoid

texts: How do we calculate the phase of a signal with a noninteger period of a signal? For example, if a signal has a frequency of 1kHz and is sampled at 8kHz, 8 samples (1 period) can be used in the FFT or Goertzel to calculate the phase. However, if a signal has a frequency of 900Hz or 1.1kHz and is sampled at 8kHz, then 8 samples (1.11 or 0.909 periods respectively) will produce incorrect phase calculation. How do we compensate for a noninteger sample period?


texts: Bilateral $\mathcal Z$-transform of exponential

texts:
We all know that $a^nu(n)$ has unilateral $\mathcal Z$-transform. But what is the $\mathcal Z$-transform of $a^n$? (bilateral) When i tried to solve, i got answer as 'zero'.

But bilateral Laplace transform of $e^t$ doesn't exist. Both are exponentials  in discrete and continuous  domain respectively. Considering the similarity between Laplace and $\mathcal Z$-transform, how to explain the above problem?

Below, this is how I got 'zero'

$$a^n=a^nu(n) + a^nu(-n-1),$$

Now taking $\mathcal Z$-transform on both sides we get
$$\frac{z}{z-a}\quad \text{and}\quad \frac{-z}{z-a}$$ respectively which add to 'zero'


texts:
I have a relatively simple question. I know if you connect multiple LTI systems in parallel, overall system's frequency response is the sum of each LTI system's frequency response. Therefore, what would be the frequency response of the overall system if two FIR Low-Pass filters are connected in parallel? Let's say one has a cut-off frequency of 50 Hz and the other one has 100 Hz. I'm guessing it would have a frequency response of 1 between 50-100 Hz and 2 below 50 Hz.

Similarly, how would you construct a band-pass filter from a pair of low-pass filters? I am totally stumped with this one.


texts: Deciding for the correct way to calculate coherence time

texts:
I have to calculate the coherence time of an indoor environment. The signal I am sending is modulated with GFSK, the carrier frequency is about 2.4 GHz and my symbol time is 1 microsecond.

During my research, I had encountered some formulas like

$T_c = \frac1{f_d}$

$T_c = \frac1{B_d} $(doppler spread)

$T_c = \frac1{2f_d} $(doppler frequency)

$T_c =\frac{9}{16f_d}$

$T_c = \frac{0.423}{f_d}$

$T_c = \frac1{2\pi f_d}$

So I am a little bit confused here. Which formula suits best for my case and how to decide the right formula in different cases?

P.S.: The Doppler shift shouldn't be that much since it is an indoor experiment and the movement with the highest velocity should be walking quickly.


texts:
I want to design a filter to remove all frequency components from a .wav file except those within $\pm 25\,\text{Hz}$ of $523\,\text{Hz}$ as well as its harmonics (up to the Nyquist frequency).

Attempt

I tried to suppress frequencies below 523-25 Hz using a high-pass filter, and suppress the ones higher than 523+25 Hz using a low pass filter. Here's my code:

[s, Fs] = wavread('x.wav');

for N=2:1:(Fs/(2*(523+25)))
    % High-pass filter
    FsNorm =  (523-25).*N / (Fs/2);
    [b,a] = butter(10, FsNorm, 'high');
    sHigh = filtfilt(b, a, s);

    % Low pass
    FsNorm = (523+25).*N / (Fs/2);
    [b,a] = butter(10, FsNorm, 'low');
    sLow = filtfilt(b, a, s);
end

wavwrite(sLow, Fs, 'x_filtered');

The code doesn't work, the resulting file is mute. My intention is to apply each filter to the output of the previous filter, but I don't know how to implement this in code.

I used $(523+25) N_{max} < \frac{F_s}{2}$ to include all the harmonics up to the Nyquist frequency.

So what do I need to do to get this filter to work?

Edit:

Here's my new code:

[s, Fs] = wavread('x.wav');

sNew = zeros(size(s));

for N=2:1:(Fs/(2*(523+25)))
    % High-pass filter
    FsNorm =  (523-25).*N / (Fs/2);
    [b,a] = butter(10, FsNorm, 'high');
    sHigh = filter(b, a, s);

    % Low pass
    FsNorm = (523+25).*N / (Fs/2);
    [b,a] = butter(10, FsNorm, 'low');
    sLow = filter(b, a, s);

    Now sum the results
    sNew = sHigh + sLow;
end

sNew

wavplay(sNew,Fs);
wavwrite(sNew, Fs, 'x_filtered');

I still get sNew=0 and the sound is still mute. What's wrong?


texts: How to discard coefficients in the frequency domain without losing too much quality?

texts:
I was reading a text that shows that one can throw away 66% of the coefficients and gets an image that is still acceptable.

I tried to replicate it using python (last example) and selectively throwing away some of the coefficients but I couldn't get a great image, mine is looking like a ghost.

How can I discard that not so important to the image?

import matplotlib.pyplot as plt
import matplotlib.image as mpimg
import numpy as np
from scipy import signal

# loading the image
img = mpimg.imread('i/super_mario_head.png')

# image in frequency domain
fft_original = np.fft.rfft2(img)

# discarding some coefficients (except initial square 80x80)
for x in range(300):
    for y in range(300):
        discard_x = x > 80
        discard_y = y > 80

        if discard_x and discard_y:
            fft_original[x,y] = 0

# reverse from frequency to spatial
i_fft_original = np.fft.irfft2(fft_original)

f, (plt1, plt2) = plt.subplots(1, 2, figsize=(15,15))

plt1.axis('off');plt1.set_title('Original');plt1.imshow(img[:,:,0], cmap='gray');
plt2.axis('off');plt2.set_title('Reversed');plt2.imshow(i_fft_original[:,:,0], cmap='gray');


texts: Does perceived ADC resolution, range and precision depend on fixed point format? What is the best way to deal with this for signal processing?

texts:
I am working on a system that is implementing a digital filter that was designed for another system that had an FPU and a 12bit ADC. The system that I am working on is using the same filter without an FPU using a 12 ADC and fixed point arithmetic.

The ADC resolution is 3.3v / 2^12 = 0.000805664062v

With an FPU each ADC integer value could be scaled with the above resolution.

On the system using fixed point where

m2^-e, e = 10

The resolution would be

2^-10 = 0.0009765625

This leads me to think that the use of fixed point is giving me slightly lower resolution.

Looking at it differently and representing a delta between ADC integers with 14bits for the fraction the resolution could be

2^(-11) + 2^(-12) + 2^(-14) = 0.000793457031

but then the range would be perceived as less than 3.3V

Is this thinking correct?


texts:
This is an image processing question.

In 2-D Euclidean space a pillbox function, $f(x)$ is given by:

$$f(x,y) = 1\text{  if } x^2 + y^2 < R, 0 \text{ otherwise } $$

$R$ is come constant.

If I try to approximate the top hat when I'm looking at an image matrix of points, a naive algorithm would make some of the image pixels 1 inside a given radius and the others 0 outside that radius.

However, if I'm looking at a 4 x 4 grid of points, the shape I'd approximate would look like a square if I tried to naively create a discrete top-hat.

How would I create a discrete top-hat function that accurately approximates it's continuous counterpart, even when my image matrix is small?

EDIT:  I believe the comment about the edge image pixels having a value equal to the area of the pixel contained within the pillbox is correct. However, implementing this algorithm seems difficult. A reference would be appreciated.

EDIT 2: Matlab has an algorithm for computing a discrete pillbox function in it's fspecial function (disk option). This is the algorithm I'm after, though I still need a reference or an explanation of the geometry, as the function is quite complicated.


texts: Compressive Sensing - Sparse in frequency example

texts:
I am learning about compressed sensing. I have a question regarding a common MATLAB "sparse in frequency" example that can be find online, for example here and here.

What confuses me in these examples is that they do not take the subset of time samples from the time signal directly, but by inverse transforming the spectrum of the signal. Is this correct? Does this not already involve knowing the frequency spectrum of the signal ? I thought that the goal was to reconstruct the time signal from a subset of time samples, without knowing the exact shape of the spectrum. Or am I missing something?


texts:
I'm writing software and the data reporting frequency is slower than the acquisition cycle, so effectively sub-sampling is being done.  Different data types require different methods of sub-sampling which I'm tempted to call a Gating method where oversampled data in one time domain goes through a Gate to another time domain. Windowing also comes to mind.

I hesitate though as although these terms imply an interval of focus/ interest, they might not necessarily imply continuity between intervals; back to back intervals where no input is ignored.

Such methods might be:

Average/Median since last report (or low pass filter)
Min/Max since last report
Sticky flags (a bit field who's bits are cleared when read)

Is there a better term than Gate or Window?

Note 1: This is borderline a programming question or code review, but I think still more for signal processing.

Note 2: Although sub-sampling causes signal artifacts such as aliasing.  This is not part of the question.


texts: I'm doing a research on the FFT method, and a term that always comes up is "frequency bin". From what I understand, this has something to do with the band created around the frequency of a given sinusoid, but I can't really figure out how. I also figured out how to go from a given bin to its related frequency, but still no intuition on what a frequency bin is.


texts: Integer and a half periods of sine wave in DFT window

texts: Why does a pure sine wave starting and ending at 0 in the sample window of a DFT but containing 1.5 cycles (or periods) produce the same spectral leakage (i.e. not more) as a pure sine wave with integer number of periods in the sample window of a DFT.


texts: Why do frequencies of analog signals range from $-\infty$ to $\infty$ while frequencies of digital signals are restricted to $[0,2\pi]$?

texts: In Fourier analysis while dealing with discrete-time signals, frequencies range from $0$ to $2\pi$ why? Intuitively how can i understand it?


texts: Neural Network learning project based on 8 wave signals over 1 second at 1 sample every 10 ms ( hence 100Hz )

texts:
I'm currentely trying to train a neural network that can decide wether a pattern produced by the movement of a hand near capacitive sensors is as expected, or random.
I have an MPR121 microchip linked with an arduino, providing me 8 signals ( ranging from 0 to 255 ), that stays at baseline ( around 135 for my current conditions ) when nothing is near or perturbating the conductances. The value lowers when something conductive approaches it. Here is some pictures to help you represent the thing.

The pattern of the eight channels together give information on the postion of the hand on one axis. The neural network is supposed to learn himself how the different channels react, in wich order, so i don't have to tell anything to the programm concerning the physical distance between two electrodes or whatever.

The fact is that i would like to learn how to do neural networks but i have no idea how to feed the data ( especially for data over time, which differs from typical image analysis tutorials often seen for neural networks over the internet ), considering that there is 100*8 values to analyse for each frame (100 values for 1 sec frame, for the eight channels ) and i don't know if i should have a network with 800 input nodes.... seems to me like i miss something...?
Can i pass a list as one input to a neural network ,hence reducing the amount of input nodes to 8 ?

And another question i had, once trained, will th algorithm be exportable on an arduino ? I read somewhere that once trained, the weight of the nodes can be gathered, alowing to make an equation to choose from inputs to output, here yes or no. Arduino should be able to perform the calculation of this equation quite complex in a short amount of time i guess ?

( If you want to know, the labelisation of the data will be performed afterward, with a button connected to the arduino, telling him to record an d labellise as positive on an SD card the past 1 sec recording of the 8 channels.
Everytime i will perform the hand pattern to be learnt, i will push this button. The randoms signals wich are labelled negative are 1 sec recordings, picked randomly over time of the day, to catch changes in baseline and hopefully some false positive noises due to people wandering around )

Thank you in advance for your help and advices on my project. =)
And by the way, sorry for my english, hope i didn't shocked anyone.


texts: Choosing the best feature extraction algorithm for die photographs without scaling and rotation

texts:
I'm working on a live stitching/"panorama" program meant to recreate the die (the wafer part, in the context of integrated circuits) of silicon chips from microscope photograph.
The first step is to use a good feature extraction algorithm adapted to my needs:

no scaling differences for the same "panorama"
very low rotations possible, generally no rotation
no perspective
almost no light change

The images are very similar in structure, like this one (20x objective)

Or this one (5x objective)

The images can be cropped automatically to remove the blurry part on the sides.
Do you have any suggestions for a fast existing or custom algorithm adapted to this kind of images?
Or more generally for images without any rotation and scaling?


texts: Check for time invariance of a _system_ using MATLAB?

texts:
I am having difficulty in developing an algorithm for checking if a given system is time invariant or not.

For example, let there be 2 systems with input output relation given as:

$$y(t)=x(t-2)$$ or $$ y(t)=tx(t-1)$$

Now the first system is time invariant whereas the second one is not. But how do I write a function in MATLAB that would take the input output relationship and find out if the system is time invariant or not.


texts: Finding a region where a function is periodic

texts:
I am looking for a method to find the region where a function is periodic. So, for example, I have a function that is periodic in some places and not in others, I would want to know the range of x values where it is periodic.
 I'm trying to do this in matlab and I tried to use autocorrelation but I can't use it to find the original range where the signal is periodic


texts: Root mean square bandwidth of a signal with flat spectrum

texts:
I am trying to derive a simple expression for root-mean-square bandwidth $B_{\rm rms}$ for a signal that has flat spectral density over its entire bandwidth $B$. The expression for the $B_{\rm rms}$ is

$$
B_{\rm rms} = \sqrt{\frac{\displaystyle\int_{-B/2}^{B/2} \lvert H(\omega)\rvert^2 \omega^2\, \mathrm{d}\omega}{\displaystyle\int_{-B/2}^{B/2} \lvert H(\omega)\rvert^2\, \mathrm{d}\omega}}
$$

where $\lvert H(\omega)\rvert^2 = A$ is the constant PSD over $B_\rm h$.

I read it at some sources that this gets simplified as follows:
$B_{\rm rms} = \frac{B}{\sqrt{3}}$.

However, when I integrate this I do not see how I can get this expression:

$$
B_{\rm rms} = \sqrt{\frac{\displaystyle\int_{-B/2}^{B/2} A \omega^2\, \mathrm{d}\omega}{\displaystyle\int\limits_{-B/2}^{B/2} A\, \mathrm{d}\omega}} = \sqrt{\frac{\displaystyle\frac{1}{3}\cdot\frac{2B^3}{8}}{\displaystyle\frac{2B}{2}}} = \frac{B}{\sqrt{12}}
.$$


texts: Why is my fundamental frequency at the 15th harmonic order

texts:
I am trying to implement a code in MATLAB to show harmonic content and order in Matlab. my fundamental is supposed to be at harmonic order 1, however the ouput plot I am getting is wrong and I am having difficulies proceeding.
I am a first time user of this forum and i have posted the same question in the matlab site
Here :
I would really apreciate any help. Thanks

[D,S,R] = xlsread('test data.xls');
v = D(:,2);
Signal = D(:,2);
Ts = 0.00005;                                               % Sampling Interval (seconds)
Fs = 1/Ts;                                                  % Sampling Frequency (Hz)
Fn = Fs/2;
%PERFOM FFT
N = length(Signal);
meanSignal = mean(Signal);                                  % ‘Signal’ Mean
FTSignal = fft(Signal-meanSignal)/N;
Fv = linspace(0, 1, fix(numel(FTSignal)/2)+1)*Fn;           % Frequency Vector
Iv = 1:numel(Fv);                                           % Index Vector
%Plotting Harmonics (THD)
har_mag = abs(FTSignal(Iv))*2;
figure;
bar(har_mag(1:20),0.4)
xlabel('Harmonic Order','fontsize',10);
ylabel('Voltage (V)','fontsize',10)
hold off

False reading of the output I am getting is this:

the output should look similar to  figure below:


texts:
I am new to seismic data processing and I really have no understanding of the term 'phase'. So, could anybody give me a simple explaination of the term?
Thank you


texts: Designing A Filter Using SPTool In Matlab

texts:
I am currently working on learning DSP and for my course we use Matlab. I have been tasked with designing a low-pass filter using SPTool with the following specification:

  A Low-Pass Kaiser Window FIR filter with a cut-off at 200hz.

Now, when i open the SPTool and choose "new" i get the filter design and analysis window. I use the following settings

Response Type = "Lowpass"
Design Method = "FIR" and "Window"
Window = "Kaiser"

I would like to know which values i should put into the "Fs" and "FPass" fields. The signal i am applying this filter to is a periodic 100hz wave.

Thanks for any feedback.


texts: In computed tomography (CT), why is 'Inverse Problem of Radon transform' studied?

texts:
As I know, the Radon transform is a very important tool in CT by Beer's law. Thus, finding $f(x)$ of Radon transform $Rf(L):=\int_{L} f(x)dl(x)$ is helpful in CT. Nowaday, the Filtered back-projection formula (FBP) for all dimensions $n$ (which can find the $f(x)$) is obtained and is used already.
My question is here. Why is people studying 'Inverse Problem' in Radon transform although FBP is obtained?


texts:
I'm writing a program in which you can synthesize waves by adding to a sound's Fourier transform, and then inverse the transform to get the modified sound. In order to do this, I need to know what to add to the DFT to synthesize a pure wave. I've tried to learn about Fourier synthesis from many sources, but they all talk about the Fourier series instead of the Fourier transform, and they all say that for a pure wave all you need is a value in the coefficient of that wave's frequency, and 0 everywhere else.
But the Fourier transform being a continuous function. And you can see here what the Fourier transform of a pure 12KHz sine wave looks like (on a logarithmic scale). As you can see, it's not just an instantaneous peak and 0 everywhere else.
So I tried to do the math myself. To put it formally, say we have a sampled signal $S_t=A cos(\frac{-2 \pi k}{N} t + \phi)$, where , $N$ is the number of samples in our signal, $0 \le t < N$, $A$ is an amplitude  $0 \le k < \frac{N}{2}$ determines the wave's frequency, and $\phi$ is the phase of the wave. The DFT of this signal would be a sequence $F_0,...,F_{N-1}$ where:
$$ F_r = \sum_{t=0}^{N-1}{{S_t}e^{\frac{-2 \pi i t r}{N}}} $$
My hope was that since our signal is a pure wave, there will be a direct formula for calculating $F_r$ that can be computed much quicker than by directly doing the math above, and that that formula will be the shape of that spike in the image I linked, because that's what I'm really interested in. For the sake of simplicity, I assumed that $\phi=0$, since once I figure out how to solve that case I'll probably be able to generalize it. I won't bore you with all the math because it's not important and you probably already know the answer. What I got in the end is that $F_k = F_{N-k} = \frac{NA}{2}$, and for all other $0 \le r < N, F_r = 0$.
Shocker. So in the end I found that it does equal to 0 everywhere except the frequency of the wave. But then I don't understand where that spike shape from the graph I linked is coming from. I'm still after it though, because just adding to the one sample which corresponds to the frequency I'm interested in doesn't work for me (I've tried).
So I guess my question is: how come plotting the Fourier transform of a pure wave doesn't actually produce zero in all samples except the one corresponding to that wave's frequency, and what's the formula for that spike that you do see.


texts:
Are Finite Impulse Response filters always linear? I looked at https://en.wikipedia.org/wiki/Nonlinear_filter and its says that if a filter does not obey linearity, it is non-linear filter. Thats understood? But does does linearity has to do with the length of impulse response?

What am I missing here?


texts:
I have an irregular shaped 3D image that I want to do frequency domain filtering on. Currently I simply zero-pad the data to make it a cube, then do an fft.

I was wondering if this was a valid approach and what effects this may have. Are there any alternate (better ways) to do this, and if so could you point me in the direction of some literature?


texts:
I try to detect some FSK signals made by audio in near real time. I know that I will get some latency.
In the start I try to detect some tones of 10Khz in a duration of 75ms.
That is achieved by a Goertzel filter.
What is the prefered technique for that issue?
I thought of making a little buffer that will always be detected by the Goertzel, but this sounds bad because I don't know when the the signal starts.
Can you give me some help or idea ? Thanks


texts: Gradient Domain Reconstruction - Scaling Problem

texts:
I am implementing reconstruction of image from gradient domain. This requires solving the following partial differential equation (a Poisson equation) on a 2D grid:

$$\nabla^{2}I=\mathbb{div} G$$

$\mathbb{div} G$ is divergence of a gradient field, which is simply an image with Laplacian operator applied (second-order derivative):

The original image have range of values 0 to 255 (8-bit image). The divergence image has values ranging from -341 to +318, which is expected.

However, when I solve the above equation for $I$, the result have very small scale. The values of reconstructed image range from -0.000038 to +0.000029.

I know the gradient domain stores information up to a constant. In this case however, the scale has been changed considerably.

When I stretch the reconstructed image to 0-255 interval, the correct result is obtained:

I am using Poisson solver based on the Full Multigrid Method (FMG) which works spectacularly in different application (HDR images). Here the values are stretched anyway so I never noticed the huge change of scale.

Is this scale change normal for gradient domain reconstruction or may there be

My ideas:

Some glitch in the numerical solution.
The $\mathbb{div} G$ image is padded with zeros so its size conforms with the solver - could this affect the scale?


texts: What principles to use, in order to create field that will be read by OCR?

texts:
I currently have an invoice, where some numbers are written by hand, and I would need to read them automatically.

My past experience is that people who write numbers will them write in so many different ways (overlapping with field name etc.) that it's really hard to isolate only digits.

Now I would like to redesign invoice, for example, to add field of some rectangular shape, that will be easily detectable (think easy to detect like a QR code is).

I'm not sure:

should I add some ordinary box, or are there any specific designs of fields/boxes/rectangles that would make finding the field and isolating everything in it (from image) simple and 99% correct?

Tech stack that would be used is:

OpenCV
Python
HTML/CSS for creating the design of invoice, that would be exported to pdf for users to fill it with the needed information


texts: How can a linear operator on DFT vector produce the same vector using only half of the DFT vector?

texts:
Suppose there is a DFT vector $\mathbf{X}$ with length N, which presents complex conjugate symmetry around its middle point, i.e., $X(1) = X(N-1)^*$, $X(2) = X(N - 2)^*$ and so forth.  $X(0)$ and $X(N/2)$ are the DC and Nyquist frequency respectively, therefore are real numbers. The remaining elements are complex.

Now, suppose there is a matrix $\mathbf{T}$, with size $N \times N$, which multiplies vector X.

\begin{align}
\mathbf{Y} = \mathbf{T}\mathbf{X}
\end{align}

Assuming the operator $\mathbf{T}$ preserves the complex conjugate symmetry in $\mathbf{Y}$ (see Conditions for precoding matrix to preserve complex conjugate symmetry on DFT vector for more details), how can I obtain the same "positive frequencies" (from $X(0)$ to $X(N/2)$), using an $(N/2 + 1) \times (N/2 +1)$ matrix $\mathbf{T}$, instead of an $N \times N$?

The motivation for this question is to reduce complexity in the operation.

Consider that $X$ is a DMT symbol and its elements are mapped from QAM constellations with different sizes. Consider also that only the positive frequencies carry useful information, the "negative" frequencies are only used for the IFFT to be real.

EDIT:
Is there a low-complexity solution other than using the $N/2+1$ upper rows of $\mathbf{T}$ multiplying the length-$N$ vector $\mathbf{X}$? This reduces the complexity to $N(N/2+1)$ CMACS plus the operations required to assemble the hermitian symmetric symbol

EDIT(2):
I suggest the following code to illustrate the problem:

N = 8;
A = rand(N,N); %must be real-valued
w = exp(-1j*2*pi/N); % twiddle factor
W = w.^(repmat(0:N-1,N,1).*repmat(0:N-1,N,1).'); % DFT matrix
T = W*A*W' % Linear operator
x =  fft(rand(8,1));

It is obviously clear that:
   T(1:N/2+1,:)x
Returns the positive frequencies from:
   Tx
However,

Using x(1:N/2+1), what linear operator (with size $(N/2 + 1) \times (N/2 + 1)$) should I use to obtain the positive frequencies?


texts: How to Prevent overflow/saturation in fixed point implementation of Fast Data Projection Algorithm

texts:
I am trying to implement an algorithm in real-time on a Fixed point DSP (The Blackfin from Analog Devices). The algorithm does a lot of stuff, but in the middle it performs an algorithm called "Fast Data Projection Method" (FDPM), which goes something like this:

Let $\mathbf{x}_{k} = [x_{1}, x_{2},\dotsb,x_{M}]^{T}$ be a random vector which contains $M$ samples from a discrete signal. In the Process we take sequentially many vectors from the sampled signal with some degree of overlapping, but that is not relevant.  We can assume that the vectors $\mathbf{x}_{k}$ come one after another.

The FDPM aims at obtaining a matrix $W \in \mathbb{R}^{M\times N}$ whose columns are the eigenvectors of the correlation matrix $R_{x} = E[\mathbf{x}\mathbf{x}^{T}]$ of the vector $\mathbf{x}$. So we initiallize the algorithm with a random matrix $W_{0}$ (Can also be the identity matrix), and perform the following steps:

1 - $\mathbf{y}_{k} = W_{k}^{T}\mathbf{x}_{k}$

2 - $\mathbf{a}_{k} = \mathbf{y}_{k} - \|\mathbf{y}_{k}\|\mathbf{e}_{1}$,  (Where $\mathbf{e}_{1} = [1,0,0,\dotsb,0]^{T}$)

3 - $G_{k+1} = I-\frac{2}{\|\mathbf{a}_{k}\|^{2}}\mathbf{a}_{k}\mathbf{a}_{k}^{T}$

4 - $W_{k+1} = Normalize\{\left[W_{k}+\mu_{k}\mathbf{x}_{k}\mathbf{x}_{k}^{T}W_{k}\right]G_{k+1}\}$

Where $Normalize\{·\}$ stand for normalizing each of the columns of the matriz individually.

So the problem is the following. I have tested this algoritm in MATLAB with double precisiòn and it works fine, i am able to obtain the eigenvectors very close to the real signal eigenvectors and the algorithms does it's job, no problem.

The thing is, when i use the fixed point toolbox to perform this operations with 16bit numbers in Q15 format, then everything goes wrong. The first one or two iterations have small error, but quickly all the vectors start to saturate (The components go to either 1 or -1) and while the matrices do not saturate, they converge to a matrix which is orthogonal to the same matrix calculated with double precision.

I am guessing that the loss of precision occurred by casting to 16bit affects too much, so i think i might have to do some signal scaling or some "re normalization" every once in a while, but i have no idea how to analize the problem and how to know where i have to modify the algorithm to make it work in 16bit.

Has anyone ever worked with this type of problem? anyone has any idea what i could do?

Thanks!!


texts:
Split from thread: Correcting Impulse Invariance Method

Back to this subject. Matt L suggested combining IIM and MZT b-coeffs to get the filter response improved at hihger frequency area.

Couldn't find suitable IIM (C++) implementation to use so, I'm still having some difficulties with calculating and combining coefficients. How's that done?


texts: Computing the derivative from a set of real-time data?

texts:
So I have an array of say 25 samples and I would want to be able to note the trends of whether it's decreasing n or increasing from those 25 sample time interval(basically 25 samples array is my buffer that is being filled by every say 1 ms). Note that it is general trend that I am looking for, not the individual derivative(as I would have obtained using finite difference or other numerical differentiation techniques). Basically I expect my data to be noisy so there might be ups and downs even after doing filtering and so on. But it's the general trend of increasing or decreasing behaviour that I am looking for. I want to integrate the increasing/decreasing behaviour in every ms to trigger some event which is more of a user interface event (blinking a LED) so it does not have to very delay of processing as long as I can detect the general trend.

Thanks in advance!


texts: What does the variable $\tau$ mean in frequency modulation?

texts:
Frequency modulation formula:
$$
x_{FM}(t) = A_{c}\cos\left(2\pi f_{c} t+ 2\pi k_{f}\displaystyle\int_0^{t}m(\tau)d\tau+\phi_{0}\right)
$$
Its change from $m(t)$ to $m(\tau)$.
What is $\tau$? I wonder what it means.


texts: Does FFT on two samples from the same audio but with a very small sampling offset return the same results?

texts: Using a single audio source, suppose I take two samples of the same length (say, 10s), with the same sample rate (say, 44100Hz), one just 1/3 sample length (7 microseconds) away from the other. Does FFT on the two samples return (roughly) the same shape? Is there a way for me to sync the two samples?


texts: 1st and 2nd Gaussian derivatives for edge detection

texts:
I have a project about edge detection, I made research on internet about filters and read many articles but there is an article talked about using Sobel and Gaussian derivatives then it mentioned "edge" and "ridge" detection successively. I found the equation for 1st and 2nd Gaussian derivatives ($3\times 3$ and $\sigma= 0.8$) and applied the before mentioned method but I noticed blurred image with slightly enhanced edges.
I still don't get the idea nor the exact role of using different Gaussian derivatives together for better edge detection, and how can I enhance curvilinear structure detection to get cleaner edges. If someone can break it down it would be appreciated.


texts:
Using Matlab's FFT function, I am getting different Amplitude values in a particular frequeny, when I use different range of Bandpass filters centering around that particular frequency. What's wrong here? How can a particular Amplitude scale be selected?

In the following code, I calculated two FFTs of the same signal over the frequency range of (59Hz to 61Hz) and (59.98Hz to 60.02Hz) after bandpassing them in the said range.

[x_noise,fs] = audioread(h);
x_noise = x_noise(1:length(x_noise)/360);
xlen = length(x_noise);
x_main=x_noise;

i=1;
%%
%Noise
f1=59;
f2=61;
fn=fs/2;
w1=f1/fn;
w2=f2/fn;

[b,a]=butter(5,[w1 w2],'bandpass');
xn_har = filter(b,a,x_noise);

%FFT of Noise
X_Noise = abs(fft(xn_har));
 p = linspace(0,1000,length(X_Noise));
 figure(1)
 plot(p,X_Noise)

%%
%Signal
f1_sig=59.98;
f2_sig=60.02;
w1_sig=f1_sig/fn;
w2_sig =f2_sig/fn;

[b,a]=butter(5,[w1_sig w2_sig],'bandpass');
x_sig_har = filter(b,a,x_main);

%FFT of Signal
X_sig = abs(fft(x_sig_har));
q = linspace(0,1000,length(X_sig));
figure(2)
plot(q,X_sig)

But after doing FFT, the Amplitude in 60Hz seems to vary a lot. As seen in the image, one has a value of around 1600, while the other has around 1.2. Why is this happening, they are at the same frequency?


texts:
I have an application that fills a buffer from my computer's microphone and does an FFT. My objective is to figure out what frequencies are in the sample.

I've tested the application by whistling as well as by playing this video from my cellphone into my computer's mic (I've tried two different mics to be sure). What I expect (taking the absolute frequency) is a smooth Gaussian-like curve with a mean at the frequency of the incoming tone, but what I'm getting is this:

What is the cause for the jaggedness of the signal's frequency spectrum? And what technique can be used to reliably identify the main frequencies present in the signal while filtering out the unwanted bumps?


texts: EEG signal processing with sampling frequency

texts:
What is the relation between sampling frequency  and the frequency of the signal after sampling? For example sampling frequency is 100 Hz,
that is the signal contain maximum frequency of 50 Hz. Then what will be the maximum frequency after sampling? Is it 100 or 50?
Or it is above it? How can we see it in MATLAB?


texts: Adding phase noise to a signal results in amplitude variation on a constellation - whats the mistake? *plots inside*

texts:
I am attempting to add phase noise to a modulated signal.

I am doing this by taking a phase noise mask (single sided 0 Hz to 1
MHz), which shows dBc values away from central value or reference
value. Phase noise is decreasing as we move in frequency away from
the reference point which is 0 Hz. [Top left image in the 4-plot
image below]
I then make this double sided.
I convolve this phase noise plot with a spectrum of a clean complex
exponential at 1 MHz, this is my oscillator for converting frequency.
I think this has imparted the phase noise onto the complex
exponential.  [Top right image in the 4-plot image below]

The spectrums are convolved  to make a oscillator with phase noise:
$\mathcal{F}(e^{j \omega t}) * \mathcal{F}( e^{j \phi(t) }) = \mathcal{F}(e^{j \omega t})* \mathcal{F}(e^{j B}) = \mathcal{F}(e^{j( \omega t + \phi(t))}) $

where $\mathcal{F}$ is fourier transform, $ \omega$ = 1 MHz, $\mathcal{F}( e^{j \phi(t) }) $ is the double sided phase noise mask on 0 Hz, then also this can be interpreted as a bandwidth $B$ of the upper and lower sidebands of the phase noise centered on 0 Hz (the double sided spectrum in the top left plot below when not in logarithmic x-axis scale). The convolution result is the top right image, symmetric phase noise translated from 0 Hz onto my a oscillator or CW or sine wave at 1 MHz.

I take the Inverse Fourier Transform of this convolution to produce
an oscillator in the time domain, my thought is that the phase noise
is on it.
I then multiply my complex modulated signal (16-APSK) by this oscillator (complex exponential with phase noise) to produce an up conversion and then downconvert with a clean complex exponential
with no phase noise. This should leave the phase noise on the
modulated signal.

The math now being done in time domain:
Up-conversion using the inverse transform of the result above
$A(m) e^{ j( \phi (m))}  e^{j( \omega t + \phi(t))} $
Down-conversion of the upconversion
$ A(m) e^{ j( \phi (m))}  e^{j( \omega t + \phi(t))} e^{-i \omega t}  $
where $A(m) e^{ j \phi (m)} $is a modulated signal

However I am showing results that are not like classical phase noise,
which would be rotation acrcs around constellation points but it
appears as an amplitude modulation?  [Final constellation image]

I know this is wrong, but I stuck why. If I subtract in the time
domain the clean sinusoid without phase noise and the oscillator with
phase noise, I get this variation across the full time domain
simulation. [Central images are the real and imag values of the oscillators in time and the bottom images are the difference]

Why is this producing amplitude variation and not phase variation?
Is there another way to do this?

--Edit Updated constellation based on replies below--


texts: Can harmonics be created or simulated a la carte

texts:
I am trying to learn about adding harmonics (audio) with a goal to make a plugin that can allow adjustment of the amplitude of each harmonic being added. I would not have the ability to write this but wanted to get an opinion about whether it is possible.  Example: to add a 2nd harmonic and adjust its level separately, add a 3rd and adjust its level etc. Ideally have a GUI that allows the ability to add harmonics and turn a knob or slider to adjust the amplitude of each harmonic up to around the 10th. The goal is to find out from some someone if it is possible and if so I’d pursue someone to write this for me.
Since the input can be any type of Audio from 20-20k, filtering is not likely the idea. I understand to some limited degree the methods of adding harmonic distortion, but have not been encouraged towards a la carte control. The end result is to custom tailor a harmonic distortion effect.
Thanks