Find peaks in 2D-Spectra and calculate peak-volume?

Share icon

Dear all,

I have some three-dimensional peaks on a plane (two-dimensional spectrum) and wondered whether there is an existing macro to automatically find these peaks and calculate their peak volumes?

Maybe someone stumbled over this in the past or is actively working with something like this?

If not, what would be the best general way to approach this problem?

Especially I wondered what would be the best (mathematically sound) way to calculate peak volumes?

Thanks for any help!

Best regards,
Peter
Well, here is a function I wrote for another customer that uses ImageAnalyzeParticles to find spots in a 2D wave. It then fits a two-dimensional Gaussian to each spot it finds. The coefficients for the 2D Gaussian are stored in a multi-column wave, one row per spot. You would then need to figure out how to compute the volume for each spot.

Note that this won't do the right thing if the spots overlap. If there are overlapping spots, the problem gets much more complicated.


Function DoParticleAnalysis(imageW)
	Wave imageW
	
	// First, find the spots
	ImageThreshold/I/M=0/T=800/Q imageW
	Wave M_ImageThresh
	ImageAnalyzeParticles /E/W/Q/M=3/A=3 stats, M_ImageThresh
	
	Wave W_ImageObjArea
	Variable nSpots = DimSize(W_ImageObjArea, 0)
	
	// Now we have an estimate of the spot locations- we can now loop through the spots and do fits
	Wave W_xmin, W_xmax, W_ymin, W_ymax
	WaveStats/Q W_ImageObjArea
	Variable maxArea = V_avg+3*V_sdev	// try to eliminate obviously bad "spots"
	Variable i
	Variable V_fitOptions=4
	Variable V_fitError
	Make/N=(nSpots, 7)/O $(nameOfWave(imageW)+"fits")
	Wave results = $(nameOfWave(imageW)+"fits")
	results = nan
	for (i = 0; i < nSpots; i += 1)
		if (W_ImageObjArea[i] > maxArea)
			continue		// this spot looks like a bad one- go around for another
		endif
		
		V_fitError = 0
		// Fit a 2D gaussian to just one spot. The -3, +3 bit is to include sufficient background to get a reliable fit
		CurveFit/N/Q Gauss2d imageW[W_xmin[i]-3, W_xmax[i]+3][W_ymin[i]-3, W_ymax[i]+3]
		Wave W_coef
		if (V_fitError == 0)
			results[i][] = W_coef[q]
		endif
	endfor
end


John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

September 8, 2015 at 11:46 am - Permalink
The image analysis expert has pointed out to me that the function has a hard-coded threshold value of 800. You will need to adjust the code for your data.

John Weeks
WaveMetrics, Inc.
support@wavemetrics.com

September 9, 2015 at 10:21 am - Permalink