Standard deviation in ROI of image

Share icon

Hello all,

I am calculating signal to noise ratio of spots in a bunch of images. I am specifying the ROI for the signal and an ROI for the noise using a mask wave where 0 is outside the ROI and 1 is inside.

Calculating the mean can be done using a nice one-liner from AG that I found in another thread. I'm curious if there is a similar one line MatrixOp for getting the standard deviation.

Here is my current code:

// snip
		m0 = 0 // m0 is mask wave for signal ROI
		MakeMaskMean(m0,xx,yy) // external function
		MatrixOP/O aa = sum(m0 * imgMat) / sum(m0) // calculates the mean
		qW[i][0] = aa[0]
		m1 = 0 // m1 is mask wave for noise ROI
		MakeMaskSD(m1,xx,yy) // external function
		MatrixOP/O bb = sum(m1 * imgMat) / sum(m1) // get mean
		Duplicate/O/FREE imgMat,m2
		m2 -=bb[0] // subtract mean from image
		m2[][] *=m1[p][q] // apply the mask
		MatrixOP/O cc = sqrt(sum(m2 * m2) / sum(m1)) // calculate SD
		qW[i][1] = cc[0]
// snip

Getting SD is 5 lines and I wonder if it would be faster to push the values in the ROI to a 1D wave and do WaveStats, or if there is some MatrixOP magic I am missing.
Thanks for any help.
Curious why use Matrixop instead of Imagestats?

Imagestats returns average and standard deviation and takes an ROI mask all in one line.

Andy

March 29, 2017 at 08:43 am - Permalink
Argh I did not realise that ImageStats took an ROI wave. Thank you!!

March 29, 2017 at 08:55 am - Permalink
First let's take a look at how you could change your code snippet:


		m0 = 0 // m0 is mask wave for signal ROI
		MakeMaskMean(m0,xx,yy) // external function
		MatrixOP/O aa = sum(m0 * imgMat) / sum(m0) // calculates the mean
		qW[i][0] = aa[0]
		m1 = 0 // m1 is mask wave for noise ROI
		MakeMaskSD(m1,xx,yy) // external function
		MatrixOP/O bb = sum(m1 * imgMat) / sum(m1) // get mean
		Variable bb0=bb[0]
		// There is no need to waste time on this duplication; MatrixOP takes care of it:
		// Duplicate/O/FREE imgMat,m2
		// m2 -=bb[0] // subtract mean from image
		// m2[][] *=m1[p][q] // apply the mask
		// MatrixOP/O cc = sqrt(sum(m2 * m2) / sum(m1)) // calculate SD
		MatrixOP/O cc=sqrt(sum(magsqr((imgMat-bb0)*m1))/sum(m1))
		qW[i][1] = cc[0]


In practice, depending on the size of your images and ROI's, ImageStats may be more efficient because in IP7 ImageStats is automatically multithreaded.

Also note that when you work with ROI in MatrixOP, included pixels correspond to a value of 1 in the mask while ImageStats requires that included pixels are marked by 0 in the ROI wave.


A.G.
WaveMetrics, Inc.

March 30, 2017 at 10:40 am - Permalink
I am now using ImageStats for this, but thank you for the input into my code. I learned two new tricks there. How to subtract a constant in MatrixOP and also the magsqr command. I'm not sure which approach would be faster, the ROIs I am using are small: tens of pixels.

Thanks!

March 30, 2017 at 01:08 pm - Permalink