Detecting peaks in images

I got a large set of infrared images of seeds, their sizes vary slightly. And I would like to find them (in thefastest possible way). Below i show zoomed in details of the images i process. After a first noise removal and blob filter this is what i have : 在这里输入图像描述

The bright white is just direct reflection of the IR lamp, white pixels never combine (stretch out) over multiple seeds.

To make it more clear i placed a letter on some individual seeds. 在这里输入图像描述

The problems i have:

  • A is a single seed (dirt on the seed) generates a slight dark line.
  • B the nearby X close to it is at its darkest intersection its still brighter as some other seeds (cannt change brightnes or remove if gray value is lower then a certain value.
  • C those are 3 seeds close to each other.
  • The smallest visible seeds above should not becomme even smaller.
  • I'm not making use of mathlab or openCV, as i work directly with locked image and memory data. i can acces pixel data by array or simple getpixel / putpixel commands. I wrote my own graphics library which is fast enough for live camera data, processing speed currently at around 13ms at around 25ms i enter stream processing lag

    I wonder how to separate those 'cloudy' blobs better. I'm thinking to find local maxima over a certain pixels range..but that should see A as one seed, while on B find that B and X are not connected. So I'm not sure here, how such a local peek function or another function should like like. Although i code in C# i looked at other c++ functions as well like dilate etc but thats not it. I also wrote a function to check for slope degree (like if it was a mountain height image) but that couldnt devide areas B and C.

    Ok well i made different slope detection code, now i dont look for some degree but just the tilting point over a small range, it works nice on X axis.. but essentially i think it should work on both X and Y here's the new result : 在这里输入图像描述 It can resolve issue A and B !!!

    However it wouldnt differentiate between seeds who are aligned in a vertical row, and and it causes small white noise (non connected lines). at places where is there is nearly nothing to detect. I'm not yet sure on how to do the same (combined) over Y axis to get the tops then erase stuff from a certain distance of the top.. (to seperate).

    Using this code just showing the pixel operations of it.

     for (int y =  raw.Height;y>5; y--)
            {slopeUP = false;
                int[] peek = new int[raw.Width];
                for (int x = raw.Width; x> 7; x--)
                {
                    int a = raw.GetPixelBleu(x, y);
                    int b = raw.GetPixelBleu(x - 1, y);
                    int c = raw.GetPixelBleu(x - 2, y);
                    int d = raw.GetPixelBleu(x - 11, y);
                    int f = raw.GetPixelBleu(x - 12, y);
    
                    if ((f + d) > (a + b))slopeUP  = true;
    
                    if ((f + d) < (a + b))
                    {
                        if (slopeUP)
                        {
                            slopeUP = false;
                            peek[x - 6] = 10;
                            //raw.SetPixel(x, y, Color.GreenYellow);
                        }
                        else peek[x - 6] = 0;
                    }
    
                }
                for (int x = raw.Width; x > 7; x--) 
                { if (peek[x-1] > 5) raw.SetPixel(x, y, Color.Lavender); }
         }
    

    So, as far as speed, I am only going off of the image you posted here... on which everything runs blazing fast because it is tiny. Note that I padded the image after binarizing and never un-padded, so you will want to either un-pad or shift your results accordingly. You may not even want to pad, but it allows detection of cut off seeds.

    Overview of pipeline: removeSaturation>>gaussian blur>>binarize>>padd>>distanceTransform>>peaks>>clustering

    That being said here is my code and results:

    void drawText(Mat & image);
    void onMouse(int event, int x, int y, int, void*);
    Mat bruteForceLocalMax(Mat srcImage, int searchRad);
    void zoomPixelImage(Mat sourceImage, int multFactor, string name, bool mouseCallback);
    Mat mergeLocalPeaks(Mat srcImage, int mergeRadius);
    Mat image;
    bool debugDisplays = false;
    
    
    int main()
    {
        cout << "Built with OpenCV " << CV_VERSION << endl;
        TimeStamp precisionClock = TimeStamp();
    
        image = imread("../Raw_Images/Seeds1.png",0);
    
    
        if (image.empty()) { cout << "failed to load image"<<endl; }
        else
        {
            zoomPixelImage(image, 5, "raw data", false);
            precisionClock.labeledlapStamp("image read", true);
    
            //find max value in image that is not oversaturated
            int maxVal = 0;
            for (int x = 0; x < image.rows; x++)
            {
                for (int y = 0; y < image.cols; y++)
                {
                    int val = image.at<uchar>(x, y);
                    if (val >maxVal && val !=255)
                    {
                        maxVal = val;
                    }
                }
            }
    
            //get rid of oversaturation regions (as they throw off processing)
            image.setTo(maxVal, image == 255);
    
            if (debugDisplays)
            {zoomPixelImage(image, 5, "unsaturated data", false);}
            precisionClock.labeledlapStamp("Unsaturate Data", true);
    
            Mat gaussianBlurred = Mat();
            GaussianBlur(image, gaussianBlurred, Size(9, 9), 10, 0);
            if (debugDisplays)
            {zoomPixelImage(gaussianBlurred, 5, "blurred data", false);}
            precisionClock.labeledlapStamp("Gaussian", true);
    
            Mat binarized = Mat();
            threshold(gaussianBlurred, binarized, 50, 255, THRESH_BINARY);
            if (debugDisplays)
            {zoomPixelImage(binarized, 5, "binarized data", false);}
            precisionClock.labeledlapStamp("binarized", true);
    
            //pad edges (may or may not be neccesary depending on setup)
            Mat paddedImage = Mat();
            copyMakeBorder(binarized, paddedImage, 1, 1, 1, 1, BORDER_CONSTANT, 0);
            if (debugDisplays)
            {zoomPixelImage(paddedImage, 5, "padded data", false);}
            precisionClock.labeledlapStamp("add padding", true);
    
            Mat distTrans =  Mat();
            distanceTransform(paddedImage, distTrans, CV_DIST_L1,3,CV_8U);
            if (debugDisplays)
            {zoomPixelImage(distTrans, 5, "distanceTransform", true);}
            precisionClock.labeledlapStamp("distTransform", true);
    
            Mat peaks = Mat();
            peaks = bruteForceLocalMax(distTrans,10);
            if (debugDisplays)
            {zoomPixelImage(peaks, 5, "peaks", false);}
            precisionClock.labeledlapStamp("peaks", true);
    
            //attempt to cluster any colocated peaks and find the best clustering count
            Mat mergedPeaks = Mat();
            mergedPeaks = mergeLocalPeaks(peaks, 5);
            if (debugDisplays)
            {zoomPixelImage(mergedPeaks, 5, "peaks final", false);}
            precisionClock.labeledlapStamp("final peaks", true);
    
            precisionClock.fullStamp(false);
            waitKey(0);
        }
    }
    
    void drawText(Mat & image)
    {
        putText(image, "Hello OpenCV",
                Point(20, 50),
                FONT_HERSHEY_COMPLEX, 1, // font face and scale
                Scalar(255, 255, 255), // white
                1, LINE_AA); // line thickness and type
    }
    
    void onMouse(int event, int x, int y, int, void*)
    {
        if (event != CV_EVENT_LBUTTONDOWN)
            return;
    
        Point pt = Point(x, y);
        std::cout << "x=" << pt.x << "t y=" << pt.y << "t value=" << int(image.at<uchar>(y,x)) << "n";
    
    }
    
    void zoomPixelImage(Mat sourceImage, int multFactor, string name, bool normalized)
    {
        Mat zoomed;// = Mat::zeros(sourceImage.rows*multFactor, sourceImage.cols*multFactor, CV_8U);
        resize(sourceImage, zoomed, Size(sourceImage.cols*multFactor, sourceImage.rows*multFactor), sourceImage.cols*multFactor, sourceImage.rows*multFactor, INTER_NEAREST);
        if (normalized) { normalize(zoomed, zoomed, 0, 255, NORM_MINMAX); }
        namedWindow(name);
        imshow(name, zoomed);
    }
    
    Mat bruteForceLocalMax(Mat srcImage, int searchRad)
    {
        Mat outputArray = Mat::zeros(srcImage.rows, srcImage.cols, CV_8U);
        //global search top
        for (int x = 0; x < srcImage.rows - 1; x++)
        {
            for (int y = 0; y < srcImage.cols - 1; y++)
            {
                bool peak = true;
                float centerVal = srcImage.at<uchar>(x, y);
                if (centerVal == 0) { continue; }
                //local search top
                for (int a = -searchRad; a <= searchRad; a++)
                {
                    for (int b = -searchRad; b <= searchRad; b++)
                    {
                        if (x + a<0 || x + a>srcImage.rows - 1 || y + b < 0 || y + b>srcImage.cols - 1) { continue; }
                        if (srcImage.at<uchar>(x + a, y + b) > centerVal)
                        {
                            peak = false;
                        }
                        if (peak == false) { break; }
                    }
                    if (peak == false) { break; }
                }
                if (peak)
                {
                    outputArray.at<uchar>(x, y) = 255;
                }
    
            }
        }
        return outputArray;
    }
    
    Mat mergeLocalPeaks(Mat srcImage, int mergeRadius)
    {
        Mat outputArray = Mat::zeros(srcImage.rows, srcImage.cols, CV_8U);
        //global search top
        for (int x = 0; x < srcImage.rows - 1; x++)
        {
            for (int y = 0; y < srcImage.cols - 1; y++)
            {
                float centerVal = srcImage.at<uchar>(x, y);
                if (centerVal == 0) { continue; }
    
                int aveX = x;
                int aveY = y;
    
                int xCenter = -1;
                int yCenter = -1;
    
                while (aveX != xCenter || aveY != yCenter)
                {
                    xCenter = aveX;
                    yCenter = aveY;
                    aveX = 0;
                    aveY = 0;
                    int peakCount = 0;
                    //local search top
                    for (int a = -mergeRadius; a <= mergeRadius; a++)
                    {
                        for (int b = -mergeRadius; b <= mergeRadius; b++)
                        {
                            if (xCenter + a<0 || xCenter + a>srcImage.rows - 1 || yCenter + b < 0 || yCenter + b>srcImage.cols - 1) { continue; }
                            if (srcImage.at<uchar>(xCenter + a, yCenter + b) > 0)
                            {
                                aveX += (xCenter + a);
                                aveY += (yCenter + b);
                                peakCount += 1;
                            }
                        }
                    }
                    double dCentX = ((double)aveX / (double)peakCount);
                    double dCentY = ((double)aveY / (double)peakCount);
                    aveX = floor(dCentX);
                    aveY = floor(dCentY);
                }
                outputArray.at<uchar>(xCenter, yCenter) = 255;
            }
        }
        return outputArray;
    }
    

    speed: 在这里输入图像描述

    debug images: 在这里输入图像描述

    results: 在这里输入图像描述

    Hope this helps! Cheers!


    In this SO answer to a similar question I applied persistent homology to find peaks in an image. I took your image, scaled it down to 50%, applied an Guassian blur with radius 20 (in Gimp), and applied the methods described in the other article (click to enlarge):

    结果

    I only show peaks with persistence (see the other SO answer) of at least 20. The persistence diagram is shown here:

    The 20-th peak would be the little peak on the top-left corner with a persistence of about 9. By applying a stronger Gaussian filter the image gets more diffuse and the peaks will be come more prominent.

    Python code can be found here.

    链接地址: http://www.djcxy.com/p/84556.html

    上一篇: 算法从n返回k个元素的所有组合

    下一篇: 检测图像中的峰值