بایگانی برچسب برای: الگوریتم کنی

الگوریتم Canny در ++C

لبه یاب کنی توسط جان اف کنی در سال 1986 ایجاد شد و هنوز یک لبه یاب استاندارد و با دقت و کیفیت بالا میباشد.الگوریتم لبه یابی کنی یکی از بهترین لبه یابها تا به امروز است. در ادامه روش کار این الگوریتم و هم چنین کد الگوریتم Canny در ++C را بررسی خواهیم کرد. این الگوریتم لبه یابی از سه بخش اصلی زیر تشکیل شده است:

  • تضعیف نویز
  • پیدا کردن نقاطی که بتوان آنها را به عنوان لبه در نظر گرفت
  • حذب نقاطی که احتمال لبه بودن آنها کم است

 

معیارهایی که در لبه یاب کنی مطرح است:
1 -پایین آوردن نرخ خطا- یعنی تا حد امکان هیچ لبه ای در تصویر نباید گم شود و هم چنین هیچ چیزی که لبه نیست نباید به جای لبه فرض شود. لبه هان پیدا شده تا حد ممکن به لبه ها اصلی
نزدیک باشند.

2 -لبه در مکان واقعی خود باشد- یعنی تا حد ممکن لبه ها کمترین فاصله را با مکان واقعی خود داشته باشند.
3 -بران هر لبه فقط یک پاسخ داشته باشیم.

4 -لبه ها کمترین ضخامت را داشته باشند- (در صورت امکان یک پیکسل).
لبه یاب کنی بخاطر توانایی در تولید لبه های نازک تا حد یک ییکسل برای لبه های پیوسته معروف شده است. این لبه یاب شامل چهار مرحله و چهار ورودی زیر است:
یک تصویر ورودی
یک پارامتر به نام سیگما جهت مقدار نرم کنندگی تصویر
یک حد آستانه بالا (Th)
یک حد آستانه پایین (Tl)

 

مراحل الگوریتم Canny:

1- در ابتدا باید تصویر رنگی را به جهت لبه یابی بهتر به یک تصویر سطح خاکسترن تبدیب کرد.

2- نویز را از تصویر دریافتی حذف کرد. بدلیل اینکه فیلتر گاوسین از یک ماسک ساده برای حذف نویز استفاده می کند لبه یاب کنی در مرحله اول برای حذف نویز آن را بکار میگیرد.

3- در یک تصویر سطح خاکستر جایی را که بیشترین تغییرات را داشته باشند به عنوان لبه در نظر گرفته می شوند و این مکانها با گرفتن گرادیان تصویر با استفاده عملگر سوبل بدست می آیند. سپس لبه های مات یافت شده به لبه های تیزتر تبدیل می شوند.

4- برخی از لبه های کشف شده واقعا لبه نیستند و در واقع نویز هستند که باید آنها توسط حد آستانه هیسترزیس فیلتر شوند.هیسترزیس از دو حد آستانه بالاتر (Th) و حد آستانه پایین تر (Tl) استفاده کرده و کنی پیشنهاد می کند که نسبت استانه بالا به پایین سه به یک باشد.

 این روش بیشتر به کشف لبه های ضعیف به درستی می پردازد و کمتر فریب نویز را می خورد و از بقیه روش ها بهتر است.

 

الگوریتم Canny    عملکرد الگوریتم Canny

 

 

 

کد الگوریتم Canny در ++C:

برای الگوریتم Canny دو کد زیر ارائه می شود که کد شماره 2 کد کاملتری است.

کد شماره  الگوریتم 1 الگوریتم Canny در ++C:

در زیر استفاده از الگوریتم کنی در ++C است. توجه داشته باشید که تصویر ابتدا به تصویر سیاه و سفید تبدیل می شود، سپس فیلتر گاوسی برای کاهش نویز در تصویر استفاده می شود. سپس الگوریتم Canny برای تشخیص لبه استفاده می شود.

 

// CannyTutorial.cpp : Defines the entry point for the console application. 
// Environment: Visual studio 2015, Windows 10
// Assumptions: Opecv is installed configured in the visual studio project
// Opencv version: OpenCV 3.1

#include "stdafx.h"
#include<opencv2/highgui/highgui.hpp>
#include<opencv2/imgproc/imgproc.hpp>
#include<string>
#include<iostream>


int main()
{

    //Modified from source: https://github.com/MicrocontrollersAndMore/OpenCV_3_Windows_10_Installation_Tutorial
    cv::Mat imgOriginal;        // input image
    cv::Mat imgGrayscale;        // grayscale of input image
    cv::Mat imgBlurred;            // intermediate blured image
    cv::Mat imgCanny;            // Canny edge image

    std::cout << "Please enter an image filename : ";     std::string img_addr;     std::cin >> img_addr;

    std::cout << "Searching for " + img_addr << std::endl;

    imgOriginal = cv::imread(img_addr);            // open image

    if (imgOriginal.empty()) {                                    // if unable to open image
        std::cout << "error: image not read from file\n\n";        // show error message on command line
        return(0);                                                // and exit program
    }

    cv::cvtColor(imgOriginal, imgGrayscale, CV_BGR2GRAY);        // convert to grayscale

    cv::GaussianBlur(imgGrayscale,            // input image
        imgBlurred,                            // output image
        cv::Size(5, 5),                        // smoothing window width and height in pixels
        1.5);                                // sigma value, determines how much the image will be blurred

    cv::Canny(imgBlurred,            // input image
        imgCanny,                    // output image
        100,                        // low threshold
        200);                        // high threshold


    // Declare windows
    // Note: you can use CV_WINDOW_NORMAL which allows resizing the window
    // or CV_WINDOW_AUTOSIZE for a fixed size window matching the resolution of the image
    // CV_WINDOW_AUTOSIZE is the default
    cv::namedWindow("imgOriginal", CV_WINDOW_AUTOSIZE);        
    cv::namedWindow("imgCanny", CV_WINDOW_AUTOSIZE);

    //Show windows
    cv::imshow("imgOriginal", imgOriginal);        
    cv::imshow("imgCanny", imgCanny);

    cv::waitKey(0);                    // hold windows open until user presses a key
    return 0;
}

 

دانلود کد فوق از طریق لینک زیر:

رمز فایل : behsanandish.com

 

 

کد شماره 2:

مرحله 1: یک blur(تار کننده) گاوسی را اعمال کنید.

اول متغیرهای ضروری اعلام شده اند و بعضی از آنها اولیه هستند. سپس Blur گاوسی اعمال می شود. برای انجام این کار یک ماسک 5×5 بر روی تصویر منتقل می شود. هر پیکسل به صورت مجموع مقادیر پیکسل در محدوده 5×5 آن ضربدر وزن گاوسی متناظر تقسیم شده توسط وزن مجموع کل ماسک تعریف می شود.

 

ماسک گاوسی

ماسک گاوسی

 

#include "stdafx.h"
#include "tripod.h"
#include "tripodDlg.h"

#include "LVServerDefs.h"
#include "math.h"
#include <fstream>
#include <string>
#include <iostream>
#include <stdlib.h>
#include <stdio.h>


#ifdef _DEBUG
#define new DEBUG_NEW
#undef THIS_FILE
static char THIS_FILE[] = __FILE__;
#endif

using namespace std;

/////////////////////////////////////////////////////////////////////////////
// CAboutDlg dialog used for App About

class CAboutDlg : public CDialog
{
public:
	CAboutDlg();

// Dialog Data
	//{{AFX_DATA(CAboutDlg)
	enum { IDD = IDD_ABOUTBOX };
	//}}AFX_DATA

	// ClassWizard generated virtual function overrides
	//{{AFX_VIRTUAL(CAboutDlg)
	protected:
	virtual void DoDataExchange(CDataExchange* pDX);    // DDX/DDV support
	//}}AFX_VIRTUAL

// Implementation
protected:
	//{{AFX_MSG(CAboutDlg)
	//}}AFX_MSG
	DECLARE_MESSAGE_MAP()
};

CAboutDlg::CAboutDlg() : CDialog(CAboutDlg::IDD)
{
	//{{AFX_DATA_INIT(CAboutDlg)
	//}}AFX_DATA_INIT
}

void CAboutDlg::DoDataExchange(CDataExchange* pDX)
{
	CDialog::DoDataExchange(pDX);
	//{{AFX_DATA_MAP(CAboutDlg)
	//}}AFX_DATA_MAP
}

BEGIN_MESSAGE_MAP(CAboutDlg, CDialog)
	//{{AFX_MSG_MAP(CAboutDlg)
		// No message handlers
	//}}AFX_MSG_MAP
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CTripodDlg dialog

CTripodDlg::CTripodDlg(CWnd* pParent /*=NULL*/)
	: CDialog(CTripodDlg::IDD, pParent)
{
	//{{AFX_DATA_INIT(CTripodDlg)
		// NOTE: the ClassWizard will add member initialization here
	//}}AFX_DATA_INIT
	// Note that LoadIcon does not require a subsequent DestroyIcon in Win32
	m_hIcon = AfxGetApp()->LoadIcon(IDR_MAINFRAME);

	//////////////// Set destination BMP to NULL first 
	m_destinationBitmapInfoHeader = NULL;

}

////////////////////// Additional generic functions

static unsigned PixelBytes(int w, int bpp)
{
    return (w * bpp + 7) / 8;
}

static unsigned DibRowSize(int w, int bpp)
{
    return (w * bpp + 31) / 32 * 4;
}

static unsigned DibRowSize(LPBITMAPINFOHEADER pbi)
{
    return DibRowSize(pbi->biWidth, pbi->biBitCount);
}

static unsigned DibRowPadding(int w, int bpp)
{
    return DibRowSize(w, bpp) - PixelBytes(w, bpp);
}

static unsigned DibRowPadding(LPBITMAPINFOHEADER pbi)
{
    return DibRowPadding(pbi->biWidth, pbi->biBitCount);
}

static unsigned DibImageSize(int w, int h, int bpp)
{
    return h * DibRowSize(w, bpp);
}

static size_t DibSize(int w, int h, int bpp)
{
    return sizeof (BITMAPINFOHEADER) + DibImageSize(w, h, bpp);
}

/////////////////////// end of generic functions


void CTripodDlg::DoDataExchange(CDataExchange* pDX)
{
	CDialog::DoDataExchange(pDX);
	//{{AFX_DATA_MAP(CTripodDlg)
	DDX_Control(pDX, IDC_PROCESSEDVIEW, m_cVideoProcessedView);
	DDX_Control(pDX, IDC_UNPROCESSEDVIEW, m_cVideoUnprocessedView);
	//}}AFX_DATA_MAP
}

BEGIN_MESSAGE_MAP(CTripodDlg, CDialog)
	//{{AFX_MSG_MAP(CTripodDlg)
	ON_WM_SYSCOMMAND()
	ON_WM_PAINT()
	ON_WM_QUERYDRAGICON()
	ON_BN_CLICKED(IDEXIT, OnExit)
	//}}AFX_MSG_MAP
END_MESSAGE_MAP()

/////////////////////////////////////////////////////////////////////////////
// CTripodDlg message handlers

BOOL CTripodDlg::OnInitDialog()
{
	CDialog::OnInitDialog();

	// Add "About..." menu item to system menu.

	// IDM_ABOUTBOX must be in the system command range.
	ASSERT((IDM_ABOUTBOX & 0xFFF0) == IDM_ABOUTBOX);
	ASSERT(IDM_ABOUTBOX < 0xF000); CMenu* pSysMenu = GetSystemMenu(FALSE); if (pSysMenu != NULL) { CString strAboutMenu; strAboutMenu.LoadString(IDS_ABOUTBOX); if (!strAboutMenu.IsEmpty()) { pSysMenu->AppendMenu(MF_SEPARATOR);
			pSysMenu->AppendMenu(MF_STRING, IDM_ABOUTBOX, strAboutMenu);
		}
	}

	// Set the icon for this dialog.  The framework does this automatically
	//  when the application's main window is not a dialog
	SetIcon(m_hIcon, TRUE);			// Set big icon
	SetIcon(m_hIcon, FALSE);		// Set small icon
	
	// TODO: Add extra initialization here

	// For Unprocessed view videoportal (top one)
	char sRegUnprocessedView[] = "HKEY_LOCAL_MACHINE\\Software\\UnprocessedView";
	m_cVideoUnprocessedView.PrepareControl("UnprocessedView", sRegUnprocessedView, 0 );	
	m_cVideoUnprocessedView.EnableUIElements(UIELEMENT_STATUSBAR,0,TRUE);
	m_cVideoUnprocessedView.ConnectCamera2();
	m_cVideoUnprocessedView.SetEnablePreview(TRUE);

	// For binary view videoportal (bottom one)
	char sRegProcessedView[] = "HKEY_LOCAL_MACHINE\\Software\\ProcessedView";
	m_cVideoProcessedView.PrepareControl("ProcessedView", sRegProcessedView, 0 );	
	m_cVideoProcessedView.EnableUIElements(UIELEMENT_STATUSBAR,0,TRUE);
	m_cVideoProcessedView.ConnectCamera2();
	m_cVideoProcessedView.SetEnablePreview(TRUE);

	// Initialize the size of binary videoportal
	m_cVideoProcessedView.SetPreviewMaxHeight(240);
	m_cVideoProcessedView.SetPreviewMaxWidth(320);

	// Uncomment if you wish to fix the live videoportal's size
	// m_cVideoUnprocessedView.SetPreviewMaxHeight(240);
	// m_cVideoUnprocessedView.SetPreviewMaxWidth(320);

	// Find the screen coodinates of the binary videoportal
	m_cVideoProcessedView.GetWindowRect(m_rectForProcessedView);
	ScreenToClient(m_rectForProcessedView);
	allocateDib(CSize(320, 240));

	// Start grabbing frame data for Procssed videoportal (bottom one)
	m_cVideoProcessedView.StartVideoHook(0);

	return TRUE;  // return TRUE  unless you set the focus to a control
}

void CTripodDlg::OnSysCommand(UINT nID, LPARAM lParam)
{
	if ((nID & 0xFFF0) == IDM_ABOUTBOX)
	{
		CAboutDlg dlgAbout;
		dlgAbout.DoModal();
	}
	else
	{
		CDialog::OnSysCommand(nID, lParam);
	}
}

// If you add a minimize button to your dialog, you will need the code below
//  to draw the icon.  For MFC applications using the document/view model,
//  this is automatically done for you by the framework.

void CTripodDlg::OnPaint() 
{
	if (IsIconic())
	{
		CPaintDC dc(this); // device context for painting

		SendMessage(WM_ICONERASEBKGND, (WPARAM) dc.GetSafeHdc(), 0);

		// Center icon in client rectangle
		int cxIcon = GetSystemMetrics(SM_CXICON);
		int cyIcon = GetSystemMetrics(SM_CYICON);
		CRect rect;
		GetClientRect(&rect);
		int x = (rect.Width() - cxIcon + 1) / 2;
		int y = (rect.Height() - cyIcon + 1) / 2;

		// Draw the icon
		dc.DrawIcon(x, y, m_hIcon);
	}
	else
	{
		CDialog::OnPaint();
	}
}

// The system calls this to obtain the cursor to display while the user drags
//  the minimized window.
HCURSOR CTripodDlg::OnQueryDragIcon()
{
	return (HCURSOR) m_hIcon;
}

void CTripodDlg::OnExit() 
{
	// TODO: Add your control notification handler code here

	// Kill live view videoportal (top one)
	m_cVideoUnprocessedView.StopVideoHook(0);
    m_cVideoUnprocessedView.DisconnectCamera();	
	
	// Kill binary view videoportal (bottom one)
	m_cVideoProcessedView.StopVideoHook(0);
    m_cVideoProcessedView.DisconnectCamera();	

	// Kill program
	DestroyWindow();	

	

}

BEGIN_EVENTSINK_MAP(CTripodDlg, CDialog)
    //{{AFX_EVENTSINK_MAP(CTripodDlg)
	ON_EVENT(CTripodDlg, IDC_PROCESSEDVIEW, 1 /* PortalNotification */, OnPortalNotificationProcessedview, VTS_I4 VTS_I4 VTS_I4 VTS_I4)
	//}}AFX_EVENTSINK_MAP
END_EVENTSINK_MAP()

void CTripodDlg::OnPortalNotificationProcessedview(long lMsg, long lParam1, long lParam2, long lParam3) 
{
	// TODO: Add your control notification handler code here
	
	// This function is called at the camera's frame rate
    
#define NOTIFICATIONMSG_VIDEOHOOK	10

	// Declare some useful variables
	// QCSDKMFC.pdf (Quickcam MFC documentation) p. 103 explains the variables lParam1, lParam2, lParam3 too 
	
	LPBITMAPINFOHEADER lpBitmapInfoHeader; // Frame's info header contains info like width and height
	LPBYTE lpBitmapPixelData; // This pointer-to-long will point to the start of the frame's pixel data
    unsigned long lTimeStamp; // Time when frame was grabbed

	switch(lMsg) {
		case NOTIFICATIONMSG_VIDEOHOOK:
			{
				lpBitmapInfoHeader = (LPBITMAPINFOHEADER) lParam1; 
				lpBitmapPixelData = (LPBYTE) lParam2;
				lTimeStamp = (unsigned long) lParam3;

				grayScaleTheFrameData(lpBitmapInfoHeader, lpBitmapPixelData);
				doMyImageProcessing(lpBitmapInfoHeader); // Place where you'd add your image processing code
				displayMyResults(lpBitmapInfoHeader);

			}
			break;

		default:
			break;
	}	
}

void CTripodDlg::allocateDib(CSize sz)
{
	// Purpose: allocate information for a device independent bitmap (DIB)
	// Called from OnInitVideo

	if(m_destinationBitmapInfoHeader) {
		free(m_destinationBitmapInfoHeader);
		m_destinationBitmapInfoHeader = NULL;
	}

	if(sz.cx | sz.cy) {
		m_destinationBitmapInfoHeader = (LPBITMAPINFOHEADER)malloc(DibSize(sz.cx, sz.cy, 24));
		ASSERT(m_destinationBitmapInfoHeader);
		m_destinationBitmapInfoHeader->biSize = sizeof(BITMAPINFOHEADER);
		m_destinationBitmapInfoHeader->biWidth = sz.cx;
		m_destinationBitmapInfoHeader->biHeight = sz.cy;
		m_destinationBitmapInfoHeader->biPlanes = 1;
		m_destinationBitmapInfoHeader->biBitCount = 24;
		m_destinationBitmapInfoHeader->biCompression = 0;
		m_destinationBitmapInfoHeader->biSizeImage = DibImageSize(sz.cx, sz.cy, 24);
		m_destinationBitmapInfoHeader->biXPelsPerMeter = 0;
		m_destinationBitmapInfoHeader->biYPelsPerMeter = 0;
		m_destinationBitmapInfoHeader->biClrImportant = 0;
		m_destinationBitmapInfoHeader->biClrUsed = 0;
	}
}

void CTripodDlg::displayMyResults(LPBITMAPINFOHEADER lpThisBitmapInfoHeader)
{
	// displayMyResults: Displays results of doMyImageProcessing() in the videoport
	// Notes: StretchDIBits stretches a device-independent bitmap to the appropriate size

	CDC				*pDC;	// Device context to display bitmap data
	
	pDC = GetDC();	
	int nOldMode = SetStretchBltMode(pDC->GetSafeHdc(),COLORONCOLOR);

	StretchDIBits( 
		pDC->GetSafeHdc(),
		m_rectForProcessedView.left,				// videoportal left-most coordinate
		m_rectForProcessedView.top,					// videoportal top-most coordinate
		m_rectForProcessedView.Width(),				// videoportal width
		m_rectForProcessedView.Height(),			// videoportal height
		0,											// Row position to display bitmap in videoportal
		0,											// Col position to display bitmap in videoportal
		lpThisBitmapInfoHeader->biWidth,			// m_destinationBmp's number of columns
		lpThisBitmapInfoHeader->biHeight,			// m_destinationBmp's number of rows
		m_destinationBmp,							// The bitmap to display; use the one resulting from doMyImageProcessing
		(BITMAPINFO*)m_destinationBitmapInfoHeader, // The bitmap's header info e.g. width, height, number of bits etc
		DIB_RGB_COLORS,								// Use default 24-bit color table
		SRCCOPY										// Just display
	);
 
	SetStretchBltMode(pDC->GetSafeHdc(),nOldMode);

	ReleaseDC(pDC);

	// Note: 04/24/02 - Added the following:
	// Christopher Wagner cwagner@fas.harvard.edu noticed that memory wasn't being freed

	// Recall OnPortalNotificationProcessedview, which gets called everytime
	// a frame of data arrives, performs 3 steps:
	// (1) grayScaleTheFrameData - which mallocs m_destinationBmp
	// (2) doMyImageProcesing
	// (3) displayMyResults - which we're in now
	// Since we're finished with the memory we malloc'ed for m_destinationBmp
	// we should free it: 
	
	free(m_destinationBmp);

	// End of adds
}

void CTripodDlg::grayScaleTheFrameData(LPBITMAPINFOHEADER lpThisBitmapInfoHeader, LPBYTE lpThisBitmapPixelData)
{

	// grayScaleTheFrameData: Called by CTripodDlg::OnPortalNotificationBinaryview
	// Task: Read current frame pixel data and computes a grayscale version

	unsigned int	W, H;			  // Width and Height of current frame [pixels]
	BYTE            *sourceBmp;		  // Pointer to current frame of data
	unsigned int    row, col;
	unsigned long   i;
	BYTE			grayValue;

	BYTE			redValue;
	BYTE			greenValue;
	BYTE			blueValue;

    W = lpThisBitmapInfoHeader->biWidth;  // biWidth: number of columns
    H = lpThisBitmapInfoHeader->biHeight; // biHeight: number of rows

	// Store pixel data in row-column vector format
	// Recall that each pixel requires 3 bytes (red, blue and green bytes)
	// m_destinationBmp is a protected member and declared in binarizeDlg.h

	m_destinationBmp = (BYTE*)malloc(H*3*W*sizeof(BYTE));

	// Point to the current frame's pixel data
	sourceBmp = lpThisBitmapPixelData;

	for (row = 0; row < H; row++) {
		for (col = 0; col < W; col++) { // Recall each pixel is composed of 3 bytes i = (unsigned long)(row*3*W + 3*col); // The source pixel has a blue, green andred value: blueValue = *(sourceBmp + i); greenValue = *(sourceBmp + i + 1); redValue = *(sourceBmp + i + 2); // A standard equation for computing a grayscale value based on RGB values grayValue = (BYTE)(0.299*redValue + 0.587*greenValue + 0.114*blueValue); // The destination BMP will be a grayscale version of the source BMP *(m_destinationBmp + i) = grayValue; *(m_destinationBmp + i + 1) = grayValue; *(m_destinationBmp + i + 2) = grayValue; } } } void CTripodDlg::doMyImageProcessing(LPBITMAPINFOHEADER lpThisBitmapInfoHeader) { // doMyImageProcessing: This is where you'd write your own image processing code // Task: Read a pixel's grayscale value and process accordingly unsigned int W, H; // Width and Height of current frame [pixels] unsigned int row, col; // Pixel's row and col positions unsigned long i; // Dummy variable for row-column vector int upperThreshold = 60; // Gradient strength nessicary to start edge int lowerThreshold = 30; // Minimum gradient strength to continue edge unsigned long iOffset; // Variable to offset row-column vector during sobel mask int rowOffset; // Row offset from the current pixel int colOffset; // Col offset from the current pixel int rowTotal = 0; // Row position of offset pixel int colTotal = 0; // Col position of offset pixel int Gx; // Sum of Sobel mask products values in the x direction int Gy; // Sum of Sobel mask products values in the y direction float thisAngle; // Gradient direction based on Gx and Gy int newAngle; // Approximation of the gradient direction bool edgeEnd; // Stores whether or not the edge is at the edge of the possible image int GxMask[3][3]; // Sobel mask in the x direction int GyMask[3][3]; // Sobel mask in the y direction int newPixel; // Sum pixel values for gaussian int gaussianMask[5][5]; // Gaussian mask W = lpThisBitmapInfoHeader->biWidth;  // biWidth: number of columns
    H = lpThisBitmapInfoHeader->biHeight; // biHeight: number of rows
	
	for (row = 0; row < H; row++) {
		for (col = 0; col < W; col++) {
			edgeDir[row][col] = 0;
		}
	}

	/* Declare Sobel masks */
	GxMask[0][0] = -1; GxMask[0][1] = 0; GxMask[0][2] = 1;
	GxMask[1][0] = -2; GxMask[1][1] = 0; GxMask[1][2] = 2;
	GxMask[2][0] = -1; GxMask[2][1] = 0; GxMask[2][2] = 1;
	
	GyMask[0][0] =  1; GyMask[0][1] =  2; GyMask[0][2] =  1;
	GyMask[1][0] =  0; GyMask[1][1] =  0; GyMask[1][2] =  0;
	GyMask[2][0] = -1; GyMask[2][1] = -2; GyMask[2][2] = -1;

	/* Declare Gaussian mask */
	gaussianMask[0][0] = 2;		gaussianMask[0][1] = 4;		gaussianMask[0][2] = 5;		gaussianMask[0][3] = 4;		gaussianMask[0][4] = 2;	
	gaussianMask[1][0] = 4;		gaussianMask[1][1] = 9;		gaussianMask[1][2] = 12;	gaussianMask[1][3] = 9;		gaussianMask[1][4] = 4;	
	gaussianMask[2][0] = 5;		gaussianMask[2][1] = 12;	gaussianMask[2][2] = 15;	gaussianMask[2][3] = 12;	gaussianMask[2][4] = 2;	
	gaussianMask[3][0] = 4;		gaussianMask[3][1] = 9;		gaussianMask[3][2] = 12;	gaussianMask[3][3] = 9;		gaussianMask[3][4] = 4;	
	gaussianMask[4][0] = 2;		gaussianMask[4][1] = 4;		gaussianMask[4][2] = 5;		gaussianMask[4][3] = 4;		gaussianMask[4][4] = 2;	
	

	/* Gaussian Blur */
	for (row = 2; row < H-2; row++) {
		for (col = 2; col < W-2; col++) {
			newPixel = 0;
			for (rowOffset=-2; rowOffset<=2; rowOffset++) {
				for (colOffset=-2; colOffset<=2; colOffset++) {
					rowTotal = row + rowOffset;
					colTotal = col + colOffset;
					iOffset = (unsigned long)(rowTotal*3*W + colTotal*3);
					newPixel += (*(m_destinationBmp + iOffset)) * gaussianMask[2 + rowOffset][2 + colOffset];
				}
			}
			i = (unsigned long)(row*3*W + col*3);
			*(m_destinationBmp + i) = newPixel / 159;
		}
	}

	/* Determine edge directions and gradient strengths */
	for (row = 1; row < H-1; row++) {
		for (col = 1; col < W-1; col++) {
			i = (unsigned long)(row*3*W + 3*col);
			Gx = 0;
			Gy = 0;
			/* Calculate the sum of the Sobel mask times the nine surrounding pixels in the x and y direction */
			for (rowOffset=-1; rowOffset<=1; rowOffset++) {
				for (colOffset=-1; colOffset<=1; colOffset++) {
					rowTotal = row + rowOffset;
					colTotal = col + colOffset;
					iOffset = (unsigned long)(rowTotal*3*W + colTotal*3);
					Gx = Gx + (*(m_destinationBmp + iOffset) * GxMask[rowOffset + 1][colOffset + 1]);
					Gy = Gy + (*(m_destinationBmp + iOffset) * GyMask[rowOffset + 1][colOffset + 1]);
				}
			}

			gradient[row][col] = sqrt(pow(Gx,2.0) + pow(Gy,2.0));	// Calculate gradient strength			
			thisAngle = (atan2(Gx,Gy)/3.14159) * 180.0;		// Calculate actual direction of edge
			
			/* Convert actual edge direction to approximate value */
			if ( ( (thisAngle < 22.5) && (thisAngle > -22.5) ) || (thisAngle > 157.5) || (thisAngle < -157.5) ) newAngle = 0; if ( ( (thisAngle > 22.5) && (thisAngle < 67.5) ) || ( (thisAngle < -112.5) && (thisAngle > -157.5) ) )
				newAngle = 45;
			if ( ( (thisAngle > 67.5) && (thisAngle < 112.5) ) || ( (thisAngle < -67.5) && (thisAngle > -112.5) ) )
				newAngle = 90;
			if ( ( (thisAngle > 112.5) && (thisAngle < 157.5) ) || ( (thisAngle < -22.5) && (thisAngle > -67.5) ) )
				newAngle = 135;
				
			edgeDir[row][col] = newAngle;		// Store the approximate edge direction of each pixel in one array
		}
	}

	/* Trace along all the edges in the image */
	for (row = 1; row < H - 1; row++) {
		for (col = 1; col < W - 1; col++) { edgeEnd = false; if (gradient[row][col] > upperThreshold) {		// Check to see if current pixel has a high enough gradient strength to be part of an edge
				/* Switch based on current pixel's edge direction */
				switch (edgeDir[row][col]){		
					case 0:
						findEdge(0, 1, row, col, 0, lowerThreshold);
						break;
					case 45:
						findEdge(1, 1, row, col, 45, lowerThreshold);
						break;
					case 90:
						findEdge(1, 0, row, col, 90, lowerThreshold);
						break;
					case 135:
						findEdge(1, -1, row, col, 135, lowerThreshold);
						break;
					default :
						i = (unsigned long)(row*3*W + 3*col);
						*(m_destinationBmp + i) = 
						*(m_destinationBmp + i + 1) = 
						*(m_destinationBmp + i + 2) = 0;
						break;
					}
				}
			else {
				i = (unsigned long)(row*3*W + 3*col);
					*(m_destinationBmp + i) = 
					*(m_destinationBmp + i + 1) = 
					*(m_destinationBmp + i + 2) = 0;
			}	
		}
	}
	
	/* Suppress any pixels not changed by the edge tracing */
	for (row = 0; row < H; row++) {
		for (col = 0; col < W; col++) {	
			// Recall each pixel is composed of 3 bytes
			i = (unsigned long)(row*3*W + 3*col);
			// If a pixel's grayValue is not black or white make it black
			if( ((*(m_destinationBmp + i) != 255) && (*(m_destinationBmp + i) != 0)) || ((*(m_destinationBmp + i + 1) != 255) && (*(m_destinationBmp + i + 1) != 0)) || ((*(m_destinationBmp + i + 2) != 255) && (*(m_destinationBmp + i + 2) != 0)) ) 
				*(m_destinationBmp + i) = 
				*(m_destinationBmp + i + 1) = 
				*(m_destinationBmp + i + 2) = 0; // Make pixel black
		}
	}

	/* Non-maximum Suppression */
	for (row = 1; row < H - 1; row++) {
		for (col = 1; col < W - 1; col++) {
			i = (unsigned long)(row*3*W + 3*col);
			if (*(m_destinationBmp + i) == 255) {		// Check to see if current pixel is an edge
				/* Switch based on current pixel's edge direction */
				switch (edgeDir[row][col]) {		
					case 0:
						suppressNonMax( 1, 0, row, col, 0, lowerThreshold);
						break;
					case 45:
						suppressNonMax( 1, -1, row, col, 45, lowerThreshold);
						break;
					case 90:
						suppressNonMax( 0, 1, row, col, 90, lowerThreshold);
						break;
					case 135:
						suppressNonMax( 1, 1, row, col, 135, lowerThreshold);
						break;
					default :
						break;
				}
			}	
		}
	}
	
}

void CTripodDlg::findEdge(int rowShift, int colShift, int row, int col, int dir, int lowerThreshold)
{
	int W = 320;
	int H = 240;
	int newRow;
	int newCol;
	unsigned long i;
	bool edgeEnd = false;

	/* Find the row and column values for the next possible pixel on the edge */
	if (colShift < 0) { if (col > 0)
			newCol = col + colShift;
		else
			edgeEnd = true;
	} else if (col < W - 1) {
		newCol = col + colShift;
	} else
		edgeEnd = true;		// If the next pixel would be off image, don't do the while loop
	if (rowShift < 0) { if (row > 0)
			newRow = row + rowShift;
		else
			edgeEnd = true;
	} else if (row < H - 1) { newRow = row + rowShift; } else edgeEnd = true; /* Determine edge directions and gradient strengths */ while ( (edgeDir[newRow][newCol]==dir) && !edgeEnd && (gradient[newRow][newCol] > lowerThreshold) ) {
		/* Set the new pixel as white to show it is an edge */
		i = (unsigned long)(newRow*3*W + 3*newCol);
		*(m_destinationBmp + i) =
		*(m_destinationBmp + i + 1) =
		*(m_destinationBmp + i + 2) = 255;
		if (colShift < 0) { if (newCol > 0)
				newCol = newCol + colShift;
			else
				edgeEnd = true;	
		} else if (newCol < W - 1) {
			newCol = newCol + colShift;
		} else
			edgeEnd = true;	
		if (rowShift < 0) { if (newRow > 0)
				newRow = newRow + rowShift;
			else
				edgeEnd = true;
		} else if (newRow < H - 1) {
			newRow = newRow + rowShift;
		} else
			edgeEnd = true;	
	}	
}

void CTripodDlg::suppressNonMax(int rowShift, int colShift, int row, int col, int dir, int lowerThreshold)
{
	int W = 320;
	int H = 240;
	int newRow = 0;
	int newCol = 0;
	unsigned long i;
	bool edgeEnd = false;
	float nonMax[320][3];			// Temporarily stores gradients and positions of pixels in parallel edges
	int pixelCount = 0;					// Stores the number of pixels in parallel edges
	int count;						// A for loop counter
	int max[3];						// Maximum point in a wide edge
	
	if (colShift < 0) { if (col > 0)
			newCol = col + colShift;
		else
			edgeEnd = true;
	} else if (col < W - 1) {
		newCol = col + colShift;
	} else
		edgeEnd = true;		// If the next pixel would be off image, don't do the while loop
	if (rowShift < 0) { if (row > 0)
			newRow = row + rowShift;
		else
			edgeEnd = true;
	} else if (row < H - 1) {
		newRow = row + rowShift;
	} else
		edgeEnd = true;	
	i = (unsigned long)(newRow*3*W + 3*newCol);
	/* Find non-maximum parallel edges tracing up */
	while ((edgeDir[newRow][newCol] == dir) && !edgeEnd && (*(m_destinationBmp + i) == 255)) {
		if (colShift < 0) { if (newCol > 0)
				newCol = newCol + colShift;
			else
				edgeEnd = true;	
		} else if (newCol < W - 1) {
			newCol = newCol + colShift;
		} else
			edgeEnd = true;	
		if (rowShift < 0) { if (newRow > 0)
				newRow = newRow + rowShift;
			else
				edgeEnd = true;
		} else if (newRow < H - 1) {
			newRow = newRow + rowShift;
		} else
			edgeEnd = true;	
		nonMax[pixelCount][0] = newRow;
		nonMax[pixelCount][1] = newCol;
		nonMax[pixelCount][2] = gradient[newRow][newCol];
		pixelCount++;
		i = (unsigned long)(newRow*3*W + 3*newCol);
	}

	/* Find non-maximum parallel edges tracing down */
	edgeEnd = false;
	colShift *= -1;
	rowShift *= -1;
	if (colShift < 0) { if (col > 0)
			newCol = col + colShift;
		else
			edgeEnd = true;
	} else if (col < W - 1) {
		newCol = col + colShift;
	} else
		edgeEnd = true;	
	if (rowShift < 0) { if (row > 0)
			newRow = row + rowShift;
		else
			edgeEnd = true;
	} else if (row < H - 1) {
		newRow = row + rowShift;
	} else
		edgeEnd = true;	
	i = (unsigned long)(newRow*3*W + 3*newCol);
	while ((edgeDir[newRow][newCol] == dir) && !edgeEnd && (*(m_destinationBmp + i) == 255)) {
		if (colShift < 0) { if (newCol > 0)
				newCol = newCol + colShift;
			else
				edgeEnd = true;	
		} else if (newCol < W - 1) {
			newCol = newCol + colShift;
		} else
			edgeEnd = true;	
		if (rowShift < 0) { if (newRow > 0)
				newRow = newRow + rowShift;
			else
				edgeEnd = true;
		} else if (newRow < H - 1) {
			newRow = newRow + rowShift;
		} else
			edgeEnd = true;	
		nonMax[pixelCount][0] = newRow;
		nonMax[pixelCount][1] = newCol;
		nonMax[pixelCount][2] = gradient[newRow][newCol];
		pixelCount++;
		i = (unsigned long)(newRow*3*W + 3*newCol);
	}

	/* Suppress non-maximum edges */
	max[0] = 0;
	max[1] = 0;
	max[2] = 0;
	for (count = 0; count < pixelCount; count++) { if (nonMax[count][2] > max[2]) {
			max[0] = nonMax[count][0];
			max[1] = nonMax[count][1];
			max[2] = nonMax[count][2];
		}
	}
	for (count = 0; count < pixelCount; count++) {
		i = (unsigned long)(nonMax[count][0]*3*W + 3*nonMax[count][1]);
		*(m_destinationBmp + i) = 
		*(m_destinationBmp + i + 1) = 
		*(m_destinationBmp + i + 2) = 0;
	}
}

الگوریتم Canny در سی پلاس پلاس قسمت 1
الگوریتم Canny در سی پلاس پلاس قسمت 2
الگوریتم Canny در سی پلاس پلاس قسمت 3
الگوریتم Canny در سی پلاس پلاس فسمت 4

الگوریتم Canny

لبه یاب کنی توسط جان اف کنی در سال 1986 ایجاد شد و هنوز یک لبه یاب استاندارد و با دقت و کیفیت بالا میباشد.الگوریتم لبه یابی کنی یکی از بهترین لبه یابها تا به امروز است. در ادامه روش کار این الگوریتم و هم چنین کد الگوریتم Canny در C را بررسی خواهیم کرد. این الگوریتم لبه یابی از سه بخش اصلی زیر تشکیل شده است:

  • تضعیف نویز
  • پیدا کردن نقاطی که بتوان آنها را به عنوان لبه در نظر گرفت
  • جذب نقاطی که احتمال لبه بودن آنها کم است

 

معیارهایی که در لبه یاب کنی مطرح می باشد:
1 -پایین آوردن نرخ خطا- یعنی تا حد امکان هیچ لبه ای در تصویر نباید گم شود و هم چنین هیچ چیزی که لبه نیست نباید به جای لبه فرض شود. لبه هان پیدا شده تا حد ممکن به لبه ها اصلی
نزدیک باشند.

2 -لبه در مکان واقعی خود باشد- یعنی تا حد ممکن لبه ها کمترین فاصله را با مکان واقعی خود داشته باشند.
3 -بران هر لبه فقط یک پاسخ داشته باشیم.

4 -لبه ها کمترین ضخامت را داشته باشند- (در صورت امکان یک پیکسل).
لبه یاب کنی بخاطر توانایی در تولید لبه های نازک تا حد یک ییکسل برای لبه های پیوسته معروف شده است. این لبه یاب شامل چهار مرحله و چهار ورودی زیر است:
یک تصویر ورودی
یک پارامتر به نام سیگما جهت مقدار نرم کنندگی تصویر
یک حد آستانه بالا (Th)
یک حد آستانه پایین (Tl)

 

مراحل الگوریتم Canny:

1- در ابتدا باید تصویر رنگی را به جهت لبه یابی بهتر به یک تصویر سطح خاکسترن تبدیب کرد.

2- نویز را از تصویر دریافتی حذف کرد. بدلیل اینکه فیلتر گاوسین از یک ماسک ساده برای حذف نویز استفاده می کند لبه یاب کنی در مرحله اول برای حذف نویز آن را بکار می گیرد.

3- در یک تصویر سطح خاکستر جایی را که بیشترین تغییرات را داشته باشند به عنوان لبه در نظر گرفته می شوند و این مکانها با گرفتن گرادیان تصویر با استفاده عملگر سوبل بدست می آیند. سپس لبه های مات یافت شده به لبه های تیزتر تبدیل می شوند.

4- برخی از لبه های کشف شده واقعا لبه نیستند و در واقع نویز هستند که باید آنها توسط حد آستانه هیسترزیس فیلتر شوند.هیسترزیس از دو حد آستانه بالاتر (Th) و حد آستانه پایین تر (Tl) استفاده کرده و کنی پیشنهاد می کند که نسبت استانه بالا به پایین سه به یک باشد.

 این روش بیشتر به کشف لبه های ضعیف به درستی می پردازد و کمتر فریب نویز را می خورد و از بقیه روش ها بهتر است.

 

الگوریتم Canny    عملکرد الگوریتم Canny

 

 

کد الگوریتم Canny در C :

برنامه زیر یک فایل BMP سیاه و سفید 8 بیت در هر پیکسل را می خواند و نتیجه را در ‘out.bmp’ ذخیره می کند.با `-lm ‘ کامپایل می شود.

 

#include <stdint.h>
#include <stdio.h>
#include <stdlib.h>
#include <float.h>
#include <math.h>
#include <string.h>
#include <stdbool.h>
#include <assert.h>
 
#define MAX_BRIGHTNESS 255
 
// C99 doesn't define M_PI (GNU-C99 does)
#define M_PI 3.14159265358979323846264338327
 
/*
 * Loading part taken from
 * http://www.vbforums.com/showthread.php?t=261522
 * BMP info:
 * http://en.wikipedia.org/wiki/BMP_file_format
 *
 * Note: the magic number has been removed from the bmpfile_header_t
 * structure since it causes alignment problems
 * bmpfile_magic_t should be written/read first
 * followed by the
 * bmpfile_header_t
 * [this avoids compiler-specific alignment pragmas etc.]
 */
 
typedef struct {
 uint8_t magic[2];
} bmpfile_magic_t;
 
typedef struct {
 uint32_t filesz;
 uint16_t creator1;
 uint16_t creator2;
 uint32_t bmp_offset;
} bmpfile_header_t;
 
typedef struct {
 uint32_t header_sz;
 int32_t width;
 int32_t height;
 uint16_t nplanes;
 uint16_t bitspp;
 uint32_t compress_type;
 uint32_t bmp_bytesz;
 int32_t hres;
 int32_t vres;
 uint32_t ncolors;
 uint32_t nimpcolors;
} bitmap_info_header_t;
 
typedef struct {
 uint8_t r;
 uint8_t g;
 uint8_t b;
 uint8_t nothing;
} rgb_t;
 
// Use short int instead `unsigned char' so that we can
// store negative values.
typedef short int pixel_t;
 
pixel_t *load_bmp(const char *filename,
 bitmap_info_header_t *bitmapInfoHeader)
{
 FILE *filePtr = fopen(filename, "rb");
 if (filePtr == NULL) {
 perror("fopen()");
 return NULL;
 }
 
 bmpfile_magic_t mag;
 if (fread(&mag, sizeof(bmpfile_magic_t), 1, filePtr) != 1) {
 fclose(filePtr);
 return NULL;
 }
 
 // verify that this is a bmp file by check bitmap id
 // warning: dereferencing type-punned pointer will break
 // strict-aliasing rules [-Wstrict-aliasing]
 if (*((uint16_t*)mag.magic) != 0x4D42) {
 fprintf(stderr, "Not a BMP file: magic=%c%c\n",
 mag.magic[0], mag.magic[1]);
 fclose(filePtr);
 return NULL;
 }
 
 bmpfile_header_t bitmapFileHeader; // our bitmap file header
 // read the bitmap file header
 if (fread(&bitmapFileHeader, sizeof(bmpfile_header_t),
 1, filePtr) != 1) {
 fclose(filePtr);
 return NULL;
 }
 
 // read the bitmap info header
 if (fread(bitmapInfoHeader, sizeof(bitmap_info_header_t),
 1, filePtr) != 1) {
 fclose(filePtr);
 return NULL;
 }
 
 if (bitmapInfoHeader->compress_type != 0)
 fprintf(stderr, "Warning, compression is not supported.\n");
 
 // move file point to the beginning of bitmap data
 if (fseek(filePtr, bitmapFileHeader.bmp_offset, SEEK_SET)) {
 fclose(filePtr);
 return NULL;
 }
 
 // allocate enough memory for the bitmap image data
 pixel_t *bitmapImage = malloc(bitmapInfoHeader->bmp_bytesz *
 sizeof(pixel_t));
 
 // verify memory allocation
 if (bitmapImage == NULL) {
 fclose(filePtr);
 return NULL;
 }
 
 // read in the bitmap image data
 size_t pad, count=0;
 unsigned char c;
 pad = 4*ceil(bitmapInfoHeader->bitspp*bitmapInfoHeader->width/32.) - bitmapInfoHeader->width;
 for(size_t i=0; i<bitmapInfoHeader->height; i++){
 for(size_t j=0; j<bitmapInfoHeader->width; j++){
 if (fread(&c, sizeof(unsigned char), 1, filePtr) != 1) {
 fclose(filePtr);
 return NULL;
 }
 bitmapImage[count++] = (pixel_t) c;
 }
 fseek(filePtr, pad, SEEK_CUR);
 }
 
 // If we were using unsigned char as pixel_t, then:
 // fread(bitmapImage, 1, bitmapInfoHeader->bmp_bytesz, filePtr);
 
 // close file and return bitmap image data
 fclose(filePtr);
 return bitmapImage;
}
 
// Return: true on error.
bool save_bmp(const char *filename, const bitmap_info_header_t *bmp_ih,
 const pixel_t *data)
{
 FILE* filePtr = fopen(filename, "wb");
 if (filePtr == NULL)
 return true;
 
 bmpfile_magic_t mag = {{0x42, 0x4d}};
 if (fwrite(&mag, sizeof(bmpfile_magic_t), 1, filePtr) != 1) {
 fclose(filePtr);
 return true;
 }
 
 const uint32_t offset = sizeof(bmpfile_magic_t) +
 sizeof(bmpfile_header_t) +
 sizeof(bitmap_info_header_t) +
 ((1U << bmp_ih->bitspp) * 4);
 
 const bmpfile_header_t bmp_fh = {
 .filesz = offset + bmp_ih->bmp_bytesz,
 .creator1 = 0,
 .creator2 = 0,
 .bmp_offset = offset
 };
 
 if (fwrite(&bmp_fh, sizeof(bmpfile_header_t), 1, filePtr) != 1) {
 fclose(filePtr);
 return true;
 }
 if (fwrite(bmp_ih, sizeof(bitmap_info_header_t), 1, filePtr) != 1) {
 fclose(filePtr);
 return true;
 }
 
 // Palette
 for (size_t i = 0; i < (1U << bmp_ih->bitspp); i++) {
 const rgb_t color = {(uint8_t)i, (uint8_t)i, (uint8_t)i};
 if (fwrite(&color, sizeof(rgb_t), 1, filePtr) != 1) {
 fclose(filePtr);
 return true;
 }
 }
 
 // We use int instead of uchar, so we can't write img
 // in 1 call any more.
 // fwrite(data, 1, bmp_ih->bmp_bytesz, filePtr);
 
 // Padding: http://en.wikipedia.org/wiki/BMP_file_format#Pixel_storage
 size_t pad = 4*ceil(bmp_ih->bitspp*bmp_ih->width/32.) - bmp_ih->width;
 unsigned char c;
 for(size_t i=0; i < bmp_ih->height; i++) {
 for(size_t j=0; j < bmp_ih->width; j++) {
 c = (unsigned char) data[j + bmp_ih->width*i];
 if (fwrite(&c, sizeof(char), 1, filePtr) != 1) {
 fclose(filePtr);
 return true;
 }
 }
 c = 0;
 for(size_t j=0; j<pad; j++)
 if (fwrite(&c, sizeof(char), 1, filePtr) != 1) {
 fclose(filePtr);
 return true;
 }
 }
 
 fclose(filePtr);
 return false;
}
 
// if normalize is true, map pixels to range 0..MAX_BRIGHTNESS
void convolution(const pixel_t *in, pixel_t *out, const float *kernel,
 const int nx, const int ny, const int kn,
 const bool normalize)
{
 assert(kn % 2 == 1);
 assert(nx > kn && ny > kn);
 const int khalf = kn / 2;
 float min = FLT_MAX, max = -FLT_MAX;
 
 if (normalize)
 for (int m = khalf; m < nx - khalf; m++)
 for (int n = khalf; n < ny - khalf; n++) {
 float pixel = 0.0;
 size_t c = 0;
 for (int j = -khalf; j <= khalf; j++)
 for (int i = -khalf; i <= khalf; i++) {
 pixel += in[(n - j) * nx + m - i] * kernel;
 c++;
 }
 if (pixel < min)
 min = pixel;
 if (pixel > max)
 max = pixel;
 }
 
 for (int m = khalf; m < nx - khalf; m++)
 for (int n = khalf; n < ny - khalf; n++) {
 float pixel = 0.0;
 size_t c = 0;
 for (int j = -khalf; j <= khalf; j++)
 for (int i = -khalf; i <= khalf; i++) {
 pixel += in[(n - j) * nx + m - i] * kernel;
 c++;
 }
 
 if (normalize)
 pixel = MAX_BRIGHTNESS * (pixel - min) / (max - min);
 out[n * nx + m] = (pixel_t)pixel;
 }
}
 
/*
 * gaussianFilter:
 * http://www.songho.ca/dsp/cannyedge/cannyedge.html
 * determine size of kernel (odd #)
 * 0.0 <= sigma < 0.5 : 3
 * 0.5 <= sigma < 1.0 : 5
 * 1.0 <= sigma < 1.5 : 7
 * 1.5 <= sigma < 2.0 : 9
 * 2.0 <= sigma < 2.5 : 11
 * 2.5 <= sigma < 3.0 : 13 ...
 * kernelSize = 2 * int(2*sigma) + 3;
 */
void gaussian_filter(const pixel_t *in, pixel_t *out,
 const int nx, const int ny, const float sigma)
{
 const int n = 2 * (int)(2 * sigma) + 3;
 const float mean = (float)floor(n / 2.0);
 float kernel[n * n]; // variable length array
 
 fprintf(stderr, "gaussian_filter: kernel size %d, sigma=%g\n",
 n, sigma);
 size_t c = 0;
 for (int i = 0; i < n; i++)
 for (int j = 0; j < n; j++) {
 kernel = exp(-0.5 * (pow((i - mean) / sigma, 2.0) +
 pow((j - mean) / sigma, 2.0)))
 / (2 * M_PI * sigma * sigma);
 c++;
 }
 
 convolution(in, out, kernel, nx, ny, n, true);
}
 
/*
 * Links:
 * http://en.wikipedia.org/wiki/Canny_edge_detector
 * http://www.tomgibara.com/computer-vision/CannyEdgeDetector.java
 * http://fourier.eng.hmc.edu/e161/lectures/canny/node1.html
 * http://www.songho.ca/dsp/cannyedge/cannyedge.html
 *
 * Note: T1 and T2 are lower and upper thresholds.
 */
pixel_t *canny_edge_detection(const pixel_t *in,
 const bitmap_info_header_t *bmp_ih,
 const int tmin, const int tmax,
 const float sigma)
{
 const int nx = bmp_ih->width;
 const int ny = bmp_ih->height;
 
 pixel_t *G = calloc(nx * ny * sizeof(pixel_t), 1);
 pixel_t *after_Gx = calloc(nx * ny * sizeof(pixel_t), 1);
 pixel_t *after_Gy = calloc(nx * ny * sizeof(pixel_t), 1);
 pixel_t *nms = calloc(nx * ny * sizeof(pixel_t), 1);
 pixel_t *out = malloc(bmp_ih->bmp_bytesz * sizeof(pixel_t));
 
 if (G == NULL || after_Gx == NULL || after_Gy == NULL ||
 nms == NULL || out == NULL) {
 fprintf(stderr, "canny_edge_detection:"
 " Failed memory allocation(s).\n");
 exit(1);
 }
 
 gaussian_filter(in, out, nx, ny, sigma);
 
 const float Gx[] = {-1, 0, 1,
 -2, 0, 2,
 -1, 0, 1};
 
 convolution(out, after_Gx, Gx, nx, ny, 3, false);
 
 const float Gy[] = { 1, 2, 1,
 0, 0, 0,
 -1,-2,-1};
 
 convolution(out, after_Gy, Gy, nx, ny, 3, false);
 
 for (int i = 1; i < nx - 1; i++)
 for (int j = 1; j < ny - 1; j++) {
 const int c = i + nx * j;
 // G = abs(after_Gx) + abs(after_Gy);
 G = (pixel_t)hypot(after_Gx, after_Gy);
 }
 
 // Non-maximum suppression, straightforward implementation.
 for (int i = 1; i < nx - 1; i++)
 for (int j = 1; j < ny - 1; j++) {
 const int c = i + nx * j;
 const int nn = c - nx;
 const int ss = c + nx;
 const int ww = c + 1;
 const int ee = c - 1;
 const int nw = nn + 1;
 const int ne = nn - 1;
 const int sw = ss + 1;
 const int se = ss - 1;
 
 const float dir = (float)(fmod(atan2(after_Gy,
 after_Gx) + M_PI,
 M_PI) / M_PI) * 8;
 
 if (((dir <= 1 || dir > 7) && G > G[ee] &&
 G > G[ww]) || // 0 deg
 ((dir > 1 && dir <= 3) && G > G[nw] &&
 G > G[se]) || // 45 deg
 ((dir > 3 && dir <= 5) && G > G[nn] &&
 G > G[ss]) || // 90 deg
 ((dir > 5 && dir <= 7) && G > G[ne] &&
 G > G[sw])) // 135 deg
 nms = G;
 else
 nms = 0;
 }
 
 // Reuse array
 // used as a stack. nx*ny/2 elements should be enough.
 int *edges = (int*) after_Gy;
 memset(out, 0, sizeof(pixel_t) * nx * ny);
 memset(edges, 0, sizeof(pixel_t) * nx * ny);
 
 // Tracing edges with hysteresis . Non-recursive implementation.
 size_t c = 1;
 for (int j = 1; j < ny - 1; j++)
 for (int i = 1; i < nx - 1; i++) {
 if (nms >= tmax && out == 0) { // trace edges
 out = MAX_BRIGHTNESS;
 int nedges = 1;
 edges[0] = c;
 
 do {
 nedges--;
 const int t = edges[nedges];
 
 int nbs[8]; // neighbours
 nbs[0] = t - nx; // nn
 nbs[1] = t + nx; // ss
 nbs[2] = t + 1; // ww
 nbs[3] = t - 1; // ee
 nbs[4] = nbs[0] + 1; // nw
 nbs[5] = nbs[0] - 1; // ne
 nbs[6] = nbs[1] + 1; // sw
 nbs[7] = nbs[1] - 1; // se
 
 for (int k = 0; k < 8; k++)
 if (nms[nbs[k]] >= tmin && out[nbs[k]] == 0) {
 out[nbs[k]] = MAX_BRIGHTNESS;
 edges[nedges] = nbs[k];
 nedges++;
 }
 } while (nedges > 0);
 }
 c++;
 }
 
 free(after_Gx);
 free(after_Gy);
 free(G);
 free(nms);
 
 return out;
}
 
int main(const int argc, const char ** const argv)
{
 if (argc < 2) {
 printf("Usage: %s image.bmp\n", argv[0]);
 return 1;
 }
 
 static bitmap_info_header_t ih;
 const pixel_t *in_bitmap_data = load_bmp(argv[1], &ih);
 if (in_bitmap_data == NULL) {
 fprintf(stderr, "main: BMP image not loaded.\n");
 return 1;
 }
 
 printf("Info: %d x %d x %d\n", ih.width, ih.height, ih.bitspp);
 
 const pixel_t *out_bitmap_data =
 canny_edge_detection(in_bitmap_data, &ih, 45, 50, 1.0f);
 if (out_bitmap_data == NULL) {
 fprintf(stderr, "main: failed canny_edge_detection.\n");
 return 1;
 }
 
 if (save_bmp("out.bmp", &ih, out_bitmap_data)) {
 fprintf(stderr, "main: BMP image not saved.\n");
 return 1;
 }
 
 free((pixel_t*)in_bitmap_data);
 free((pixel_t*)out_bitmap_data);
 return 0;
}

 

دانلود کد فوق از طریق لینک زیر:

Canny in C

رمز فایل : behsanandish.com

الگوریتم Canny

لبه یاب کنی توسط جان اف کنی در سال 1986 ایجاد شد و هنوز یک لبه یاب استاندارد و با دقت و کیفیت بالا میباشد.الگوریتم لبه یابی کنی یکی از بهترین لبه یابها تا به امروز است. در ادامه روش کار این الگوریتم و هم چنین کد الگوریتم Canny در #C را بررسی خواهیم کرد. این الگوریتم لبه یابی از سه بخش اصلی زیر تشکیل شده:

  • تضعیف نویز
  • پیدا کردن نقاطی که بتوان آنها را به عنوان لبه در نظر گرفت
  • حذب نقاطی که احتمال لبه بودن آنها کم است

 

معیارهایی که در لبه یاب کنی مطرح است:
1 -پایین آوردن نرخ خطا- یعنی تا حد امکان هیچ لبه ای در تصویر نباید گم شود و هم چنین هیچ چیزی که لبه نیست نباید به جای لبه فرض شود. لبه هان پیدا شده تا حد ممکن به لبه ها اصلی
نزدیک باشند.

2 -لبه در مکان واقعی خود باشد- یعنی تا حد ممکن لبه ها کمترین فاصله را با مکان واقعی خود داشته باشند.
3 -بران هر لبه فقط یک پاسخ داشته باشیم.

4 -لبه ها کمترین ضخامت را داشته باشند- (در صورت امکان یک پیکسل).
لبه یاب کنی بخاطر توانایی در تولید لبه های نازک تا حد یک ییکسل برای لبه های پیوسته معروف شده است. این لبه یاب شامل چهار مرحله و چهار ورودی زیر است:
یک تصویر ورودی
یک پارامتر به نام سیگما جهت مقدار نرم کنندگی تصویر
یک حد آستانه بالا (Th)
یک حد آستانه پایین (Tl)

 

مراحل الگوریتم Canny:

1- در ابتدا باید تصویر رنگی را به جهت لبه یابی بهتر به یک تصویر سطح خاکسترن تبدیب کرد.

2- نویز را از تصویر دریافتی حذف کرد. بدلیل اینکه فیلتر گاوسین از یک ماسک ساده برای حذف نویز استفاده می کند لبه یاب کنی در مرحله اول برای حذف نویز آن را بکار میگیرد.

3- در یک تصویر سطح خاکستر جایی را که بیشترین تغییرات را داشته باشند به عنوان لبه در نظر گرفته می شوند و این مکانها با گرفتن گرادیان تصویر با استفاده عملگر سوبل بدست می آیند. سپس لبه های مات یافت شده به لبه های تیزتر تبدیل می شوند.

4- برخی از لبه های کشف شده واقعا لبه نیستند و در واقع نویز هستند که باید آنها توسط حد آستانه هیسترزیس فیلتر شوند.هیسترزیس از دو حد آستانه بالاتر (Th) و حد آستانه پایین تر (Tl) استفاده کرده و کنی پیشنهاد می کند که نسبت استانه بالا به پایین سه به یک باشد.

 این روش بیشتر به کشف لبه های ضعیف به درستی می پردازد و کمتر فریب نویز را می خورد و از بقیه روش ها بهتر است.

 

الگوریتم Canny    عملکرد الگوریتم Canny

 


 

کد الگوریتم Canny در #C:

الگوریتم در 5 مرحله جداگانه اجرا می شود:

1. صاف کردن: تار شدن تصویر برای حذف نویز. پیکربندی توسط فیلتر گاوسی با اندازه مشخص هسته (N) و پارامتر پوشش گاوسی سیگما. پوشاننده فیلتر گاوسی توسط تابع زیر تولید می شود:

private void GenerateGaussianKernel(int N, float S ,out int Weight)
{

float Sigma = S ;
float pi;
pi = (float)Math.PI;
int i, j;
int SizeofKernel=N;

float [,] Kernel = new float [N,N];
GaussianKernel = new int [N,N];
float[,] OP = new float[N, N];
float D1,D2;

D1= 1/(2*pi*Sigma*Sigma);
D2= 2*Sigma*Sigma;

float min=1000;

for (i = -SizeofKernel / 2; i <= SizeofKernel / 2; i++)
{
for (j = -SizeofKernel / 2; j <= SizeofKernel / 2; j++)
{
Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] = ((1 / D1) * (float)Math.Exp(-(i * i + j * j) / D2));
if (Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] < min)
min = Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j];

}
}
int mult = (int)(1 / min);
int sum = 0;
if ((min > 0) && (min < 1))
{

for (i = -SizeofKernel / 2; i <= SizeofKernel / 2; i++)
{
for (j = -SizeofKernel / 2; j <= SizeofKernel / 2; j++)
{
Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] = (float)Math.Round(Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] * mult, 0);
GaussianKernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] = (int)Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j];
sum = sum + GaussianKernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j];
}

}

}
else
{
sum = 0;
for (i = -SizeofKernel / 2; i <= SizeofKernel / 2; i++)
{
for (j = -SizeofKernel / 2; j <= SizeofKernel / 2; j++)
{
Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] = (float)Math.Round(Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] , 0);
GaussianKernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j] = (int)Kernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j];
sum = sum + GaussianKernel[SizeofKernel / 2 + i, SizeofKernel / 2 + j];
}

}

}
//Normalizing kernel Weight
Weight= sum;

return;
}

 

زیر روال ذیل نویز را توسط فیلتر گوسی حذف می کند.

 

private int[,] GaussianFilter(int[,] Data)
        {
            GenerateGaussianKernel(KernelSize, Sigma,out KernelWeight);

            int[,] Output = new int[Width, Height];
            int i, j,k,l;
            int Limit = KernelSize /2;

            float Sum=0;

 Output = Data; // Removes Unwanted Data Omission due to kernel bias while convolution

            for (i = Limit; i <= ((Width - 1) - Limit); i++)
            {
                for (j = Limit; j <= ((Height - 1) - Limit); j++)
                {
                    Sum = 0;
                    for (k = -Limit; k <= Limit; k++)
                    {

                       for (l = -Limit; l <= Limit; l++)
                        {
                            Sum = Sum + ((float)Data[i + k, j + l] * GaussianKernel [Limit + k, Limit + l]);                        

                        }
                    }
                    Output[i, j] = (int)(Math.Round(Sum/ (float)KernelWeight));
                }

            }

            return Output;
        }

 

2. پیدا کردن شیب ها: لبه ها باید مشخص شوند، جایی که شیب های تصویر بزرگ می شوند.

ماسک های سوبل  X و Y برای تولید گرادیان های تصویر X و Y استفاده می شود؛ تابع بعدی تمایز را با استفاده از فیلتر ماسک sobel اعمال می کند.

 

private float[,] Differentiate(int[,] Data, int[,] Filter)
        {
            int i, j,k,l, Fh, Fw;

            Fw = Filter.GetLength(0);
            Fh = Filter.GetLength(1);
            float sum = 0;
            float[,] Output = new float[Width, Height];

            for (i = Fw / 2; i <= (Width - Fw / 2) - 1; i++)
            {
                for (j = Fh / 2; j <= (Height  - Fh / 2) - 1; j++)
                {
                  sum=0;
                   for(k=-Fw/2; k<=Fw/2; k++)
                   {
                       for(l=-Fh/2; l<=Fh/2; l++)
                       {
                          sum=sum + Data[i+k,j+l]*Filter[Fw/2+k,Fh/2+l];


                       }
                   }
                    Output[i,j]=sum;

                }

            }
            return Output;

        }

 

3. توقیف غیر حداکثر: فقط حداکثرهای محلی باید به عنوان لبه ها مشخص شود.

ما جهت گرادیان را پیدا می کنیم و با استفاده از این جهت، ما توقیف غیر حداکثر را انجام می دهیم (“پردازش تصویر دیجیتال- آموزش توسط گنزالس-پیرسون ” را بخوانید)

 

// Perform Non maximum suppression:
           // NonMax = Gradient;

            for (i = 0; i <= (Width - 1); i++)
            {
                for (j = 0; j <= (Height - 1); j++)
                {
                    NonMax[i, j] = Gradient[i, j];
                }
            }
     
            int Limit = KernelSize / 2;
            int r, c;
            float Tangent;

                for (i = Limit; i <= (Width - Limit) - 1; i++)
            {
                for (j = Limit; j <= (Height - Limit) - 1; j++)
                {

                    if (DerivativeX[i, j] == 0)
                        Tangent = 90F;
                    else
                        Tangent = (float)(Math.Atan(DerivativeY[i, j] / DerivativeX[i, j]) * 180 / Math.PI); //rad to degree



                    //Horizontal Edge
                    if (((-22.5 < Tangent) && (Tangent <= 22.5)) || ((157.5 < Tangent) && (Tangent <= -157.5)))
                    {
                        if ((Gradient[i, j] < Gradient[i, j + 1]) || (Gradient[i, j] < Gradient[i, j - 1]))
                            NonMax[i, j] = 0;
                    }

                    //Vertical Edge
                    if (((-112.5 < Tangent) && (Tangent <= -67.5)) || ((67.5 < Tangent) && (Tangent <= 112.5)))
                    {
                        if ((Gradient[i, j] < Gradient[i + 1, j]) || (Gradient[i, j] < Gradient[i - 1, j]))
                            NonMax[i, j] = 0;
                    }

                    //+45 Degree Edge
                    if (((-67.5 < Tangent) && (Tangent <= -22.5)) || ((112.5 < Tangent) && (Tangent <= 157.5)))
                    {
                        if ((Gradient[i, j] < Gradient[i + 1, j - 1]) || (Gradient[i, j] < Gradient[i - 1, j + 1]))
                            NonMax[i, j] = 0;
                    }

                    //-45 Degree Edge
                    if (((-157.5 < Tangent) && (Tangent <= -112.5)) || ((67.5 < Tangent) && (Tangent <= 22.5)))
                    {
                        if ((Gradient[i, j] < Gradient[i + 1, j + 1]) || (Gradient[i, j] < Gradient[i - 1, j - 1]))
                            NonMax[i, j] = 0;
                    }

                }

            }

 

4. آستانه دوگانه: لبه های بالقوه توسط آستانه تعیین می شود.

5. ردیابی لبه توسط هیسترسیس: لبه های نهایی توسط توقیف تمام لبه هایی که به یک لبه بسیار قطعی (قوی) متصل نیستند، مشخص می شوند.

 

private void HysterisisThresholding(int[,] Edges)
        {

            int i, j;
            int Limit= KernelSize/2;


            for (i = Limit; i <= (Width - 1) - Limit; i++)
                for (j = Limit; j <= (Height - 1) - Limit; j++)
                {
                    if (Edges[i, j] == 1)
                    {
                        EdgeMap[i, j] = 1;

                    }

                }

            for (i = Limit; i <= (Width - 1) - Limit; i++)
            {
                for (j = Limit; j <= (Height  - 1) - Limit; j++)
                {
                    if (Edges[i, j] == 1)
                    {
                        EdgeMap[i, j] = 1;
                        Travers(i, j);
                        VisitedMap[i, j] = 1;
                    }
                }
            }




            return;
        }

//Recursive Procedure 
private void Travers(int X, int Y)
        {

            
            if (VisitedMap[X, Y] == 1)
            {
                return;
            }

            //1
            if (EdgePoints[X + 1, Y] == 2)
            {
                EdgeMap[X + 1, Y] = 1;
                VisitedMap[X + 1, Y] = 1;
                Travers(X + 1, Y);
                return;
            }
            //2
            if (EdgePoints[X + 1, Y - 1] == 2)
            {
                EdgeMap[X + 1, Y - 1] = 1;
                VisitedMap[X + 1, Y - 1] = 1;
                Travers(X + 1, Y - 1);
                return;
            }

           //3

            if (EdgePoints[X, Y - 1] == 2)
            {
                EdgeMap[X , Y - 1] = 1;
                VisitedMap[X , Y - 1] = 1;
                Travers(X , Y - 1);
                return;
            }

           //4

            if (EdgePoints[X - 1, Y - 1] == 2)
            {
                EdgeMap[X - 1, Y - 1] = 1;
                VisitedMap[X - 1, Y - 1] = 1;
                Travers(X - 1, Y - 1);
                return;
            }
            //5
            if (EdgePoints[X - 1, Y] == 2)
            {
                EdgeMap[X - 1, Y ] = 1;
                VisitedMap[X - 1, Y ] = 1;
                Travers(X - 1, Y );
                return;
            }
            //6
            if (EdgePoints[X - 1, Y + 1] == 2)
            {
                EdgeMap[X - 1, Y + 1] = 1;
                VisitedMap[X - 1, Y + 1] = 1;
                Travers(X - 1, Y + 1);
                return;
            }
            //7
            if (EdgePoints[X, Y + 1] == 2)
            {
                EdgeMap[X , Y + 1] = 1;
                VisitedMap[X, Y + 1] = 1;
                Travers(X , Y + 1);
                return;
            }
            //8

            if (EdgePoints[X + 1, Y + 1] == 2)
            {
                EdgeMap[X + 1, Y + 1] = 1;
                VisitedMap[X + 1, Y + 1] = 1;
                Travers(X + 1, Y + 1);
                return;
            }


            //VisitedMap[X, Y] = 1;
            return;

        } 
          
        //Canny Class Ends

    }

 

این کار با یک تابع بازگشتی انجام می شود که آستانه دوگانه را با دو آستانه بالا (Threshold (TH و (Low Threshold (TL و تجزیه و تحلیل 8-اتصال انجام می دهد.

 

دانلود کد فوق از طریق لینک زیر:

Canny Edge Detection C#

رمز فایل : behsanandish.com