package ij.process;
import ij.measure.Calibration;
import java.util.Arrays;
public class FloatStatistics extends ImageStatistics {
public FloatStatistics(ImageProcessor ip) {
this(ip, AREA+MEAN+MODE+MIN_MAX, null);
}
public FloatStatistics(ImageProcessor ip, int mOptions, Calibration cal) {
this.width = ip.getWidth();
this.height = ip.getHeight();
setup(ip, cal);
double minT = ip.getMinThreshold();
double minThreshold,maxThreshold;
boolean limitToThreshold = (mOptions&LIMIT)!=0;
if (!limitToThreshold || minT==ImageProcessor.NO_THRESHOLD) {
minThreshold=-Float.MAX_VALUE;
maxThreshold=Float.MAX_VALUE;
} else {
minThreshold=minT;
maxThreshold=ip.getMaxThreshold();
}
if (limitToThreshold)
saveThreshold(minThreshold, maxThreshold, cal);
getStatistics(ip, minThreshold, maxThreshold);
if ((mOptions&MODE)!=0)
getMode();
if ((mOptions&ELLIPSE)!=0 || (mOptions&SHAPE_DESCRIPTORS)!=0)
fitEllipse(ip, mOptions);
else if ((mOptions&CENTROID)!=0)
getCentroid(ip, minThreshold, maxThreshold);
if ((mOptions&(CENTER_OF_MASS|SKEWNESS|KURTOSIS))!=0)
calculateMoments(ip, minThreshold, maxThreshold);
if ((mOptions&MEDIAN)!=0)
getMedian(ip, minThreshold, maxThreshold);
if ((mOptions&AREA_FRACTION)!=0)
calculateAreaFraction(ip);
}
void getStatistics(ImageProcessor ip, double minThreshold, double maxThreshold) {
double v;
float[] pixels = (float[])ip.getPixels();
nBins = ip.getHistogramSize();
histMin = ip.getHistogramMin();
histMax = ip.getHistogramMax();
histogram = new int[nBins];
double sum = 0;
double sum2 = 0;
byte[] mask = ip.getMaskArray();
double roiMin = Double.MAX_VALUE;
double roiMax = -Double.MAX_VALUE;
for (int y=ry, my=0; y<(ry+rh); y++, my++) {
int i = y * width + rx;
int mi = my * rw;
for (int x=rx; x<(rx+rw); x++) {
if (mask==null || mask[mi++]!=0) {
v = pixels[i];
if (v>=minThreshold && v<=maxThreshold) {
if (v<roiMin) roiMin = v;
if (v>roiMax) roiMax = v;
}
}
i++;
}
}
min = roiMin; max = roiMax;
if (histMin==0.0 && histMax==0.0) {
histMin = min;
histMax = max;
} else {
if (min<histMin) min = histMin;
if (max>histMax) max = histMax;
}
binSize = (histMax-histMin)/nBins;
double scale = nBins/(histMax-histMin);
int index;
pixelCount = 0;
for (int y=ry, my=0; y<(ry+rh); y++, my++) {
int i = y * width + rx;
int mi = my * rw;
for (int x=rx; x<(rx+rw); x++) {
if (mask==null || mask[mi++]!=0) {
v = pixels[i];
if (v>=minThreshold && v<=maxThreshold && v>=histMin && v<=histMax) {
pixelCount++;
sum += v;
sum2 += v*v;
index = (int)(scale*(v-histMin));
if (index>=nBins)
index = nBins-1;
histogram[index]++;
}
}
i++;
}
}
area = pixelCount*pw*ph;
mean = sum/pixelCount;
umean = mean;
calculateStdDev(pixelCount, sum, sum2);
}
void getMode() {
int count;
maxCount = 0;
for (int i = 0; i < nBins; i++) {
count = histogram[i];
if (count > maxCount) {
maxCount = count;
mode = i;
}
}
dmode = histMin+mode*binSize;
if (binSize!=1.0)
dmode += binSize/2.0;
}
void calculateMoments(ImageProcessor ip, double minThreshold, double maxThreshold) {
float[] pixels = (float[])ip.getPixels();
byte[] mask = ip.getMaskArray();
int i, mi;
double v, v2, sum1=0.0, sum2=0.0, sum3=0.0, sum4=0.0, xsum=0.0, ysum=0.0;
for (int y=ry,my=0; y<(ry+rh); y++,my++) {
i = y*width + rx;
mi = my*rw;
for (int x=rx; x<(rx+rw); x++) {
if (mask==null || mask[mi++]!=0) {
v = pixels[i]+Double.MIN_VALUE;
if (v>=minThreshold && v<=maxThreshold) {
v2 = v*v;
sum1 += v;
sum2 += v2;
sum3 += v*v2;
sum4 += v2*v2;
xsum += x*v;
ysum += y*v;
}
}
i++;
}
}
double mean2 = mean*mean;
double variance = sum2/pixelCount - mean2;
double sDeviation = Math.sqrt(variance);
skewness = ((sum3 - 3.0*mean*sum2)/pixelCount + 2.0*mean*mean2)/(variance*sDeviation);
kurtosis = (((sum4 - 4.0*mean*sum3 + 6.0*mean2*sum2)/pixelCount - 3.0*mean2*mean2)/(variance*variance)-3.0);
xCenterOfMass = xsum/sum1+0.5;
yCenterOfMass = ysum/sum1+0.5;
if (cal!=null) {
xCenterOfMass = cal.getX(xCenterOfMass);
yCenterOfMass = cal.getY(yCenterOfMass, height);
}
}
void getCentroid(ImageProcessor ip, double minThreshold, double maxThreshold) {
float[] pixels = (float[])ip.getPixels();
byte[] mask = ip.getMaskArray();
double count=0.0, xsum=0.0, ysum=0.0, v;
int i, mi;
for (int y=ry,my=0; y<(ry+rh); y++,my++) {
i = y*width + rx;
mi = my*rw;
for (int x=rx; x<(rx+rw); x++) {
if (mask==null||mask[mi++]!=0) {
v = pixels[i];
if (v>=minThreshold && v<=maxThreshold) {
count++;
xsum+=x;
ysum+=y;
}
}
i++;
}
}
xCentroid = xsum/count+0.5;
yCentroid = ysum/count+0.5;
if (cal!=null) {
xCentroid = cal.getX(xCentroid);
yCentroid = cal.getY(yCentroid, height);
}
}
void calculateAreaFraction(ImageProcessor ip) {
int sum = 0;
int total = 0;
float t1 = (float)ip.getMinThreshold();
float t2 = (float)ip.getMaxThreshold();
float v;
float[] pixels = (float[])ip.getPixels();
boolean noThresh = t1==ImageProcessor.NO_THRESHOLD;
byte[] mask = ip.getMaskArray();
int i, mi;
for (int y=ry,my=0; y<(ry+rh); y++,my++) {
i = y*width + rx;
mi = my*rw;
for (int x=rx; x<(rx+rw); x++) {
if (mask==null||mask[mi++]!=0) {
v = pixels[i];
total++;
if (noThresh) {
if (v!=0f) sum++;
} else if (v>=t1 && v<=t2)
sum++;
}
i++;
}
}
areaFraction = sum*100.0/total;
}
void getMedian(ImageProcessor ip, double minThreshold, double maxThreshold) {
if (pixelCount==0) {
median = Double.NaN;
return;
}
float[] pixels = (float[])ip.getPixels();
float[] pixels2 = new float[pixelCount];
byte[] mask = ip.getMaskArray();
int i, mi;
float v;
int count = 0;
for (int y=ry,my=0; y<(ry+rh); y++,my++) {
i = y*width + rx;
mi = my*rw;
for (int x=rx; x<(rx+rw); x++) {
if (mask==null||mask[mi++]!=0) {
v = pixels[i];
if (v>=minThreshold && v<=maxThreshold) {
if (count==pixels2.length) {
median = Double.NaN;
return;
}
pixels2[count++] = v;
}
}
i++;
}
}
Arrays.sort(pixels2);
int middle = pixels2.length/2;
if ((pixels2.length&1)==0) median = (pixels2[middle-1] + pixels2[middle])/2f;
else
median = pixels2[middle];
}
}