Last active
March 22, 2019 13:26
-
-
Save eudoxos/b88f36e27eb782c6c60f25481c642aaf to your computer and use it in GitHub Desktop.
sobel gradient filter kernel (along X), integral and floating-point (c++17, Eigen)
This file contains hidden or bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
| #include<Eigen/Core> | |
| #include<iostream> | |
| #include<numeric> | |
| // based on and verified with https://stackoverflow.com/a/46262628/761090 | |
| // works for size 11 with int and 15 with long int (on amd64) | |
| template<typename Int> | |
| std::tuple<Eigen::Matrix<Int,Eigen::Dynamic,Eigen::Dynamic>,Int> sobelGradFilterX_int(int sz){ | |
| assert(sz>0); | |
| assert(sz%2==1); | |
| Int sz2=sz/2; | |
| Int denom=1; | |
| for(Int i=0; i<=sz2; i++){ | |
| for(Int j=0; j<=sz2; j++){ | |
| Int d=(Int)(i*i+j*j); | |
| if(d==0) continue; | |
| denom=std::lcm(denom,d/std::gcd(d,i)); | |
| } | |
| } | |
| Eigen::Matrix<Int,Eigen::Dynamic,Eigen::Dynamic> ret; ret.setZero(sz,sz); | |
| Int ww=0; | |
| for(int i0=0; i0<sz; i0++){ | |
| for(int j0=0; j0<sz; j0++){ | |
| int i=i0-sz2, j=j0-sz2; | |
| Int d=(Int)(i*i+j*j); | |
| if(i==0) continue; | |
| assert((denom*i)%d==0); | |
| Int w=denom*i/d; | |
| ww+=i*w; | |
| ret(j0,i0)=w; | |
| } | |
| } | |
| return std::make_tuple(ret,ww); | |
| } | |
| template<typename Scalar> | |
| Eigen::Matrix<Scalar,Eigen::Dynamic,Eigen::Dynamic> sobelGradFilterX(int sz){ | |
| assert(sz>0); | |
| assert(sz%2==1); | |
| int sz2=sz/2; | |
| int denom=1; | |
| for(int i=0; i<=sz2; i++){ | |
| for(int j=0; j<=sz2; j++){ | |
| int d=(int)(i*i+j*j); | |
| if(d==0) continue; | |
| denom=std::lcm(denom,d/std::gcd(d,i)); | |
| } | |
| } | |
| Eigen::Matrix<Scalar,Eigen::Dynamic,Eigen::Dynamic> ret; ret.setZero(sz,sz); | |
| Scalar ww=0; | |
| for(int i0=0; i0<sz; i0++){ | |
| for(int j0=0; j0<sz; j0++){ | |
| int i=i0-sz2, j=j0-sz2; | |
| Scalar d=i*i+j*j; | |
| if(i==0) continue; | |
| Scalar w=denom*i*1./d; | |
| ww+=i*w; | |
| ret(j0,i0)=i*1./d; | |
| } | |
| } | |
| return denom*ret/ww; | |
| } | |
| int main(void){ | |
| auto [K,w]=sobelGradFilterX_int<int>(11); | |
| std::cerr<<K<<std::endl<<w<<std::endl; | |
| auto K2=sobelGradFilterX<double>(13); | |
| std::cerr<<K2<<std::endl; | |
| } |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment
Example output from
g++ mat2.cc -std=c++17 -Wall -O3 -o mat2 -I/usr/include/eigen3 && ./mat2