texts: How do I identify repeated subimages within a larger image?

texts:
I have a large image of a microscopic die of an integrated circuit. The pixel dimensions are something like 15000 x 10000.

Example image:
Commodore Amiga Paula Image http://tinytransistors.net/images/chips/Amiga_A500/Paula_total.jpg

The image contains a number of repeated structured subimages. Some large blocks might even repeat 3 or 4 times because there was a need for duplicated functionality. (See good candidates on the left of the example image)

I'd like to identify the largest(in pixel dimensions) repeating subimage block, and then highlight the locations of those subimages. I would repeat this process looking for the next largest repeating subimage block, highlight those, and so on.

To be clear, I don't have the subimages in advance. The problem isn't only to find the location of certain subimages within an image. That's Step two. :)

While the image is in fact a color image, it could be easily turned into a gray-scale image without loss of content.

There is obviously some "noise" to an image like this, and there'd have to be some adjustable "fudge factor" on what to consider a match.

I'm not entirely sure where to start on this, or if this is even the correct SE. Would FFT be a potential solution? Using MatLab? Specific one-off image analysis tools?

I could write some simple image analysis software that breaks an image up into smaller blocks, and then searches the image for those blocks. But this sounds incredibly slow to me, and there's got to be a more efficient way to do this.


texts: Proof of the convolution property of Fourier Series in continuous time

texts:
I am facing problem in understanding the proof of Convolution property of Fourier Series (FS) in continuous time CT;
that is:  $$\mathrm{FS} \big\{x_1(t)\star x_2(t)\big\}=T\sum_{n=-\infty}^{\infty}c_n d_ne^{j \Omega_0 n t}$$ where $x_1(t)$ and $x_2(t)$ are periodic with $T$, and $c_n$ and $d_n$ are the Fourier series coefficient of $x_1(t)$ and $x_2(t)$.

Below I am showing the proof and my confusion,

PROOF:
\begin{align}
x_1(t)\star x_2(t) &= \int_T x_1(\tau)x_2(t-\tau)d\tau\\
\therefore \mathrm{FS}\big\{x_1(t)\star x_2(t)\big\}&=\int_T FS[x_1(\tau)] FS[x_2(t-\tau)]d\tau\\ &=\int_T \sum_{n=-\infty}^{\infty}c_n e^{j \Omega_0 n \tau} \sum_{p=-\infty}^{\infty}d_p e^{j \Omega_0 p (t- \tau)} d\tau\\
&=\sum_{n=-\infty}^{\infty} \sum_{p=-\infty}^{\infty} c_n d_p \int_T e^{j \Omega_0 (n-p) \tau} e^{j \Omega_0 p t} d\tau\\
&=\sum_{n=-\infty}^{\infty} c_n d_n T e^{j \Omega_0 n t}
\end{align}
CONFUSION:
Now I am confused  how we are getting
$$ \sum_{n=-\infty}^{\infty} c_n d_n T e^{j \Omega_0 n t}  \quad\text{from}
\quad\sum_{n=-\infty}^{\infty} \sum_{p=-\infty}^{\infty} c_n d_p \int_T e^{j \Omega_0 (n-p) \tau} e^{j \Omega_0 p t} d\tau \quad ?$$

Please explain.


texts: Predicting Execution Time of DSP programs

texts:
I've been tasked with predicting the performance of our DSP code, for the obvious reasons (reduce time to market, add certainty in choosing chips, etc.)

Consider a system composed of many subsystems with timing results $P_{i}$ (in seconds) (filters, gains, dynamics processing etc.) with a clockrate $C$ (seconds/sample) and sampling rate $Fs$ (samples/second)

Our initial assumption is that:

$$\text{CPU utilization} = 100\% * \sum_{i}^{i=N}{P_{i}*Fs/C} \stackrel{!}{=}
\texttt{top -d 1}$$

i.e. CPU utilization estimate should equal the output of top -d 1`.

However, this doesn't hold up to real data! It's always underestimated.

We measure individual subsystems by wrapping each of them in a test harness and timing their execution. Our system runs in a single thread and on a single core.

Depending on the platform, we can get within 5 to 20%, but always underestimating top.

Why is this?

I can find nothing in the performance literature that lines up with our assumption, but the explanations are so numerous that it's almost like "we don't really know why we can't predict performance."

I am beginning to think that the solution is monitoring like a hawk during early development and scaling back the system complexity as we uncover new solutions.

However, it would be pretty awesome if we could design a system with the performance requirements ahead of time. Our requirements have largely been determined by "this is better than the old thing the vendor had, and what we had before worked with our old system."


texts: What is impulse response (or equivalently the frequency response) of the channel ?

texts: Please allow me to ask this silly question. I googled, but I still not clear about this. What is impulse response of a channel (wireless communication channel). Someone can explain me with a lot of pictures? I know this question is too ez, but I really need to understand it fully. Thank you again.


texts:
Which one of the following is the condition of aliasing?

(a) Tails of the replicas enter into the Nyquist
interval

(b) The tails of the replicas enter into the Nyquist
interval and add to the original spectrum, distorting it.

According to my understanding, for aliasing condition $(a)$ is enough and
hence the overlapping condition specified in condition $(b)$ (and hence
of the textbook) is not necesary.


texts: Does delta-sigma ADC also reduce Gaussian noise on input signal to ADC or just quantization noise?

texts:
The motive for posing this questions arises from a difference of analysis between a colleague an myself.

Our general environment is in the construction of an analog front-end which takes in signals to an instrumentation amp in the low microvolt range. followed by a gain stage which then feeds a 24 bit delta-sigma converter.

My understanding is that we can achieve a certain number of noise-free bits from the ADC, where the noise is the quantization noise. For any Gaussian noise on the input to ADC, this noise will not be reduced by the ADC itself. And in this regard, we need to do some additional averaging of the sampled data points for reducing the Gaussian noise as an additional processing step.

My colleague thinks that the delta-sigma ADC will also reduce the Gaussian noise on the input signal, and therefore additional averaging is not needed.

As a second question, is it the case that we should arrange the gain so that as much of the signal fills the dynamic range of the ADC so as to improve the SNR of the signal due to having a large signal quantized as a function of quantization noise. In other words, the quantization noise is a smaller fraction of the signal?

Could you help us out by commenting on this question?


texts: Is interpolation (interp1) better than FIR filtering when rational integers are close to 1?

texts:
Question
I've been attempting to resample a GPS signal in MATLAB. I've built a few FIR filters using fvatool and from handmade transfer functions (punched out with a HP-35s). Most are kaiser-windowed LPFs, but some are least-squared LPFs. All of the filters I've built have nominal responses that should adequately filter the aliasing out of the interpolations. I've also built a FFT-Resampling tool that resamples well and very quickly. However, I've hit a strange anomoly, at least as I see it. A simple 'nearest' neighbor interpolation using interp1 seems to much more accurately resample my signal than any FFT-Resampling I have done, or any upsample, filter, decimate tool I've used, to include dsp.FIRRateConverter, resample,upfirdn, and even straight upsample, filter,decimate. Why is this? Can anyone explain why this occurs? My theory is that when the rational integers of your sampling ratio are very close to one, straight interpolation will provide more accurate results because you reduce occilations and/or other noise distortion.

Background
My signal is a standard collected GPS signal, pulled from the air, sampled at 25e6 sps and saved in a binary file using signed 16-bit integers. When pulled into MATLAB (has to be done in chunks due to size) the signal is a complex row vector.

I'm upsampling to 26.25e6 sps, therefore my rational integers p and q are 21 and 20, respectively. To use interp1 I simply build two time vectors:

 t0 = (0:1:numSamples - 1)*1/sampleRate;
 t1 = (0:1:numResamples - 1)*1/resampleRate;

where

 sampleRate= 25e6;
 resampleRate = 26.25e6;
 samplingRatio = resampleRate/sampleRate;
 numSamples = length(s0) % length of my original signal, sampled at 25e6 samples
 numResamples = length(s0)*samplingRatio % length of my sampled signal; should equal length(s1)
 s1 = interp1(t0,s0,t1,'nearest','extrap');

When using a FIR filter, I have used many methods. But the simpliest is to call resample like this (resample is a least squares FIR filter with a kaiser window called through firls that gets pushed through upfirdn):

 s1 = resample(s0,p,q);

where

 [p,q] = rat(resampleRate/sampleRate,1e-12);


texts: Bilinear Interpolation Algorithm for up-sampling 2D images

texts:
In keras it is possible to use UpSampling2D layer to up-sample an image. You can use Bilinear Interpolation.
Given an image ${h\times w}$ it is possible to increase its size in ${h*k\times w*l}$, where $k$ and $l$ are factors greater than or equal to $1$.
I understand the Bilinear Interpolation formula applied to $4$ neighbouring pixels (${2\times 2}$ sub-matrix). But I just can't figure out how to apply it to perform up-sampling. What is the actual algorithm being implemented? How is Bilinear Interpolation mathematically defined for up-sampling 2D images?
Nowhere have I found the real algorithm, no one explains how to construct the new matrix correctly and consistently. It is always explained in an incomprehensible way, please help me.


texts: How many data points are required for one revolution?

texts:
I have a question about an FFT. I have recorded the vibration of a gearbox with a sampling rate of 48 kHz. The speed of the input is 50 Hz, the speed of the output is 12.5 Hz. At a sampling rate of 48 kHz I have 48000 data points per second and after the output rotates 12.5 times per second I have (48000/12.5) 3840 points per rotation.

Does it make a difference if I use only 3840 points for the FFT (i.e. 1 revolution) or if I use 48000 points and thus have 12.5 revolutions?

If someone has literature on this or can explain to me whether there are rules for this, I would be very grateful.


texts: Convolution output signal is very large, do I need to scale inputs or outputs?

texts:
I am writing some simple matlab code to refresh my convolution knowledge (it has been awhile) in regards to DSP.
I have an input signal x and a filter -
x = [1 3 6 2 8];
h = [2 9 5 3];
y=conv(x,h);

which results in
y = [2 15 44 76 73 100 46 24];

I understand these values are large because convolution has lots of multiplication and addition, however
if this was an audio signal I was trying to filter, these values seem very large.
Do I need to scale the input filter or output signal some how?  If so which one and by how much?


texts:
Consider the discrete-time time-invariant system with input $x[n]$ and
output $y[n]$ satisfying
$$y[n] = \sum_{k=1}^5{x[n-k]}$$
Consider approximating the desired system with a second-order IIR
system with system function
$$H^{'}(z)= \frac{1}{1+a_1z^{-1}+a_2z^{-2}}$$
Use the following error criterion:
where $h_d$ is the desired impulse response.
$$E = \sum_{n=-\infty}^{\infty}\left\lvert h_d[n]+a_1h_d[n-1]+a_2h_d[n-2]\right\rvert^2$$

How can the particular error function be useful in solving the system of equations or the desired impulse response?
Do I need to take derivative w.r.t both $a_1$ and $a_2$ and make it 0 to get the system of equations?
Now this is what I am getting:
$$0 = \sum_{n=-\infty}^{\infty}( h_d[n]h_d[n-1]+a_1h_d[n-1]h_d[n-1]+a_2h_d[n-2]h_d[n-1])$$
and
$$0 = \sum_{n=-\infty}^{\infty}( h_d[n]h_d[n-2]+a_1h_d[n-1]h_d[n-2]+a_2h_d[n-2]h_d[n-2])$$
Now how can I solve these system of two equations and how can I approximate the desired impulse response from that?


texts: image reconstruction using a set of kernels

texts:
I try to be brief. For what I understood Convolutional Neural Network CNN for style-transfer extract/learn the main features of source Image. I think these features basically are Kernels (of different sizes depending on Levels). What I ask is how to reconstruct image once we obtained the Kernels. To simplify let's assume I have N-kernels of the same size. How should I reconstruct the Image using only the Kernels ? If question is not clear I'll try to clarify it. (I try using and weighting their responses, but not succeed) Thanks.

PS I wrongly posted to stackoverflow the question... but I think this DSP is the right place. link to Stackoverflow: https://stackoverflow.com/questions/41833641/image-reconstruction-using-kernels


texts:
The transfer-function:
  $$G(s) = \frac{\beta s}{s^2 + \beta s + \omega_0^2}$$

  is to be used in an application that requires the magnitude of the frequency response to be of the band-pass form.

  Show that this property is present by considering the cases when $\omega = 0, \omega = 1$ and show also that $\lvert G(j\omega_0)\rvert = 1$.
  Determine also the equation that defines the frequencies where the filter magnitude is $−3\textrm{ dB}$.

I use $G(s)=G(j\omega)|_{s=j\omega}$, then
$$
\lvert G(j\omega_0)\rvert= \frac{j\beta \omega_0}{-\omega_0^2 + j\beta\omega_0 + \omega_0^2}=1
$$
when $\omega=0$, the numerator equal to $0$ so $\omega=0$ is be filtered out.

My questions are:

Could I use the $G(s)=\lvert G(j\omega)\rvert$ directly or $G(s)=G(j\omega)=\lvert G(j\omega)\rvert e^{-jk\omega}$?
I have no idea about the case when $\omega=+\infty$?
The last question, I use $-20\log\lvert G(j\omega)\rvert=-3\textrm{ dB}$, so we can solve the $\omega$ from this equation, is that correct?


texts: Need help understanding 16QAM constellation diagram

texts:
Is the identification of each point on the constellation diagram arbitrary as long as the sender and receiver know the mapping?

I have seen some diagrams which order in an outward spiral, while others go left to right. Does it matter?


texts:
For my FPGA programming lab final project, I want to implement an audio spectrum analyzer. Since I won't be recreating any signals, I only need the frequency magnitude, and not the phase. I'm hoping that ignoring the phase will spare me the trouble of using complex arithmetic, and make the algorithm leaner. I found some references to the discrete hartley transform, which has "no intrinsic involvement of complex numbers" according to wikipedia. I also found the "real DFT" which splits the complex numbers into sine and consine parts. I was wondering, though, how the radix-2 fft butterfly method would be modified for the real DFT or the hartley dft. Or, maybe there's another method to get just the magnitude information that's better for hardware.

Thank you for any advice.


texts: How can a RLS algorithm utilise Wiener filter as FIR (M-tap) block?

texts:
I'm currently working with a dataset of $5000$ pulses of $N=15000$ samples each. I managed to implement the RLS algorithms with a FIR M-Tap filter such that $M\leq 15000$ ($150$ seems to achieve the best results).

I would like to initialise the coefficients of my filter to the values achieved via Wiener (optimal filter of $N$ coefficients).


texts: Black-box system identification procedure

texts:
I am fairly new to system identification (some background in DSP), and I am trying to use a linear black-box system identification procedure. I am specifically looking for recommendations of sources with a well-described procedure on how to implement black-box system ID and criteria for a properly identified black-box system. I am currently using the method described in "Nonlinear System Identification: NARMAX Methods in the Time, Frequency, and Spatio-Temporal Domains" by Billings, but I am not having much luck. Basically any good starting place would be helpful!

More possibly relevant info about my project: I have a system that I am reasonably confident is linear. I only have data for a step input and I do not have a way of generating other kinds of responses. Signal to noise ratio is extremely high. The manufacturer of the system has purported a simple model. It is clearly an erroneously low order, but is this enough to use a gray-box system ID? I know even less about that, so if I have it all backwards, just ignore this part of the question.

Thanks for the help!


texts: Given the absolute value of a function, find hilbert pairs that make it up

texts:
This might be a solved problem (or a known unsolved problem) because it is related with finding a causal signal that has a desired magnitude response.

$$F(f)=R(f)+iI(f)$$

$$F(f)=R(f)+i \int^\infty_{-\infty}\frac{R(u)}{f-u}du$$

$$F(f)F^*(f)=R^2(f)+\bigg(\int^\infty_{-\infty}\frac{R(u)}{f-u}du\bigg)^2$$

We only know $F(f)F^*(f)$

Is it possible to solve this, is the solution unique?


texts: So I am interested in finding delay between two noisy narrowband signals,for example, two arterial pulse waveform. I can easily get a correct estimate of delay using cross-correlation method, but it fails in low SNR cases. I know Generalized Cross Correlation (GCC) method is supposed to give better results in these cases, but as I am working with narrowband signals, and also due to incorrect estimates of power spectral densities, GCC is very unstable, and the results are worse than normal cross-correlation. So my question is are there any methods for finding time delay between signals in low SNR and are more accurate than the above 2 mentioned methods?


texts: Plotting impulse functions on the same time axis with other functions in MATLAB

texts:
I am trying to plot the following function:

I have to plot another function along with it all on the same plot. Got the other function working, but I can't seem to get this one working. The code below is my attempt at plotting the function:

t = 0:1:10;

for t=>1 && t<=1

    n = 0:1;

end

plot(t,n)

Need help


texts: In this question, it has been proved that $x(t) = \sin(2\pi f_0 t)$ and $x_k(t) = \sin(2\pi (f_0 + k f_s) t)$ have the same sample points. So sinusoids with frequencies $f_0$ and $f_0 + kf_s$ will identify as the same signals. In other words, we keep sampling frequency constant and increase frequency by $kf_s$. Apparently, only frequencies $f_0 + kf_s$ lead to aliasing. Why this frequency is special? I tried to see this result in the frequency domain but it wasn't useful. I know by sampling in time domain we get shifted replicas of the original spectrum which is two delta functions at $f_0$ and $-f_0$ but how this corresponds to $f_s$ and aliasing at $f_0 + kf_s$? What I'm trying to get here is a visual explanation of mentioned result. In the case of the sampling theorem it's really easy to see aliasing occurs when there is an overlap.


texts: Why entropy is undefined in low contrast image

texts:
I am trying to measure contrast of image by entropy of histogram of image.

Code for computing entropy

float measureContrast_inImage(Mat imagel) {

    Mat hist;
    /// Establish the number of bins
    int histSize = 256;
    /// Set the ranges ( for B,G,R) )
    float range[] = { 0, 256 };
    const float* histRange = { range };
    bool uniform = true; bool accumulate = false;
    /// Compute the histograms:
    calcHist(&imagel, 1, 0, Mat(), hist, 1, &histSize, &histRange, uniform, accumulate);
    hist /= imagel.total();

    Mat logP;
    cv::log(hist, logP);

    cv::imshow("logP", logP);

    //Performs an element-wise multiplication or division of the two matrices
    float entropy = -1 * sum(hist.mul(logP)).val[0];

    cout << entropy << endl;

    //float entropy = 0;

    waitKey(0);
    return entropy;
}

Lena image

But entropy is -nan(ind)

similarly, entropy is undefined for another Image.

Entropy is 5.29657 for another image.

Is there any problem in my code.

I wonder that why in low contrast image, entropy is undefined. Kindly help me to figure out.

Edited

I computed histogram and get entropy by the following code
and my entropy is 256

float entropy(Mat seq, Size size, int index)
{

    int cnt = 0;

    float entr = 0;

    float total_size = size.height * size.width; //total size of all symbols in an image

    for (int i = 0; i<index; i++)
    {

        float sym_occur = (float)seq.at<uchar>(0, i); //the number of times a sybmol has occured

        if (sym_occur>0) //log of zero goes to infinity

        {

            cnt++;

            entr += (sym_occur / total_size)*(log2(total_size / sym_occur));

        }
    }

    cout << "cnt: " << cnt << endl;

    return entr;
}

void imhist(Mat image, int histogram[])
{
    // initialize all intensity values to 0
    for (int i = 0; i < 256; i++)
    {
        histogram[i] = 0;
    }

    // calculate the no of pixels for each intensity values
    for (int y = 0; y < image.rows; y++)
        for (int x = 0; x < image.cols; x++)
            histogram[(int)image.at<uchar>(y, x)]++;
}

Mat histDisplay(int histogram[], const char* name)
{
    int hist[256];
    for (int i = 0; i < 256; i++)
    {
        hist[i] = histogram[i];
    }
    // draw the histograms
    int hist_w = 512; int hist_h = 400;
    int bin_w = cvRound((double)hist_w / 256);

    Mat histImage(hist_h, hist_w, CV_8UC1, Scalar(255, 255, 255));

    // find the maximum intensity element from histogram
    int max = hist[0];
    for (int i = 1; i < 256; i++) {
        if (max < hist[i]) {
            max = hist[i];
        }
    }

    // normalize the histogram between 0 and histImage.rows

    for (int i = 0; i < 256; i++) {
        hist[i] = ((double)hist[i] / max)*histImage.rows;
    }

    // draw the intensity line for histogram
    for (int i = 0; i < 256; i++)
    {
        /*line(histImage, cvPoint(bin_w*(i), hist_h),
        cvPoint(bin_w*(i), hist_h - hist[i]),
        Scalar(0, 0, 0), 1, 8, 0);*/
        line(histImage, cvPoint(bin_w*(i), hist_h),
            cvPoint(bin_w*(i), hist_h - hist[i]),
            Scalar(0, 0, 0));
    }

    // display histogram
    //namedWindow(name, CV_WINDOW_AUTOSIZE);
    //imshow(name, histImage);
    return histImage;
}

int main () {
    Mat image = imread("lena.png", CV_LOAD_IMAGE_GRAYSCALE);
    int hist[256];
    imhist(image, hist);
    Mat histImg = histDisplay(hist, "hist");

    float en = entropy(histImg, image.size(), 256);
}

I believe that it is still wrong.
Please correct me .


texts: Hilbert Transform on split datastream produces phase flip

texts:
this is my first post here!

I'm trying to do phase correction on a signal from a linear CCD (2048 pixels, 10bit per pixel, 250kHz linerate) in an FPGA. Therfore I need to generate an analytical signal from the real sensor data. The sensor outputs 8 pixels in parallel to achieve the high datarate of 512 megapixels per second.

The signal is first k-linearized and normalized and then fed into a real FFT, the imaginary part is set to zero. As a tradeoff between speed and resources I buffer the datastream and perform the FFT on 4 parallel channels at double speed. The first channel gets the 0th,4th,8th pixel, the second channel the 1st,5th,9th pixel, etc.

Then the upper half of each spectrum is then zeroed and an inverse FFT is performed. The 4 resulting analytic signal channels are then combined according to the algorithm from this page (https://www.dsprelated.com/showarticle/63.php).

This all works fine if the main frequency of the original signal lies between 0 and 256 or between 513 and 768. But when the frequency is 257-512 or 769-1024 the phase flips by 180° (the imaginary part has the opposite sign). These numbers correspond to the size of the four segments (256 pixels per segment.)

This doesn't happen when the same calculation is done with only one channel.
But due to speed and resource constraints, the FFTs have to be parallelized.

Is there any way to prevent this?

Unfortunately, I don't have much experience with the mathematical side of things, I just convert C algorithms into parallelized, pipelined FPGA code.
But I have built the whole processing in Matlab and verified the problem is not just a bug in the VHDL implementation.

I'd be grateful for any input you might have! Thanks!


texts: Is there a way to get output in second order sections (SOS) format from Scipy's signal.firwin

texts:
I am using Scipy's signal.firwin to implement a lowpass filter. Using the rule of thumb formula(s) to estimate the order of the filter, I am getting a very high number. I looked into this question, and as it turns out I need to break up the filter into second-order sections.

The linked question mentions how to get output from Butterworth filter in 'sos' format. But I couldn't find how to get output in second-order sections (SOS) format from Scipy's signal.firwin.

Is there an easy way to do this?


texts:
What are the differences between random binary sensing matrix  and
random Gaussian sensing matrix?
What the advantages and disadvantages of each matrix?

How can I choose the suitable matrix for a certain signal?


texts: Confusion about subtle difference between discrete-time and continuous-time

texts:
In Alan Oppenheim's book Signals and Systems a comparison is made between the properties of discrete-time and continuous-time complex exponential signals in section 1.3 pg. 26.  Specifically it says:

  The continuous-time complex exponential $e^{j\omega_0t}$ has two
  properties:
      1.  the larger the magnitude of $\omega_0$, the higher is the rate of oscillation in the signal
      2. $e^{j\omega_0t}$ is periodic for any value of $\omega_0$

  Consider the discrete-time signal $e^{j(\omega_o+2\pi)n}= e^{j2\pi n}
e^{j\omega_0 n}=e^{j\omega_0 n}$. From this we see that the
  exponential at the frequency $\omega_0+2\pi$ is the same as that at
  frequency $\omega_0$.Thus we have a very different situation from the
  continuous-time case, in which signals $e^{j\omega_ot}$ are all distinct for
  distinct values of $\omega_0$

The part I don't understand is the one in bold . It seems to me I can do exactly the same thing for continuous-time complex exponential by replacing $n$ with $t$:

$e^{j(\omega_o+2\pi)t}= e^{j2\pi t}e^{j\omega_0 t}=e^{j\omega_0 t}$

So there we have it.  I have produced the same signal in the continuous-time complex exponential.  So what are they talking about? I see no difference.


texts: How can i define the vector to variable scalar?

texts:
I have A,B,C,D,E values all in size 1x1155.
A=[1.0001, 1.0005, 1.00010,...,1.0009]. similar to B, C, D and E.
I need (A(1)+B(1))/(C(1)+D(1)E(1)) + (A(2)+B(2))/(C(2)+D(2)E(2)) + ... + (A(1155)+B(1155))/(C(1155)+D(1155)E(1155)). Therefore i need each individual variable using loop rather than defining one by one from 1 to 1155. Is there a way to do that in Matlab?


texts:
I'm trying implement an IIR filter. I know the IIR filter come from that equation:

and i found that for high pass filter the coefficients are:

where fc = cut freq.

here my matlab code:

fs = 1000;
fc=0.3;
r=exp(-2*pi*fc);
a = [(1+r)/2 -(1+r)/2];
b= [1 -r];
n = 0:1/(fs):20;
n1 = 1024;
x1 = cos(0.3*2*pi.*n);
x2 = cos(0.8*2*pi.*n);
%plot(n,x1);
%hold on
%plot(n,x2);
%hold on
%plot(n,x1+x2);
%aux = fft((x1+x2),n1);
%aux = aux(1:n1/2);
%f = (0:n1/2-1)*(fs/n1);
%length(f)
%legend('x1[n]','x2[n]','x1+x2')
%plot(f,abs(aux));
%hold on

y = filter([0.5759 -0.5759],[1 -0.1518] , x1+x2);
plot(n,y,'r-',n,x2,'b--',n,x1,'g-.');
grid
%axis([0 100 -1.2 4]);
ylabel('Amplitude');
xlabel('Time index n');
legend('y[n]','x2[n]','x1[n]')

the commented parts are from FFT part of program and it's working well but the filter function isn't working. Help me plz and sorry about my english.


texts: Determine the system function H(s) of a system and find out the differential equation

texts:
I have created the following system for practice purposes. From this system I want to determine the system function H(s).

In the picture I have worked with auxiliary (dummy) variables, which should simplify the whole calculation (see blue marking).

How did I approach the matter now? Well, first I have drawn in my auxiliary variables, then I have looked at the summing elements.

I know the relation $H(s)=\frac{Y(s)}{X(s)}$, so I have to transfer my equations into this form.

So these are my equations:

$$\text{Eq 1:}\quad X(s) + G(s) = P(s)$$
$$\text{Eq 2:}\quad 3P(s)=Y(s),$$
$$\text{Eq 3:}\quad G(s) = s^{-1}Y(s)+2Y(s)$$
Eq 3 in 1 results in:
$$X(s) + s^{-1}Y(s)+2Y(s) = P(s)$$
Using Eq. 2 results in:
$$Y(s) = 3X(s) + 3s^{-1}Y(s) + 6Y(s)$$
Therefor:
$$Y(s)(-5-3s^{-1})= 3X(s)$$
$$H(s)=\frac{Y(s)}{X(s)}=\frac{3}{-5-3s^{-1}}$$

My first question is, if you see the system function $H(s)$, is it correct so far?

Then I would be interested in the conversion into a differential equation, here I have already prepared something:

$$Y(s)(-5-3s^{-1})= 3X(s)$$

This should now be handled with the inverse Laplace transformation:

$$-5y(t)+\mathcal{L}^{-1}(-3s^{-1})= 3x(t) \\-5y(t)-3\int_{0}^{\tau} x(\tau) d\tau = 3x(t)$$

My problem at this point is that when I look into my Laplace transform table, I find that $s^{-1}Y(s)$ transformed back results in an integral. Now it is said that these systems can also be given in a differential equation (in terms ot the z-trafo this would be comparable to the difference equation there). How would one come to this differential equation? I have a hunch but I would be interested in your answers :)

I hope my question is so far understandable. Thank you very much.


texts: What is the difference between Antenna port and Physical Antenna in LTE?

texts:
I have just started working on the 5G NR physical layer and I have difficulty understanding the concept of Antenna Port at the downlink.

In release 10 there supposedly 15-22 antenna ports for CSIRS (Channel State Information Reference Signal) and different transmission modes have different antenna ports assigned to them how are these antenna ports actually mapped to physical antennas.

As I understand each port has a resource grid for itself how does this antenna ports resource grid map onto a physical antenna?

Is it possible to have two antenna ports but only one physical antenna or one antenna port but two physical antenna's? If yes then how does the Resource elements map onto the physical antenna's?


texts: Strange results from my FM synthesis algorithm

texts:
It seems that I haven't studied enough about FM Synthesis.
My current algorithm looks like this:

modulator = sin( 2 * pi * fm_freq * i/sample_rate )
carrier = sin( 2 * pi * ( freq + modulator ) * i/sample_rate )

Instead of generating one stable waveform, it generates a sound warping from a pure sine wave to something awfully distorted. What am I doing wrong?


texts:
My question has to do with integrating gaussian noise.

Let us assume we have samples of discrete gaussian white noise  with mean $\mu = 0$ and variance $\sigma_{th}^2$. These noise samples are passed through the system shown in the Figure (a cascade of two integrators with outputs $y_1[n]$ and $y_2[n]$, respectively).

What will be the mean and variance of $y_1$ and $y_2$ (let us say after $N$ cycles)?


texts: Sharp EL-5500III (PC-1403) cassette tape encoding

texts:
I have an old audio cassette tape of digital data I converted to .wav and am hacking to recover the original data. I demodulated the 4K/2K FSK modulation with minimodem. Then I get stuck with insufficient understanding of the NRZI bit-stuffing algorithm. Short-cut: does anybody know of sox-like or minimodem-like program that does NRZI decoding through Linux pipes?

Many references explain that bits are stuffed into NRZI data streams after 5 consecutive bits of one polarity. But to be precise, I haven't found a clarification. Does this mean:

a bit stuff happens after 5 bits, or
a bit stuff happens in order to prevent a run of 6 bits?

By example, if the original data is ...11111001... does it get transmitted as is, or does a bit stuff occur to send the data ...111110001... ? If the bit-stuff is only to cause a transition after 5 bits, then the stuff is not necessary. But the docs say "after 5 bits" which would create unnecessary stuffs if the original data had inherent runs of exactly 5 bits.

Here's the bit stream I'm trying to decode.

root@kali:~# minimodem --rx -R 44100 --mark 4000 --space 2000 --binary-raw 26 -f monoshort8008.wav 500
CARRIER 500 @ 4000.0 Hz
11111111111111111111111111
11111111111111111111111111
11111111111111111110111010
00001111100000100000111110
00001000001111100000100000
11111000001000001111100000
10000011111000001000001111
10000010000011111010101011
11111110100010001011111000
00100000110000010010111000
00101100110101110111111000
10101000110000010101111000
00100000110100010001011000
00100101110101110101011000
10101001110110010101111000
10101000110101110000011011
00100000110101110100011001
00101011110110010100011000
00101011110000010000011010
00100111110000010110011010
11101011110001010100111000
00101011110000010000011001
00100001110000010110011010
11101001110001010100111000
00101011110000010000011011
00100100110000010110011000
11100010110001010100011000
00101011110111110111111011
11101111110000110111011110
NOCARRIER ndata=29 confidence=150.653 ampl=0.126 bps=500.04 (0.0% fast)

Notice the sets of 5 bits near the beginning. Notice also the periodicity of 13 bytes, visible as vertical column runs of 1 or 0 through all the data.

Andy suggested Octave check the frequencies of the modulation.  Here is a plot of the data - x axis is an unscaled count of data points. Data was collected at 44,100 Hz.  Working the scale factor shows 2KHz and 4KHz, 2msec bit time. Matches Audacity waveform view and Minimodem FSK demodulation.


texts: Wave-time trace to Amplitude and PSD spectrum

texts:
For a project I’m currently working on, I composed a irregular ocean wave (function of time) by the summation of n regular waves. For the validation and further processing I’d like to produce the following graphs:

Power Spectral Density (PSD) spectrum of the wave (single sided) with:
Units: x-axis: [Hz], y-axis: [m^2/Hz]
Units: x-axis: [rad/s], y-axis: [m^2/(rad/sec)]

Question 2 In order to produce the single sided PSD graphs I found some code online. The peaks are shown at the right period or frequency however, the information on the y-axis has to be scaled some way (changing the amount of samples per second (i.e. parameter "Fs") or lenght of the signal (i.e. parameter "StopTime") results in other values on the y-axis). Is this value shown in "bins"? How to change this correctly when the periodogram function is used?

CODE 2

Fs          = 6;                     % samples per second
dt          = 1/Fs;                  % seconds per sample
StopTime    = 30*3600;               % seconds
t           = (0:dt:StopTime-dt)';   % seconds
N           = size(t,1);

W = 5*cos(t) + 10*cos(8*t);
p = 5;  %number of periodograms
        %p=10;
        %p=5;
        %W is divided into 20 samples/periodograms 18000/20=900 samples/periodogram
        %W is divided into 10 periodograms 18000/10=1800 samples/peridogram
        %W is divided into 5 periodograms 18000/5=3600 samples/periodogram
        %Very long records needed!!! 90 min for 10 periodograms or 45 for 5.

%Periodogram samples
%x=[];
for i=1:p
    x(:,i)=W((1+(i-1)*N/p):((N/p)*i));
end

%Next power of 2 greater or equal to length of the sample vector to
%calculate fft. FFT needs a length power of 2 to be efficient. For a
%shorter record, zeros will be added to reach this value.
nfft=2^(nextpow2(length(x(:,1))));

%Calculate PSD
Pxx = [];
for i=1:p
    [Pxx(:,i),f]=periodogram(x(:,i),[],nfft,Fs);
end

mW=sum(Pxx,2)/p;    %Average periodograms

ft=1./f; %Wave period vector for the plot. f represents the frequency vector

%Wave frequency spectral density plot
figure;
plot(f,mW,'-')
xlabel('Frequency (Hz)')
title('Power wave Spectrum')
grid;

%Wave Spectral density plot with wave period
figure;
plot(ft,mW);
xlabel('Wave period (s)')
xlim([0 25])
title('Power wave Spectrum')
grid;


texts: Does the system function H(z) for a filter always have symmetry above and below the real x-axis?

texts:
When I look at the Pole/Zero diagrams in the DSP books it appears that whatever
pole or zero is above the real x-axis from:

$$\omega=0 \rightarrow \pi$$

is always reflected as a mirror image below the x-axis from:

$$\omega=0 \rightarrow -\pi$$.

(And whatever pole or zero is on the real axis with $Im\{H(z)\}=0$ is not reflected anywhere.)

For instance, if you have a zero at $e^{j\omega_0}$ and ($\omega_0\ \ne \pi$, and  $\omega_0\ \ne 0$) then you are always going to get a zero that is below the x-axis at $e^{-j\omega_0}$ regardless what type of filter or function the system function H(z) is implementing.

I'm wondering what law or rule this relationship comes from?  and are there any cases where this is not true for a realizable system function of H(z)?


texts:
I am trying to reconstruct the Zone Plate image and am struggling to remove the last remaining aliasing.

ShaderToy: https://www.shadertoy.com/view/wdGGWK

In the shader above you can see that I am sampling each pixel independently 1000 times with [-3, 3] jitter, and filtering those samples with (currently) Lanczos2. Why is there still aliasing?

My understanding is that the aliasing just at the edge of the visible pattern is from the filter and can't be removed without a different choice of filter, but the subtle aliasing in the middle of the axis seems like an error, and one that I have been unable to remove.

Any ideas?


texts: Why can't realisable IIR filters have linear phase?

texts: I am studying IIR Filter Design and came across this arbitrary statement in my textbook which says that 'physically realisable and stable IIR filters can not have linear phase'. Would really appreciate if anyone can elaborate on this.


texts: Finding Duplicate Pixel/Objects along Image seams

texts:
I have a seam of two photos joined by a seam alignment algorithm that on occasion generates duplicated pixels/objects/artifacts on both sides of the seam in the final, aligned photo. The images are large and currently the seams are checked manually for duplicate artifacts/objects. Is there an algorithm/software/programming package that will allow me to more quickly check along the seams for duplicates, preferably automatically? The algorithm should be able to spot arbitrary duplicates on both sides of the seam. What I am dealing with is sketched below as:

 |
 |
 | <---- Seam
 |
o|o <---- Duplicated object
 |
x|x <---- Duplicated object
 |
 |


texts:
Bandwidth is the difference between the upper and lower frequencies in
  a continuous band of frequencies. A key characteristic of bandwidth is
  that any band of a given width can carry the same amount of
  information, regardless of where that band is located in the frequency
  spectrum.

So I understood bandwidth to be a measure of how much space it takes up in the frequency domain. Am I right in saying that because of symmetry about the Y axis, the same information is repeated for negative frequency values, this is why the "bandwidth" of the bandlimited signal given below is not $2B = B - (-B)$ but $B$.

1) Since the left is a mirror image of the right part only for real signals, if the signal was not real, can we still say bandwidth is $B$?

2) Why do we have two different definitions for passband and baseband bandwidth? Is there any other reason?

  Passband bandwidth is the difference between the upper and lower cutoff
  frequencies of, for example, a band-pass filter, a communication
  channel, or a signal spectrum. Baseband bandwidth applies to a
  low-pass filter or baseband signal; the bandwidth is equal to its
  upper cutoff frequency.

3) Am I right in understanding that in case of either type of signals, alias-free signal sampling is possible at $2B$ samples per second?


texts: Calculate impulse response when output contributes in input

texts:
I have an exercise in which I need to find the impulse response for this given system:
$$y(n)=\frac{1}{2}y(n−1)+x(n−1)+x(n)$$
As per my knowledge, I need to find the homogenous solution. My homogenous solution gives me $$y(n)=\frac{1}{2^n} u(n)$$ But the actual answer shows $$y(n)=\frac{1}{2^n} u(n) + \frac{1}{2^{n-1}}u(n-1)$$ I don’t get where the second part came from.
How to solve this?
Thanks!


texts: How to plot noise shaped spectrum of First order Incremental Sigma Delta ADC's output?

texts:
A little background:

This question is about viewing the output spectrum of incremental delta sigma modulators (IDSM). IDSM is a subcategory of general delta sigma modulators where the integrators are reset for every conversion cycle with the advantage of being multiplexable along several channels. (https://pdfs.semanticscholar.org/1f62/af5ec0ba627e8ef5ebf81a479307b059a116.pdf).

I am aware that a general sigma-delta noise shaping architecture has a noise shaping characteristic of a differentiator as shown below:

1,2 and 3 indicate the orders of noise shaping and sigma delta architectures too.

For example, if we consider the following continuous-time sigma-delta architecture (first order), the noise shaping is governed by the equation:

NTF(z)=(1-z^(-1)), which translates to the curve (1) of the first figure.

Simulation background for IDSM:

When I simulate an incremental sigma delta converter (IDSM), I run it for (N*(1/fs)*OSR), where it runs for OSR clock cyles for single conversion. I run the output through the counter and then scale it by (1/M) for every reset (Decimation) and then view the spectrum with N points, I am able to get the desired output SNR for an IDSM.

One conversion cycle is shown below (M=OSR).

Core question:

I designed a continuous-time first-order sigma delta modulator in matlab and the output spectrum is shown below:

Here I am able to see the first order noise shaping (-20dB/decade noise increment).

Now, how do I simulate the output from an IDSM to view the result of noise shaping?.

I saw some articles, plotting noise shaped spectrum for IDSM (in free running DSM mode). I am really not sure how to run an IDSM in this mode, because logically it would require a reset for each conversion and therefore I am stuck.

(Again, I could decimate the output of IDSM and then view the spectrum - which gives me the correct SNR. But I'm interested in viewing the noise shaping behavior.)

If what i am expecting is not achievable, please let me know as well - that would be helpful too.


texts: What's the meaning of a complex zero/pole?

texts:
I have been studying signal processing and control for a while now, and I use Laplace and Fourier transforms almost everyday. Also another tools such as Nyquist or Bode plots.

However, I had never thought of this until today: what is the physical meaning of a complex number when dealing with frequencies?

This may sound silly, but I was asked this question and I didn't know what to answer. Why do we talk about $j\omega$ and not just $\omega$ in, for example, Fourier transforms and Bode or Nyquist plots? What is the physical sense of the real and imaginary part of a zero or a pole in the Laplace domain?


texts: Why will the FDMA usually be used in transmitting analog signal,not in digital signal?

texts:
I saw this in the WiKi:https://en.wikipedia.org/wiki/Frequency-division_multiple_access

FDMA can be used with both analog and digital signal but it generally used with analog signal

Why?Don't we usually transmit the digital signal now?why don't we usually use FDMA to transmit digital signal but analog?


texts: MATLAB, How to filter a discrete signal?

texts:
I have a 1-by-10000 matrix of double`s stored in a file named "fecg.mat". The matrix represents the magnitude of a recorded FECG signal.

I've plotted it against time (from 0 to 9999):

For removing the baseline I wonder if I can use a high-pass filter. How do I design a proper filter?

P.S
Signal processing isn't my field of study. I have no idea how to filter a discrete time-domain signal.


texts: Low-pass filters not producing smooth output

texts:
I'm relatively new to this material, so apologies in advance if I'm missing something obvious.
I'm trying to smooth some very noisy signals. What I observe is that if I apply an ewma, even with a large center of mass, the resulting signal is still relatively noisy (in the sense of jittering from step to step). If however, I apply a filter with relatively low weights at the front end then the result is much smoother.
This effect is illustrated in the examples below.
At the level of impulse response this sort of makes sense to me, as there is large difference in filter weight between lag 0 and lag 1. Is it possible to understand this effect at the level of frequency response? I would have expected to see the ewma crossover killing high frequency noise much more than the ewma filters. However, this does not appear to be the case.
Replicating code
import matplotlib.pyplot as plt
import pandas as pd
import numpy as np
from scipy import signal

def apply_filter(X, f):
    result = []
    for i in range(len(X)):
        result += [(X.iloc[:i + 1] * f[:i + 1][::-1]).sum() ]

    return pd.Series(result, index=X.index)

X = pd.Series(np.random.normal(size=1200))  # .diff()

indicator = pd.Series([0] * 1000 + [1] + [0] * len(X))
filter = {'raw': indicator,
          'ewm20': indicator.ewm(20).mean(),
          'ewm100': indicator.ewm(100).mean(),
          'mean5ewm20': indicator.rolling(5).mean().ewm(20).mean(),
          '5-15crossover': indicator.cumsum().ewm(5).mean() - indicator.cumsum().ewm(15).mean(),
          'firwin': pd.Series([0] * 1000 + signal.firwin(80, 0.1, window=('kaiser', 8)).tolist() + [0] * len(X))}
filter = {k: (v / v.abs().sum()).iloc[1000:1001 + len(X)].values for k, v in filter.items()}
smoothed = {k: apply_filter(X, v) for k, v in filter.items()}
freq = {k: pd.Series(signal.freqz(v)[1], index=signal.freqz(v)[0]) for k, v in filter.items()}

fig, axs = plt.subplots(len(filter), 3, figsize=(16, 8))
pd.DataFrame(smoothed).loc[:, filter.keys()].plot(subplots=True, ax=axs[:, 0])
pd.DataFrame(filter).loc[:, filter.keys()].iloc[:400].plot(subplots=True, ax=axs[:, 1])
np.log10(pd.DataFrame(freq).abs().loc[:, filter.keys()]).plot(subplots=True, ax=axs[:, 2])
axs[0, 0].set_title('Data')
axs[0, 1].set_title('Impulse Response')
axs[0, 2].set_title('Frequency Response')
axs[-1, 0].set_xlabel('Time')
axs[-1, 1].set_xlabel('Lag')
axs[-1, 2].set_xlabel('Frequency [rad/sample]')
plt.tight_layout()


texts: How does averaging increase the Signal to Noise ratio?

texts:
When attempting to spot activities on an EEG waveform given a stimuli, it is recommended that several trials be taken and their average is done to increase the SNR.

Intuitively, since noise is not stationary or localized, averaging different trials will have a smearing effect on the noise.

But how exactly does it increase the SNR since the noise is still present (but averaged out)? Can we more rigorously define the effect of averaging on the SNR?


texts: Event detection in a running person's acceleration data

texts:
I am trying to detect some events on this signal :

This is the acceleration of someone running. By eyes we can detect many blocs, sometimes separated by a little interval. I would like to be able to detect automatically all those blocs.  With this post, it helped me detect those blocks if the signals is less noisy. But if the signal is really noisy and intervals between blocks are really short I can not use the idea of that post (as known force value above a threshold to be 1 and take parts when we find 100 consecutive values of 1 for exemple).

Any documentation or info could be useful... :)


texts: Total bandwidth in a network calculation

texts:
I have data on the number of devices in a wired network and the data rates on links among these devices. I also know the positions of these devices from which I can compute the distances between the devices.

The idea is that if I were replace the wired network with wireless infrastructure, how much bandwidth am I looking for. Would let's say 200MHz be sufficient for such a case.

Assuming that the BER is very good, how should I go about in calculating the total bandwidth in such a network.


texts: Where can I find "recipes" for synthesizing natural sounds?

texts:
For fun, I recently wrote some code in C to generate waveforms by hand and output them. I can play pure sines, sawtooth waves and box waves with simple ADSR envelopes. It's nice.
I'd like to get into synthesizing more "natural" instruments. Classic drum machine sounds like snares and kick drums, piano or vibraphone sounds, that kind of thing. I've found several books with general information on audio synthesis techniques: frequency modulation, cutoffs, general techniques like that. But what I'd be really interested is a kind of "cookbook" for generating sounds: specific recipes on how to create waveforms that sound like a piano, a snare, etc. Some actual case studies in how to use these techniques to do useful work.
Is this type of information available anywhere in a systematic format?


texts: Okay so i understand the first det, which is why i got the correct solution. However for the second det, im not sure how to solve it since because of the different coefficents (0,2,4) and i tried of different solutions. If anyone can give me a clear way to solve these it would be appreciated. Thank You.


texts: How Can I Resample a Signal with an Arbitrary Factor (For Example - 128000Hz to 16000.1Hz) in MATLAB?

texts:
I need to simulate the sampling of a continuous (fsCtu = 128000Hz), acoustic signal with two microphones that have a slight offset in sampling rate (fsMic1 = 16000, fsMic2 = 16000.1) in Matlab. What is the best way to do this?

Things I tried:

The Matlab "resample" command only works for resampling to 16000, not to 16000.1
"interp1" doesn't seem to be an option because I think I need to use bandlimited interpolation for a correct simulation. (Is this assumption correct?)
I tried to write my signal to a wav file and resample it via a system call using this software, then load the processed file with wavread. I'm not sure if this is a good solution. A quick test revealed that this method doesn't give the same result as the "resample" command for resampling to 16000Hz, which I find strange.

Any ideas or suggestions?


texts: even symmetry of magnitude and odd symmetry of phase

texts:
I'll appreciate it if any of you guys could help me with this question:
Suppose that x(t) is a real signal, prove that the magnitude of its Fourier transform has even symmetry and its phase has odd symmetry.
Thanks.


texts: Computing of speech signals modulation entropy

texts:
I got a question regarding to modulation entropy of audio signal. There is some literature I found about this topic, http://www.irit.fr/recherches/SAMOVA/FeaturesExtraction.html#me.

I've read sample .wav file. How do I calculate entropy of this sampled signal of lenght N samples? How should I understand, in this case (in paper equation from above), propability of event?


texts:
I am trying to understand the workings of an MRI scanner and wanted to figure out how to simulate it by solving Bloch equations in the rotating frame of reference.

I have a small query regarding the apparent Larmor frequency in the rotating frame of reference. So, if $\omega_0$ is the Larmore frequency in the laboratory frame of reference and in the rotating frame the coordinate system is rotating about the $z$ axes with a frequency $-\omega_r$. In this case, now the apparent Larmor frequency would be: $w_{eff} = \omega_0 - \omega_r$.

Now, I have an RF pulse which is entered on the Larmor frequency and has some bandwidth which means that it will excite spins with resonant frequencies +-100 Hz around the Larmor frequency and say I have a slice selection gradient which achieves this frequency difference.

Now, it is clear that for the spins which are in the centre, they will have an apparent resonant frequency of 0 Hz and the other spins would have an apparent resonant frequency  +- 100 Hz depending on their position along the z-axes. So, now the time dependent effective field experienced by a spin due to the RF pulse, would this be something like:

$$
B_x(t) = rf(t) \cos(2 \pi \omega_{eff} t + \phi_{rf}) \\
B_y(t) = rf(t) \sin(2 \pi \omega_{eff} t + \phi_{rf})
$$

Here $rf(t)$ is the amplitude of the rf waveform at time $t$. $w_{eff}$ is the apparent resonant frequency in the rotating frame of reference for a given spin and $\phi_{rf}$ is the phase of the applied RF field. I am not sure if this is correct and wanted to verify.

I am not sure if this should belong to the physics forum but I never get any answers there regarding MRI. But please feel free to migrate if this is not appropriate here.


texts: Using Chebyshev polynomials in wavetables

texts:
How are Chebyshev polynomials used in wavetables for audio synthesis? I dont understand how to store these inside wavetables?

How does the synthesized audio result in matching to the correct index in the wavetable and what part of the Chebyshev polynomials are actually stored in the wavetable?


texts: Statistical Signal Analysis - what to read

texts:
Possible Duplicate:
  Introduction to Statistical Signal Processing

Currently I am doing an (undergrad) class on statistical signal analysis, but I find the reader at hand (no book, buh) very hard to read and understand. The bad thing is, my professor did not include a read list of any kind, and this is not really a topic I have favourite authors in.

Do you have any books that you can recommend me to read, to finally understand this subject? I have looked at the MIT opencourseware link, and there is a list of books there, but any of these seems to cover only a single chapter in my reader. If there is a single book that covers all of them it would be great.

The chapters in my reader are: Modeling of deterministic signals and systems, Modelling of Stochastic Signals, Linear filtering of stochastic processes, Estimation, Optimal Filtering and Estimation of correlation functions and power spectra.


texts: Compressive sensing based sparse vector estimation

texts:
I am newbie in compressive sensing (CS), I read about compressive sensing and its use for sparse vector estimation. As I understood CS can be used either in time or frequency domain.  For me, The part I couldn't understand in my case is how to build the measurement matrix for my vector required to be estimated. (In my case the sparse vector $h$ needs to be estimated)
I have a vector $X$ where some values of it are known (let's say the values located at $1:4:end$), called pilots, that vector was convoluted with a sparse vector $h$, to have:
$y = X*h$
$*$ denotes the convolution operation.
My goal is to estimate the sparse vector $h$ based on the known values (pilots) in $X$.  How to build the measurement matrix in that case?


texts: Given a log-plot of frequency-magnitude-phase how to apply it as an EQ curve to a signal?

texts:
I have several hundred data points that represents an EQ correction curve, where each point contains frequency, magnitude, and phase (-1..+1). The frequencies are distributed exponentially, not linearly. Magnitude is expressed as a multiplier. A log-plot of this data is very smooth in appearance.

What is the best way to apply this correction curve to digital audio samples in real time (as a software plugin or implemented on DSP hardware)?


texts: GNU Radio - Bandwidth and Band Reject Filter

texts:
I am currently trying to generate a psk modulated carrier with gnu radio.

Attach is an image of my flow graph.

As can be seen I am generating a bit stream using the Random Source module and implementing the modulation using the Constellation Modulator module. I then Low Pass Filter it and move forward to the conversion from baseband to passband generating Fc +- Fbb. I would like to utilize a bandpass filter to deal with the lower carrier frequency. It seems though when I implement the filter I still see both of the carrier frequency.

Why is this? I can't wrap my head around it.

EDIT: I figured out my bandpass filter question. All I had to do was set the FIR Type in the filter to "Complex->Complex(Complex Taps)(Decim)" from "Complex->Complex(Real Taps)(Decim)"

Now for the bandwidth question...

This is is fairly general question (I am assuming), but how does one define bandwidth in this context?

Typically BW with PSK is established at BB with the highest frequency component representing the BW.


texts: Check if I has right: Is this the Extended Kalman Filter

texts:
I have learn the Kalman-Buncy filter for the LQG controllers. I know that this is a signal processing forum and not robotics not math forum. But Extended Kalman Filters are daily discussed here.

First of all! Let's say that I have model a nonlinear state space model:
$$ \dot{x} = f(x, u, t) + d(t) $$
$$ y_m = g(x, u, t) + n(t)$$

Where:

$ \dot{x}$ is the derivative of the state vector of the system
$x$ is the state vector of the system
$u$ is the input signal vector
$t$ is time
$d(t)$ is the disturbance vector
$f(x, t)$ is the standard form of nonlinear state space model
$g(x, t)$ is the standard form of the output of a nonlinear state space model
$n(t)$ is the noise vector
$y_m$ is the measured output vector

My goal is to compute the Extended Kalman Filter gain $K$ fo the system so I can build my LQG controller.

Extended Kalman Filter uses jacobian linearizations to make the model $f(x, u, t)$ linear on equilibrium points.

So I build new linear state space matrix $A(t)$ and a output matrix $C(t)$ with respect to time. The linear state space model is:

$$ \dot{x} = Ax + Bx$$
$$ y = Cx + Du$$

But in this case I going to build:
$$\delta \dot{x} = A(t) \delta x + B(t) \delta u$$
$$\delta y = C(t) \delta x + D(t) \delta u $$

Where:

$\delta x = x - x_0$
$\delta y = u - u_0$

The $x_0$ and $u_0$ is the equilibrium point and $u$ is the input vector.

I linearize $f(x, u, t)$ and $g(x, u, t)$  and get:

$$A(t) = \left [ \frac{\partial f({x}, u, t)}{\partial {x}} \right ]_{x = \hat{x}}$$
 $$C(t) = \left [ \frac{\partial g({x}, u, t)}{\partial {x}} \right ]_{x = \hat{x}}$$

Where:

$\hat{x}$ is the estimated state vector form the LQE-observer

And of course I get $B(t)$ matrix:
 $$B(t) = \left [ \frac{\partial f({x}, u, t)}{\partial {u}} \right ]_{x = u_0}$$

But to finding Extended Kalman Filter gain matrix $K$, matrix $B(t)$ is not important.

So now I have my LQE-observer:
$$\dot{\hat{x}} = f(\hat{x} , t) + K[y_m - g(\hat{x}, u, t)]$$

And I find the Extended Kalman Filter gain $K$ from this Riccati equation:

$$\dot{P} = A(t)P + PA^T(t) - PC^T(t) R^{-1}C(t)P + Q$$

Where:

$\dot{P}$ representing the state estimation error covariance matrix derivative
$P$ is the state estimation error covariance matrix
$Q$ and $R$ is the weighing matrices

Solve $P$ from the Riccati equation and put in in here:

$$K = PC^T(t)R^{-1}$$

Now we have our Extended Kalman Filter gain $K$.

Now I just plug $K$ into my LQG controller diagram:

Questions:

So all I need to inclue more are linearizations with the jacobian matrix?
Have I done it right?
Which comes first, the linearizations or the estimations? I need the estimated vector $\hat{x}$ to get $K$ and I need $K$ to get $\hat{x}$.
Should $\dot{P} = 0$ when I solving the Riccati equation?


texts: Gaussian White Noise - Relation Between Distribution and Correlation

texts:
Im a beginner in signal processing so my question may be obvious.

A white noise has the property to have its autocorrelation function that is equal to

$$\mathbb{E}[f(t+\tau)f(t)]=\sigma^2 \delta(\tau).$$

We say that random signals can be Gaussian white noise, but for me if the signal has a Gaussian distribution it is necesseraly a white noise because as soon as we know the distribution then $\mathbb{E}[f(t+\tau)f(t)]$ is automatically determined.

Am I wrong ?


texts:
I used xcorr function in Matlab to compute the autocorrelation of my signal. I obtained a graph, which is shown below. How to interpret this? I read that autocorrelation helps in finding the patterns in the data. And i have went through some of the properties of the auto-correlation.

But i am finding difficulty in obtaining information regarding how my input data is correlated with itself. Do i need to compute any other things in order to get this information? How to understand the autocorrelation graph?


texts: The Distribution of Filtered Gaussian White Noise

texts:
Suppose I have some time series $s(t)$ which contains Gaussian white noise generated by a distribution $N(0,\sigma^2)$

Then I apply a filter to s(t) with a frequency response $H(\omega)$, giving me $s_H(t)$

What distribution does $s_H(t)$ have?

I'm particularly interested in cases where $s_H(t)$ is band limited Gaussian white noise.


texts: Why does this parallel IIR filter have such a weird form?

texts:
Referring to "Parallel 1st-order IIR filters" on this wikipedia page, the structure of a 4-parallel IIR filter is discussed. However, the format in which it is implemented seems very odd.

Why is it necessary to iterate the IIR equation N parallel times? Could you not just parallelize it in a direct form. Refer to my chicken-scratch block diagram below.


texts: Calculating Beam pattern, Directivity of microphone array

texts:
As answered in Relationship of Fourier Transform and Spatial Directivity Pattern the relation between directivity pattern and weights is given.  But I am bit confused, I have Beamformer weights (complex valued, for each frequency bin) for a particular steering vector. The relationship mentioned in the above link needs the values to be calculated for each theta, but the beamformer weights are for a particular steering angle theta.
Do I need to calculate beamformer weights for each of the theta, or I have a wrong understanding of the directivity pattern.
Suggest methods for calculating directivity using beamformer weights.
Thanks


texts: Compute the Second Order Derivative of Digital Image with Finite Differences

texts:
I was looking for how to compute second order derivative of an image and came across the question kernels to Compute Second Order Derivative of Digital Image. In the top voted answer, it gives an example to get $3\times 3$ kernel $I_{xx}$ and $I_{xy}$. The 1st order derivative filter is defined as
$$
\left(
\begin{array}{cc}
 -1 & 1 \\
 -1 & 1 \\
\end{array}
\right)
$$
Then

$$
 I_{xx} = I_x \cdot I_x = \left(
\begin{array}{cc}
 -1 & 1 \\
 -1 & 1 \\
\end{array}
\right)\cdot\left(
\begin{array}{cc}
 -1 & 1 \\
 -1 & 1 \\
\end{array}
\right)=\left(
\begin{array}{ccc}
 1 & -2 & 1 \\
 2 & -4 & 2 \\
 1 & -2 & 1 \\
\end{array}
\right)
$$

and

$$
I_{xy} = I_y \cdot I_x = \left(
\begin{array}{cc}
 -1 & -1 \\
 1 & 1 \\
\end{array}
\right)\cdot\left(
\begin{array}{cc}
 -1 & 1 \\
 -1 & 1 \\
\end{array}
\right)=\left(
\begin{array}{ccc}
 1 & 0 & -1 \\
 0 & 0 & 0 \\
 -1 & 0 & 1 \\
\end{array}
\right)
$$

The computation of $I_{xx}$ makes sense to me, however, for $I_{xy}$, shouldn't it be defined as below?
$$
I_{xy} = I_x \cdot I_y = \left(
\begin{array}{cc}
 -1 & 1 \\
 -1 & 1 \\
\end{array}
\right)\cdot\left(
\begin{array}{cc}
 -1 & -1 \\
 1 & 1 \\
\end{array}
\right)=\left(
\begin{array}{ccc}
 -1 & 0 & 1 \\
 0 & 0 & 0 \\
 1 & 0 & -1 \\
\end{array}
\right)
$$


texts:
Hello I am not a domain expert in signal processing but I need to align two audio signals.

I have seen the following page
https://librosa.github.io/librosa_gallery/auto_examples/plot_music_sync

I have implemented DTW and so many other pages. Most of these are using DTW but the problem is that none of them are showing WHAT TO DO AFTER I PERFORM DTW!. DTW gives you a cost matrix and wrapping path (D, wp) but nowhere it is shown how to use these two to align the actual signals.

I have come to a github repo
https://github.com/chausies/voice_transmogrifier/blob/master/sync_audio.py

But most of the parts of the code seems difficult to understand. Like why he used variable_phase_vocoder or what close_points is doing. Why to use PchipInterpolator ? In fact why interpolation is done at all in audio?

Please help me in the project. I am new to audio processing.

Thanks


texts: Derivation of the LMMSE (Linear Minimum Mean Squared Error) Estimate and the MMSE Under Gaussian Prior

texts:
I am learning estimation theory through Steven M. Kay - Fundamentals of Statistical Signal Processing, Volume 1: Estimation Theory. In the Chapter 12 (Linear Bayesian Estimator), Theorem 12.1 (Bayesian Gauss-Markov Theorem) gives the LMMSE estimation of the signal based on the linear noisy measurement under the Gaussian prior assumption:

If the data are described by the Bayesian linear model form
\begin{equation}
\boldsymbol{x}=\boldsymbol{H\theta}+\boldsymbol{w} \tag{12.25}
\end{equation}
where $\boldsymbol{x}$ is an $N \times 1$ data vector, $\boldsymbol{H}$ is a known $N\times p$ observation matrix, $\boldsymbol{\theta}$ is a $p \times 1$ random vector of parameters whose realization is to be estimated and has mean $E(\boldsymbol{\theta})$ and covariance matrix $\boldsymbol{C}_{\theta\theta}$, and $\boldsymbol{w}$ is an $N \times 1$ random vector with zero mean and covariance matrix $\boldsymbol{C}_w$ and is uncorrelated with $\boldsymbol{\theta}$ (the joint PDF $p(\boldsymbol{w},\boldsymbol{\theta})$ is otherwise arbitrary), then the LMMSE estimator of $\boldsymbol{\theta}$ is
\begin{align}
\hat{\boldsymbol{\theta}} & = E(\boldsymbol{\theta})+\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \tag{12.26} \\
& = E(\boldsymbol{\theta})+(\boldsymbol{C}_{\theta\theta}^{-1}+\boldsymbol{H}^T\boldsymbol{C}_w^{-1}\boldsymbol{H})^{-1}\boldsymbol{H}^T\boldsymbol{C}_w^{-1}(\boldsymbol{x}-\boldsymbol{H}E(\boldsymbol{\theta})) \tag{12.27}
\end{align}
The performance of the estimatior is measured by the error $\boldsymbol{\epsilon}=\boldsymbol{\theta}-\hat{\boldsymbol{\theta}}$ whose mean is zero and whose covariance matrix is
\begin{align}
\boldsymbol{C}_\boldsymbol{\epsilon} &= E_{\boldsymbol{x},\boldsymbol{\theta}}(\boldsymbol{\epsilon}\boldsymbol{\epsilon}^T) \\
& = \boldsymbol{C}_{\theta\theta} - \boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T(\boldsymbol{H}\boldsymbol{C}_{\theta\theta}\boldsymbol{H}^T+\boldsymbol{C}_w)^{-1}\boldsymbol{H}\boldsymbol{C}_{\theta\theta} \tag{12.28} \\
& = (\boldsymbol{C}_{\theta\theta}^{-1}+\boldsymbol{H}^T\boldsymbol{C}_w^{-1}\boldsymbol{H})^{-1} \tag{12.29}
\end{align}

Since the prior of $\boldsymbol{\theta}$ is Gaussian, the LMMSE estimate $\hat{\boldsymbol{\theta}}_{LMMSE}$ is equivalent to the MMSE estimate $\hat{\boldsymbol{\theta}}_{MMSE}$, and $\hat{\boldsymbol{\theta}}_{MMSE}$ is equal to the posterior mearn $E(\boldsymbol{\theta}|\boldsymbol{x})$. Since the prior and likelihood are both Gaussian, the posterior distribution $p(\boldsymbol{\theta}|\boldsymbol{x})$ is also Gaussian.
Here I am trying to derive $\hat{\boldsymbol{\theta}}_{MMSE}$ and $\boldsymbol{C}_\boldsymbol{\epsilon}$ from the perspective of PDF multiplication, that is, calculate $p(\boldsymbol{\theta}|\boldsymbol{x}) \propto p(\boldsymbol{x}|\boldsymbol{\theta})p(\boldsymbol{\theta})=\mathcal{N}(\boldsymbol{x};\boldsymbol{H\theta},\boldsymbol{C}_{w})\mathcal{N}(\boldsymbol{\theta};E(\boldsymbol{\theta}),\boldsymbol{C}_{\theta\theta})$, and formulate the quadratic and firse-order terms of $\boldsymbol{\theta}$ at the exponential to form a Gaussian PDF. The covariance matrix of $p(\boldsymbol{\theta}|\boldsymbol{x})$ I got matches 12.29, but the posterior mean is the following form:
\begin{equation}
E(\boldsymbol{\theta}|\boldsymbol{x}) = \boldsymbol{C}_{\boldsymbol{\epsilon}}(\boldsymbol{H}^T\boldsymbol{C}_w^{-1} \boldsymbol{x}+\boldsymbol{C}_{\theta\theta}^{-1}E(\boldsymbol{\theta})) \tag{q1}
\end{equation}
So my question is, is the posterior mean I got in q1 equal to the $\hat{\boldsymbol{\theta}}$ given in 12.26 and 12.27? If so, how can I reach that?
By the way, I can't find the way from 12.26 to 12.27 (12.28 to 12.29 either). So can someone give me a hint?


texts: Hilbert transform of a function $g(t)$ which is defined in time domain, would result in another function in time domain. Is there any other transformation like Hilbert that the results be in time domain?


texts:
Why is binary polar signalling not feasible for a data rate of 2.0Mbps within a bandwidth not exceeding 1MHz?

Is it due to the fact that binary signalling can only transfer 1 bit of data per symbol thus, a 1MHz bandwidth would not be large enough for the given data rate?

Then a 4-PAM signal would carry 2 bits of information per symbol which would be more appropriate?


texts: Bit error rate curves, theoretical and simulation

texts:
The BER is the average number of errors that would occur in a sequence of
$n$ bits.  When $n = 1$, we can
think of the BER as the probability that any given bit will be received in error. Basically, it lies between 0 and 1 --

when BER = 0, we can say no error in the received bits.
If at least one error is seen, we can say with 100
percent certainty that the BER is definitely not zero.
as long as at least one bit was successfully received, we can state that the BER is most definitely less than 1.

When we measure the performance of a digital communication system with and without equalization, with that of the theoretical BER for a particular kind of modulation, we prefer the simulation to be exact match with the theoretical.
My Questions are :

(1) Is theoretical BER the upper or lower performance bound?

(2) Is it desired that the BER tends to zero? Is it reasonable to say that  since this may not happen in real communication, so there is always a quest to come up with an equalizer designed with an objective that the BER tends to zero?

(3) What if the simulation curve with and without equalizer falls below the theoretical? What can be said about the equalization method playing a role and if it is desired that the simulation curve - with and without equalizer is below the theoretical.

(4) What can be said about performance of the system with and without equalizer if the simulated BER is above the theoretical?

(5) Does an equalizer essentially help to improve the BER? By improving -how much improvement is desired? Tradeoff? Do we wish the simulated curve to be an exact match to the theoretical and then say that the equalizer is functioning properly in removing intersymbol interference and noise?

What is the interpretation for the bit error performance?


texts: How can I find expansion coefficients of the a vector in a given basis?

texts:
How can I find the coefficient of the vector $\mathbf y$? And how can the inner product be done on these vectors?

Let $\mathbf y = \begin{bmatrix}1\\2\\0\\1\end{bmatrix}$
What are the expansion coefficients of $\mathbf y$ in the basis $\left\{\mathbf v^{(0)}, \mathbf v^{(1)}, \mathbf v^{(2)}, \mathbf v^{(3)}\right\}$ where
$$\mathbf v^{(0)} = \begin{bmatrix}\frac 12\\\frac 12\\\frac 12\\\frac 12\end{bmatrix}, \mathbf v^{(1)} = \begin{bmatrix}\frac 12\\\frac 12\\-\frac 12\\-\frac 12\end{bmatrix}, \mathbf v^{(2)} = \begin{bmatrix}\frac 12\\-\frac 12\\\frac 12\\-\frac 12\end{bmatrix}, \mathbf v^{(3)} = \begin{bmatrix}\frac 12\\-\frac 12\\-\frac 12\\\frac 12\end{bmatrix}, \quad $$


texts:
In section 1.3 from Signals, Systems and Inference, Alan V. Oppenheim, George C. Verghese. Finite-action signals, which are also called absolutely summable signals, are defined by the condition
\begin{equation}
 \int_{-\infty}^{\infty}|x(t)|dt<\infty
\end{equation}
whereas for discrete time signal, its as
\begin{equation}
 \sum_{k = -\infty}^{\infty}|x[k]|<\infty
\end{equation}
The integration and sum on the left are called the action of the signal. Therefore also known as finite-action signals.
So I would like to know what does finite action signal means? Or is there intuitive explanation about it.


texts: Incoherence: Compressed Sensing (CS) vs Matrix Completion (MC)

texts:
I am seeking a clarification of the concept of Incoherence within the MC framework.

Specifically,

1) the literature mentions the application of a "strong incoherence" given a set of assumptions. Intuitively, how is the strong incoherence different than "incoherence".

2) Does the Restricted Isometry Property have anything to do with low rank Matrix Completion of natural images, or does it only pertain to the reconstruction of signals in CS (i.e. y = Ax).

Background

There is an incoherence property $\mu$ used in CS to measure the correlation between the sensing and representation bases; whereby low correlation enables signal reconstruction of sparse signals with fewer samples and high correlation makes it nearly impossible regardless of the sparsity of the signal.

For MC, there is also an incoherence property $\mu$: it measures the correlation between the singular vectors of the matrix (U and V from SVD) and the standard basis vectors (i.e. $\ e_i$ = vector of zeros except for a 1 at index i). As with CS, low rank (i.e. singular value-sparse) matrices exhibiting a low correlation are more likely to be reconstructed than those with high correlation.

Both of these incoherence properties are used to stipulate the minimum number of samples required to reconstruct a signal (1D, 2D, respectively). However, I cannot see where the matrix incoherence property addresses the nature of the sampling; I.e. it only tells me which matrices are more likely to be reconstructed. I believe this is because the "nature of the sampling" that is assumed (by CS community doing MC) for low rank MC problems is strictly "uniform random sampling".

Now, restricting our attention only to those low rank matrices with low $\mu$, is there an additional MC incoherence property that characterizes the nature of the sampling to say whether or not a matrix can be reconstructed? (I'm looking for a link back to CS where they compare the relationships between different representation/sensing bases to see if they work well together to reconstruct the signal.)

Thank you.


texts: I ran into a strange bug while implementing a complex band-pass filter. The filter design process is based on the windowing method and seems to be working well, but when i apply the coefs of the filter on my time domain input samples and perform DFFT the carrier frequency is doubled. I can´t really figure out why is this happening but my guess is that i'm doing something bad during the convolution. The first plot shows the magnitude response of the complex band pass filter. The second plot shows FFT of my input 5kHz carrier signal, the third one shows the same data but with BP filter applied. Everything is happening in the complex plain.


texts:
What's the best method for calculating phase of complex nonsinusoidal signals?

I've tried using Hilbert transform, but I've found that my results change drastically depending on the amount of data I pass through the Hilbert transform. What I'm seeing is that the phase angle calculated depends on the maximum amplitude of the original signal. I don't understand how to deal with this.

I have motion capture data - tracking one marker on a foot in 3D space. During one part of my experiment, the foot moves less and at a lower velocity. During another part of the experiment, the amplitude and velocity increase. If I calculate phase for all of my data (both parts of the experiment) I get completely different results than if I calculate phase for each part of the experiment separately. When calculating phase for the entire trial, phase fluctuations are TINY (hovering near +pi or near -pi) except for the time when the foot amplitude and velocity reach their overall maximum - at that time phase fluctuates between pi and -pi. The maximum velocity/amplitude is about 1.8 times greater than the average during the first part of the experiment, but there is still sizeable motion during the first part of the experiment. If I calculate phase for each part of the experiment separately, I get more normal looking phase results, fluctuating between +- pi .

I ultimately want to calculate synchronicity between two such signals.
Can anyone explain this or provide advice?


texts: Direct sequence spread spectrum modulation

texts:
When the receiver gets the signal $m_t$ sent by the transmitter, the signal will have an additive interference $i_t$. Thereby signal received will be : $$r_t = m_t + i_t$$

I know that to get back the data from the signal, we have to multiply by the pseudo-noise sequence (PN).

That will take out the pseudo-noise from the signal $m_t$ allowing us to get the data.

What do we do with the noise $i_t$?

It has been multiplied by the PN but I supposed we cant take out the data from the $r_t$ until we eliminate the noise. I don't really see what kind of signal do we get multiplying the additive interference by the PN. I've read that we use an low pass filter but i dont really undestand why.


texts: How to calculate FFT of a non-uniformly sampled signal?

texts:
Is it possible to calculate the FFT of a signal that was not sampled uniformly? I have an MPU6050 connected to an Arduino uno.
I write the samples at the default frequency, approximately 20 milliseconds,
but it's not always constant (see the image). Also see this related question.

This is my first time I approach to the FFT, and I have a civil engineer background, so is my first time with signal analysis, never studied Signal Theory before.

My interest for the FFT is to define the best low pass filter (example apply a butterworth but I don’t know how to choose the filter order and cutoff frequency.

Thanks for your time, thanks for you patience. I tried to solve this with the code below, but I think is only rubbish.

M=csvread('DATALOG.TXT');
n=length(M);
cleaning=500;
n=(n-cleaning);

MatFilt= M(cleaning:n,:);
m=length(MatFilt);
A= MatFilt(2:m,1);
B=MatFilt(1:m-1,1);
C=A-B;
meanrate= mean(C);
f= 1000/meanrate;
T= sum(MatFilt);
TemporalSample= MatFilt(:,1);
Ax=MatFilt(:,2);
Axs=Ax*19.6/32768; %I converted the data in m/s^2
SignalLength= (MatFilt(m,1)-MatFilt(1,1))/1000;
plot(TemporalSample,Axs);
pause
title('Signal Corrupted with Zero-Mean Random Noise')
xlabel('time (milliseconds)');
NFFT = 2^nextpow2(m); % Next power of 2 from length of y
Y = fft(Axs,NFFT)/m;
f = f/2*linspace(0,1,NFFT/2+1);


texts: Derive the Forward Euler substitution for transfer function

texts:
In the book "The control handbook. Volume 1 " by Levine, the author shows that the transfer function:

can be aproximated and discretized in the transfer function:

using the forwar euler integral aproximation. In this way Levine shows that if i have a transfer function $G(s)=\frac{1}{s}$ i can simply substitute $s=\frac{z-1}{T}$ to apply the Forward Euler aproximation.
At this point the author doesn't show how to generalize this particular case to a generic transfer function $H(s)$. Are you able to show me ?


texts: Inconsistency Between Analysis and Simulation Results of Fourier Transform of a Sinc Function

texts:
First of all, hello.
This question is about a problem that I’ve faced during an attempt to obtain both time and frequency responses of a sinc function in MATLAB. The problem is an inconsistency between my calculations that I’ve done on the paper and simulation results that I’ve observed in MATLAB.
I’ve picked up $sinc(10 \cdot t)$ function and at first, I’ve taken its Fourier transform on paper. Afterwards, in order to have an experience in coding, I implemented this operation into MATLAB and plotted the functions both in the time and frequency domains.
However, in the amplitude spectrum in the frequency domain, the amplitude of the resultant rectangular pulse didn’t match with my results on the paper. According to my analysis results, the amplitude of the rectangular pulse is 0.1 but in MATLAB, I’ve got 10 instead of it.
My MATLAB code for simulation is given below:
t = -10: 0.01: 10;
Ts = 0.01;  % Sampling period in seconds
fs = 100;  % Sampling frequency in Hertz
function_1 = sinc(10 * t);
number_of_samples_1 = length(function_1);
number_of_samples_2 = power(2, nextpow2(number_of_samples_1));
f = linspace(-fs / 2, fs / 2, number_of_samples_2);
function_2 = fft(function_1, number_of_samples_2);
function_3 = fftshift(function_2);
function_4 = abs(function_3);
subplot(2, 1, 1);
plot(t, function_1, "r");
xlabel("Time (s)");
ylabel("Function-1 (t)");
title("Sinc(10.t) Function in the Time Domain");
subplot(2, 1, 2);
plot(f, function_4, "r");
xlabel("Frequency (Hz)");
ylabel("|Function-4 (f)|");
title("Amplitude Spectrum of Sinc(10.t) Function in the Frequency Domain");

My question is this: Why has this result occurred in contrast to my expectations and analysis?
I will be happy if you give feedback to my problem and question.
Thank you.


texts: Specifications (N) FIR Low Pass Resamplig, Hamming Window

texts:
I'm learning about DSP and I've doubts.

The resampling process can be do it with only one filter.

$L$: Increment Sample Rate Integer

$M$: Decrement Sample Rate Integer

$f_S$: Frequency Sampling

For Low Pass Filter:

The frequencies allowed are in the range: from $0$ until $f_{pass}$,

The frequencies denied are in the range: from $f_{stop}$ until $f_S/2$

$f_{pass}$: Highest pass frequency

$f_{stop}$: Lowest cut frequency

I'm simplifying my calculations as well:

I'm defining: $ 0 < \delta < 1$

$$f_{stop} = f_S/2\\
f_{pass} = \delta f_{stop} = \delta f_S/2$$

According to Known formula:

$$\Delta F = (f_{stop}-f_{pass})/f_S\\
f_c = (f_{stop}+f_{pass})/2
$$

For Hamming Window:

$$N = 3.3/\Delta F$$

Questions

Assuming $f_S = 24000$, $\delta = 0.97$ and the most demanding requirements for FIR filter calculus.

If $L=3$ and $M=2$; How can I to calculate (step by step) $N=?$
and $f_c = ?$
If $L=2$ and $M=3$; How can I to calculate (step by step) $N=?$
and $f_c = ?$

EDIT:

Case 1

$$f_{S1} = 24000\\f_{stop1} = f_{S1}/2 = 12000\\f_{pass1} = 0.97*f_{S1}/2 = 11640$$

$L=3$
$$f_{S2} = 72000\\f_{stop2} = f_{S2}/2 = 36000\\f_{pass2} = 0.97*f_{S2}/2 = 34920$$

$M=2$
$$f_{S3} = 36000\\f_{stop3} = f_{S3}/2 = 18000\\f_{pass3} = 0.97*f_{S3}/2 = 17460$$

Case 2

$$f_{S1} = 24000\\f_{stop1} = f_{S1}/2 = 12000\\f_{pass1} = 0.97*f_{S1}/2 = 11640$$

$L=2$
$$f_{S2} = 48000\\f_{stop2} = f_{S2}/2 = 24000\\f_{pass2} = 0.97*f_{S2}/2 = 23280$$

$M=3$
$$f_{S3} = 16000\\f_{stop3} = f_{S3}/2 =  8000\\f_{pass3} = 0.97*f_{S3}/2 =  5173$$


texts: Is this waterfall characteristic of QAM or OFDM (or something other) modulation?

texts:
I recorded this signal at 2.45GHz, it is coming from a camera.

Is this QAM or OFDM ?

Based on those side-lobes I would guess OFDM but I am not sure.

Also there seems to be an interesting preamble at the top.


texts: How can I find the transfer function of the following block diagram?

texts:
I've the following image and I want to find the transfer function from input $x(t)$ to output $y(t)$.
I know that I have to apply Laplace Transform, so the integrator becomes $\dfrac{1}{s}$, but I don't know what to do with the numbers $a$ and $b$. Should it still be the same?
If a and b still the same, then I found that $H(s) = \dfrac{Y(s)}{X(s)} = s^2-as-b$.


texts: how to get the location long and lat out of the raw measurements?

texts:
I'am making a research about GNSS
I used GNSS logger to collect measurements
but couldn't understand how these measurements turn into location of long and lat?
the picture below of the measurements attributes.


texts: Can a pole at the origin be considered a dominant pole?

texts:
I am studying control syetems, and I have studied that the dominant poles of a system are the ones closer to the imaginary axis.

Now, consider a system with a pole at the origin. Can this pole at the origin be considered a dominant pole?


texts:
I am working on a photo stitching application that uses multi-band blending. I need to get rid of unpleasant edges appearing at some places:

Here is the area of overlap (left - new image added to the mosaic, right - current mosaic contaning pixels of new image on background pixels to improve blending, middle - blending mask):

If I just compute weighted average between left and right image according to mask, the result is of course OK. However, this would leave a visible seam as two images have usually slightly different exposure.

So all three images need to be successively blurred to build a Gaussian pyramid - here is how one level of the pyramid looks like:

You can see that the top part of blending mask "touches" the border. The Gaussian blurring filter reflects on image borders and this causes inaccuracy in lower band.

I colored the images to make the problem more visible:

I am not sure about how to blur the mask so that it would deal with overlap area edges nicely.

Few suggestions:

change behavior of the blurring filter (how?)
extend area of blending mask so that one part never "touches" the border
update weighting masks so that resulting blend mask is more "edge-aware"

Any other suggestions/hints?


texts: Use Recursive Least Square as an "Actual" Filter?

texts:
I am currently studying recursive least square (RLS) and as far as I understand, the setup is (for example) that we are given a process with uncertain parameters/noise and we want to adapt parameters of a model of that process online.

To me, that sounds like a Luenberger observer, where the parameters of the underlying model are adapted online, depending on the model quality.

If this is true, a side question would be why this technique isn't usually called an observer, but a filter instead? However, my main question is:

Question: Can RLS also be used instead of a standard low-pass filter?

So, given only a signal from which I want to remove high-frequency components (like noise), can a RLS filter be used to do that as well? And if so, would it have any advantage over normal low-pass filters?


texts: Adaptive Line noise removal using IIR notch filters

texts: I'm trying to filter out line noise(60Hz and 120Hz) from a EEG signals received over a bluetooth link. I'm proposing to use a IIR notch filter to filter out the line-noise which varies w.r.t to distance between the bluetooth transmitter and receiver. I'm not able to establish a clear relationship between the Notch filter parameters and the magnitude of line noise so as to avoid over-attenuation or under-attenuation to avoid dips and peaks at these frequencies in the spectrum. The problem is Q-factor and stop-band Bandwidth which define the notch filter don't directly relate to the stop band attenuation. If this is not possible is there any other equally efficient adaptive stop-band filtering method to achieve the same ? How can this be implemented w.r.t line noise measurement and corresponding filter parameter adjustments ?


texts: Use of the Sobel Filter for Image Gradient [MATLAB]

texts:
I am using the Sobel filter for an RGB image. I have found two different ways to do that and the results look a little different. What is the difference between these two methods?

Method 1

 sbl=fspecial('Sobel');
 sbl_red=imfilter(Image(:,:,1));

Method 2

 [magred,dirred]=imgradient(Image(:,:,1),'Sobel')

Is there an advantage over one of the two methods? I am trying to find a structure of the image that is parallel to the x axis.


texts: How to Formulate a Constraint Which Ensures All Variables Have the Same Sign

texts:
I'm trying to include a constraint in my problem (to be solved by any convex optimization solver). Let {a,b,c,d ...} be a finite set of continuous variables. How to formulate a constraint which ensure all of these variables, in the same time interval, take either "positive or zero" or "negative or zero" values.

sincerely,

Mohan


texts:
I currently have written a code which takes in a raw file as input and does some audio processing and writes it to another different raw file.
The way I am currently inputting is
.\my_code_binary < input.raw > output.raw

as you can see, I am making the input.raw as stdin and output.raw as stdout for the execution my program.
fread(tmp, sizeof(short), channels * size_of_frame, stdin); // the way I am using the input.raw
fwrite(tmp, sizeof(short), channels * FRAME_SIZE, stdout); // the way I am using the output.raw

Now I want to make my program run real-time, as in, take my mic input as stdin and mic output as stdout. Any resources or code snippets will help me out, I am a beginner at audio processing in C.
I am using a Raspberry Pi 4


texts: Generating a signature to identify RF modulation types

texts:
Thanks to the Realtek 2832U I'm just discovering the world of Software Defined Radio.  There are a lot of signals I'm able to pick up with this device, but identifying them is not always easy and cannot be done automatically.

When scanning for radio transmissions, I can look at a waterfall plot and see whether a broadcast is AM, FM or a handful of digital modes I have learnt to identify, but this is slow and error prone.  I am wondering whether there is a way to automate this.

What I am after is an algorithm to examine the block of RF spectrum that has been digitised, to identify any signals along with their transmission modes.  For example, in a particular part of the spectrum it might identify three broadcasts, two that are FM and one that is a POCSAG data stream.

I am thinking that for this you would need to generate some sort of signature for each detected signal, and then compare that to a list of known modulation types.  Is it possible to generate a signature like this, taking into account the different amplitudes (signal strength) that result from varying distances to the transmitter?  The signature (i.e. an identifying number) would need to be the same for all instances of that signal - for instance, each consumer FM radio station would have the same signature.

I am not really sure where to start, or how to deal with the bandwidth ranges involved (some signals are narrow and close together, others are very wide.)  I'm also unsure whether you would have to pick out strong signals first and only process those (losing weaker ones) or whether you would examine the spectrum uniformly regardless of the power level at that particular frequency, discarding any frequencies that matched the signature for "no signal present."

The reason I am asking is that I would like to write a program that can tune the radio to a given frequency and then display all the transmissions it picks up, decoding them all automatically where possible.  In order to do this, it must be able to identify signals in enough detail to demodulate them.


texts:
What are the statistics of the cepstrum of gaussian white noise?

\begin{align}\newcommand{\Nfft}{ {N_{\mathrm{FFT}} }}\DeclareMathOperator{\FFT}{FFT}\DeclareMathOperator{\IFFT}{IFFT}
x_i &\sim \mathcal{N}(0, \sigma^2)& i&\in\{0,\ldots,\Nfft-1\}\\[1.2em]
%
X   &= \FFT(x)\\
X_i &\sim \mathcal{N}\left(0, \frac{\sigma^2}{\Nfft}\right)\\[1.2em]
%
\left\lvert X_i\right\rvert^2 &\sim \chi^2\\[1.2em]
%
C   &= \IFFT\left(\log{\left\lvert X\right\rvert^2}\right)\\
C_i &\sim ~\ ???
\end{align}

What if I apply a window function?

\begin{align}
X &= \FFT(w \cdot x)\\
  &= \FFT(w)*\FFT(x) \\[1.2em]
%
C &= \IFFT\left(\left\lvert X\right\rvert^2\right)\\
C_i &\sim ???
\end{align}

Even an approximation would be nice. I'm working on a detector and understanding the noise structure is important for my approach.


texts: How to calculate the Fourier Transform of a solvable chaos waveform?

texts:
Recently I am stucking in frequency estimation of a solvable chaos waveform.
Its local analytic expression in time domain is
$$
z(t)=s_m(u_m-s_m)e^{\beta(t-mT)}\cos(\omega_0 t+\varphi),mT\leq t<(m+1)T
$$
where $u_m \sim U(-1,1)$ , $s_m=sgn[u_m]$ and $s_m\in[+1,-1]$.
$T=2\pi/\omega_0$ and $\beta=T \mathrm{ln}2$. And $m = 0,1,2,...$.

And its waveform is shown in following ($\omega_0 = 200\pi$)

Its spectrum shows spikes at base frequency $f_0 =2\pi/\omega_0$ and its harmonics (here we set the base frequency as 100Hz). Besides, the spectrum has regularly repeated shape which implies some kind of determined solvable expression of Fourier Transform.

To begin with, I tried to calculate the FT directly, as
$$
Z(\omega)=\int_{-\infty}^\infty z(t)e^{j\omega t}\mathrm{d}t\\
=\sum_{m=-\infty} ^{\infty}s_m(u_m-s_m) \int_{mT}^{(m+1)T} e^{\beta(t-mT)}\cos(\omega_0 t+\varphi) e^{j\omega t}\mathrm{d}t\\
=\frac 1 2 \sum_{m=-\infty} ^{\infty}s_m(u_m-s_m) \int_{mT}^{(m+1)T} e^{\beta(t-mT)} \{e^{j[(\omega_0+\omega) t+\varphi]}+e^{-j[(\omega_0-\omega) t+\varphi]}\}\mathrm{d}t
$$

However, the integeration is limted in a window with width of $T$ and it is hard to calculate the FT from it.

Then I tried another way to prove the existence of spikes at base frequency and its harmonics. Here, I calcute the integration of $z(t)\sin(n\omega_0 t)$ as
$$
a_0(n\omega_0)=\int_{-\infty}^\infty z(t)\sin (n\omega_0 t)\mathrm{d}t\\
=\sum_{m=-\infty} ^{\infty}s_m(u_m-s_m) \int_{mT}^{(m+1)T} e^{\beta(t-mT)}\cos(\omega_0 t+\varphi) \sin (n\omega_0 t)\mathrm{d}t\\
=\sum_{m=-\infty} ^{\infty}s_m(u_m-s_m) \int_{0}^{T} e^{\beta t}\cos(\omega_0 t+\varphi) \sin (n\omega_0 t)\mathrm{d}t\\
=\int_{0}^{T} e^{\beta t}\cos(\omega_0 t+\varphi) \sin (n\omega_0 t)\mathrm{d}t\sum_{m=-\infty} ^{\infty}s_m(u_m-s_m)
$$

where $n=1,2,3...$ .
And
$$
b_0(n\omega_0)=\int_{-\infty}^\infty z(t)\cos (n\omega_0 t)\mathrm{d}t\\
=\int_{0}^{T} e^{\beta t}\cos(\omega_0 t+\varphi) \cos (n\omega_0 t)\mathrm{d}t\sum_{m=-\infty} ^{\infty}s_m(u_m-s_m)
$$

And the amplitude of spectrum at $n\omega_0$ is
$$
c(n\omega_0)=\sqrt{a_0(n\omega_0)^2+b_0(n\omega_0)^2}
$$

Still, it is difficult to prove the existence of delta funcions in the spectrum from the amplitude at base frequency.

Hence, how to prove the frequency properties shown in Figure 2 ? Is there any possiblity to calculate the FT of $z(t)$?


texts: How to verify if two images have the same object/content?

texts:
I am truly new to computer vision field, but it is fascinating me! I have now a challenge in my hands and I am seeking for mentors/advisers to give me some follow up.

My project is: From a picture of a video-game cover, search that pic through a video-game cover database, and if there is a very good match the app will return a string with the name of the videogame and the platform.

Problem example:

Take a photo of a cover similar to this:

Cover match to this in the database:

App gives string: "Fifa 12 Playstation 2"

In my preliminary research I found that I should save in my cover database, the name of the game, the platform, the URL for the cover, and the image features (keypoints and descriptors).

I am using SURF features detector/extractor.

There are some concepts that I am still confused about... I am not looking for similarity right? I just need to observe if there are some sort of good keypoints-paired/matched, right? Because in the examples images mentioned before I get "img1 - 1087 features (query image), img2 - 1755 features - 30 % - SIMILARITY - 321/333 inliers/matched" What are the inliers? My calculation of similarity seems wrong to me... I would say that these two images are like 70 % look alike...

ps: I am using Python so...

Results so far:

Thanks for your time and help. Sorry If I could not explain any better my problems/concerns.


texts:
Let $X(t)$ and $Y(t)$ be two orthogonal processes with power spectral densities
$$S_{xx}(f) = S_{yy}(f)=\begin{cases}
1-\lvert f\rvert, & \lvert f\rvert<1
\\[1ex]
0,& \text{otherwise}
\end{cases}$$

Define a new process $Z(t) = Y(t)-X(t-1)$. Determine and sketch the power spectral density $S_{zz}(f)$.

The solution is the following:
The autocorrelation function of the process is
$$R_{zz}= E \bigg\{\big[Y(t-\tau)-X(t+\tau-1)\big]\big[Y(t)-X(t-1)\big] \bigg\} $$
And now comes the part which I do not understand at all.
\begin{align}
R_{zz}&=R_{yy}(\tau)-R_{yx}(\tau+1)-R_{xy}(\tau-1)+R_{xx}(\tau)\\
&=R_{yy}(\tau) + R_{xx}(\tau)\ \quad\text{since}\quad R_{yx} = R_{xy} = 0 \quad\text{from orthogonality.}
\end{align}

Therefore, $S_{zz}(f)=2S_{yy}(f)=2S_{xx}(f)$ as shown below.

Can anyone explain it like I am five?


texts: why is the plot of the PSD of this autocorrelation function oscillating?

texts:
I was doing an exercise trying to compute the PSD of the following discrete autocorrelation function:

\begin{equation}
r_x(m) =(1/2)^{|m|}
\end{equation}

In this simple  matlab code that I came up with to compute the PSD, I'm getting something that I don't quite understand, why is my frequency response oscillating?.If I plot the abs(Sx) i get a smooth plot but , the dft of a real and even sequence should be real and even as well in my understanding . Why I'm not getting a real frequency response and instead this oscillating response appears? Also, matlab is ignoring the imaginary part of the PSD computed as it says in the command window. What is going on  here?

N=100;
m=-N:1:N;
x=(0.5).^abs(m);
subplot(2,1,1)
stem(m,x)
Nfft=512;
Sx=fft(x,Nfft);
subplot(2,1,2)
plot((0:1:Nfft-1)/Nfft,Sx);

Here are  the figures:


texts:
I have just started working on FFT and its related phase factors.

I have read for an N=2 point FFT, We have two phase factors, one on the positive axis and one on negative axis.

similarly for an 4 point FFT, we have 4 phase factors as theta=(pi/2), at

 0'  W 40=1
-90' W 41=-j
180' W 42=-1
 90' W 43= j

here we can see the values that are 180 degree apart are negative to each other.

Here I wanted to know how it works for an 32 point FFT.

As theta=(pi/16) do we have 16 phase factors on the positive axis and 16 phase factors on the negative axis with the interval of (pi/16) ?

so will it follow the pattern as I mentioned below?

 0'       W16 0 =1
-11.25'   W16 1 =0.98 -j0.195
-22.5'    W16 2 =0.92 -j-0.38
and so on
....
....
till
-168.75'  W16 15=
and the negative parts of the above values from **W16 16** to **W16 32**

Am I correct? please verify the above details and let me know how the phase factors works for an 32 point FFT

Thanks in advance


texts:
I must strengthen my knowledge in the mathematical properties of signals, and would like to know the following thing:

if I have a complex signal, written in magnitude and phase, it's like if I have 1 signal that holds 2 signal himself? I have this question because I have to do the inverse transformation of a signal in the frequency domain , that is complex and I have the magnitude and phase plot. Is it right that my phase, when I do the inverse- transform retards my signal ( the inverse transform of the magnitude) in the time domain?
Thank you so much.. I have not on focus how the phase works yet

Thank you for the answer, I know what you have written... my question is:
If I have one complex signal in the magnitude and phase form, do I have to consider magnitude and phase as two different signals? Look

The signal is written in  magnitude with a complex exponential with the phase... when I do the inverse transform the exponential is in the time domain a delay ... and this it feels strange that one signal bring with him a delay ... it's like if phase is another signal , that retards the inverse transform of the magnitude.
I hope I have explained my question in the right way