Can someone verify my code for calculating a Gaussian kernel?

[eemans]

Here's the paper that discusses the convolution method.

[altenbach]

I cannot load the arithmetic" express VI, It seems there are some signalexpress components (CG .......vi) missing.

Where did you get that one? What does it do?

[eemans]

I found the arithmetic expresss vi. by going to Functions>>Signal Express>>Processing>>Analog Signals>>Arithmetic. 

LabVIEW help had this to say about the function:

Arithmetic (Time Domain) (LabVIEW SignalExpress)

Performs arithmetic operations on two signals. The available operations change depending on the type of the input signals you select.

To select the correct type of operation, select the first signal to process from the Input signal 1 pull-down menu. The Input signal 2 pull-down menu displays only the list of compatible signals, and the step displays the available operations.

Dialog Box Options

Parameter	Description
Input Signals	Displays the two input signals.
Resulting Signal	Displays the signal that results from the arithmetic operation.
Input	Contains the following options: Input signal 1:—Specifies the first input signal. Input signal 2:—Specifies the second input signal.
Configuration	Contains the following option: Operation Setup—Contains the following options: Operation—Specifies to add, subtract, multiply, or divide the signals or compute the RMS sum of the signals. The default is to add the signals. Ignore timestamp—Ignores eventual differences in timestamps by forcing the second signal timestamp to be equal to the first signal timestamp. Interpolate if needed—Resamples the signals to align the samples. Optimized for single record—Optimizes the resampling operation for a single record. Interpolation type—Contains the following options: Coerce—Sets each output sample value to equal the input sample value that is closest to it in time. Linear—Sets each output sample value to be a linear interpolation between the two input samples that are closest to it in time. Spline—Uses the spline interpolation algorithm to compute the resampled values. Filter based—Uses an interpolation method based on the convolution of the signal with an FIR filter. Resulting interval—Specifies if the resulting signal covers the Common or Global time interval of the signals.
Filter Setup	Available when you select Filter based from Interpolation type. Contains the following option: FIR Filter Specification—Contains the following options: Normalized bandwidth—The normalized cut-off frequency of the FIR filter to use. The default is 0.4000. Alias rejection (dB)—The minimum stopband attenuation of the FIR filter to use. The default is 80 dB.

The operation I'm using is "add", I am NOT ignoring time stamps, I'm using linear interpolation, optimized for a single record, and the two inputs are time waveforms

[altenbach]

OK, signalexpress != LabVIEW and I don't have it. 😉

Still I undestand what you are trying to do, and you are doing this orders of magnitude more complicated than needed. You create tons of waveforms, one for each spike, and then later add them all together. This is a lot of work, memory inefficient and it is slow.

Here's the trivially simple solution using plain convolution.

You simply create an array containing the spikes, then convolute it with the Gaussian kernel. Compare the speed to your solution and you'll be amazed. (Since I create a real histogram for the time points, you could even reduce the time resolution drmatically and still get about the right answer).

Let me know if you have any questions.

Message Edited by altenbach on 05-29-2009 04:04 PM

[eemans]

The convolution method is about 10x faster than the addition method.  However, the addition method is arguably more accurate given that the kernels are already aligned at the peak times.  I say arguably because the convolution method could be made more accurate (at the cost of speed) by decreasing the delta t.  Anyway, I've chosen to use the convolution method becasue of it's already widely accepted use.  Thank you for showing me how to set that up. 

I think I finally came up with a method of generating the kernel function that we can both agree upon.   It calculates the kernel without any loops or special scaling at the end.  The area under the curve always comes out to 1 (.9999) no matter what the width or resolution is set to.    Could you please take a look at the attached VI and let me know what you think.  

[eemans]

Could someone who is familiar with Gaussian kernel functions please take a look at my code and let me know if I'm calculating it correctly?  I've attached my VI for generating the kernel.  You'll have to set the width and resolution on the Front Panel before it can be run.

[altenbach]

Yes, it looks good.

You probably should keep the FP closed during run, else the waveform graph update would add to the cost of running it. You don't need the waveform graph except for debugging.

If this subVI is used often with similar or identical inputs during a run, you might want to add a feedback node to cache the result. Now simply check if the inputs are the same from the previous run and use the cached ouput if they are. (see image, the false case simply connects the output of the SR to the output tunnel of the case structure).

Since your data is symmetric, you need in principle to only calculate one half of the curve, then invert and merge. It's probably not worth for this simple equation, though. It would also complicate the code, because it is slightly more involved for odd numbered outputs.

Message Edited by altenbach on 06-08-2009 12:25 PM

[Ken_Brooks]

Hi, Altenbach. I happened across this thread, searching for information on Gaussian kernels myself. And the thing you are talking about, "width at half height", would be just the thing for me. If it's handy, would you mind showing me a version that takes that for a parameter?

I'm actually needing to compute a 2-D Gaussian convolution kernel -- but if I understand right, that's just the outer product of a 1D kernel with itself. True?

Thanks,
Ken