WebSVN – seema-scanner – Diff – /src/algorithm/algorithmtools.h


#ifndef ALGORITHMTOOLS_H
#define ALGORITHMTOOLS_H
 
#include <opencv2/opencv.hpp>
#include <opencv2/core/ocl.hpp>
#include <QTime>
 
// Convert an unsigned binary number to reflected binary Gray code.
// Source: http://en.wikipedia.org/wiki/Gray_code
inline unsigned int binaryToGray(unsigned int num) {
    return (num >> 1) ^ num;
}
 
// Convert a reflected binary Gray code number to a binary number.
// Source: http://en.wikipedia.org/wiki/Gray_code
inline unsigned int grayToBinary(unsigned int num){
    unsigned int mask;
    for(mask = num >> 1; mask != 0; mask = mask >> 1)
        num = num ^ mask;
    return num;
}
 
// Return the Nth bit of an unsigned integer number
inline bool getBit(int decimal, int N){
 
    return decimal & 1 << (N-1);
}
 
// Return the position of the least significant bit that is set
inline int leastSignificantBitSet(int x){
  if(x == 0)
      return 0;
 
  int val = 1;
  while(x>>=1)
      val++;
 
  return val;
}
 
// Compute the power of a number (where the exponent is an integer)
inline unsigned int powi(int num, unsigned int exponent){
 
    if(exponent == 0)
        return 1;
 
    float res = num;
    for(unsigned int i=0; i<exponent-1; i++)
        res *= num;
 
    return res;
}
 
// Compute the power of a 2 (where the exponent is an integer)
inline unsigned int twopowi(unsigned int exponent){
    return 1 << exponent;
}
 
inline cv::Vec3b getColorSubpix(const cv::Mat& img, cv::Point2f pt){
    assert(!img.empty());
    assert(img.channels() == 3);
 
    int x = (int)pt.x;
    int y = (int)pt.y;
 
    int x0 = cv::borderInterpolate(x,   img.cols, cv::BORDER_REFLECT_101);
    int x1 = cv::borderInterpolate(x+1, img.cols, cv::BORDER_REFLECT_101);
    int y0 = cv::borderInterpolate(y,   img.rows, cv::BORDER_REFLECT_101);
    int y1 = cv::borderInterpolate(y+1, img.rows, cv::BORDER_REFLECT_101);
 
    float a = pt.x - (float)x;
    float c = pt.y - (float)y;
 
    uchar b = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[0] * a) * (1.f - c)
                           + (img.at<cv::Vec3b>(y1, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[0] * a) * c);
    uchar g = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[1] * a) * (1.f - c)
                           + (img.at<cv::Vec3b>(y1, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[1] * a) * c);
    uchar r = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[2] * a) * (1.f - c)
                           + (img.at<cv::Vec3b>(y1, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[2] * a) * c);
 
    return cv::Vec3b(b, g, r);
}
 
inline cv::Mat computePhaseVector(unsigned int length, float phase, float pitch){
 
    cv::Mat phaseVector(length, 1, CV_8UC3);
    //phaseVector.setTo(0);
 
    const float pi = M_PI;
 
    // Loop through vector
    for(int i=0; i<phaseVector.rows; i++){
        // Amplitude of channels
        float amp = 0.5*(1+cos(2*pi*i/pitch - phase));
        phaseVector.at<cv::Vec3b>(i, 0) = cv::Vec3b(255.0*amp,255.0*amp,255.0*amp);
    }
 
    return phaseVector;
}
 
inline cv::Mat getPhase(const cv::Mat I1, const cv::Mat I2, const cv::Mat I3){
 
    cv::Mat_<float> I1_(I1);
    cv::Mat_<float> I2_(I2);
    cv::Mat_<float> I3_(I3);
 
    cv::Mat phase;
 
    // One call approach
    cv::phase(2.0*I1_-I3_-I2_, sqrt(3.0)*(I2_-I3_), phase);
    return phase;
 
}
 
// Phase unwrapping by means of a phase cue
inline cv::Mat unwrapWithCue(const cv::Mat &up, const cv::Mat &upCue, float nPhases){
 
    const float pi = M_PI;
 
    // Determine fringe order
    cv::Mat P = (upCue*nPhases-up)/(2.0*pi);
 
    // Round to integers
    P.convertTo(P, CV_8U);
    P.convertTo(P, CV_32F);
 
    // Add to phase
    cv::Mat upUnwrapped = up + P*2*pi;
 
    // Scale to range [0; 2pi]
    upUnwrapped *= 1.0/nPhases;
 
    return upUnwrapped;
}
 
// Absolute phase and magnitude from N frames
inline std::vector<cv::Mat> getDFTComponents(const std::vector<cv::Mat> &frames){
 
    unsigned int N = frames.size();
 
//    std::vector<cv::Mat> framesReverse = frames;
//    std::reverse(framesReverse.begin(), framesReverse.end());
 
    // DFT approach
    cv::Mat I;
    cv::merge(frames, I);
    unsigned int w = I.cols;
    unsigned int h = I.rows;
    I = I.reshape(1, h*w);
    I.convertTo(I, CV_32F);
 
    QTime time;
    time.start();
 
    cv::Mat fI;
    cv::dft(I, fI, cv::DFT_ROWS + cv::DFT_COMPLEX_OUTPUT);
 
    std::cout << "elapsed: " << time.restart() << std::endl;
 
    fI = fI.reshape(N*2, h);
 
    std::vector<cv::Mat> fIcomp;
    cv::split(fI, fIcomp);
 
    return fIcomp;
 
}
 
#endif // ALGORITHMTOOLS_H
 

Rev 207	Rev 236
1	`#ifndef ALGORITHMTOOLS_H`	1	`#ifndef ALGORITHMTOOLS_H`
2	`#define ALGORITHMTOOLS_H`	2	`#define ALGORITHMTOOLS_H`
3		3
4	`#include <opencv2/opencv.hpp>`	4	`#include <opencv2/opencv.hpp>`
-		5	`#include <opencv2/core/ocl.hpp>`
-		6	`#include <QTime>`
5		7
6	`// Convert an unsigned binary number to reflected binary Gray code.`	8	`// Convert an unsigned binary number to reflected binary Gray code.`
7	`// Source: http://en.wikipedia.org/wiki/Gray_code`	9	`// Source: http://en.wikipedia.org/wiki/Gray_code`
8	`inline unsigned int binaryToGray(unsigned int num) {`	10	`inline unsigned int binaryToGray(unsigned int num) {`
9	`return (num >> 1) ^ num;`	11	`return (num >> 1) ^ num;`
10	`}`	12	`}`
11		13
12	`// Convert a reflected binary Gray code number to a binary number.`	14	`// Convert a reflected binary Gray code number to a binary number.`
13	`// Source: http://en.wikipedia.org/wiki/Gray_code`	15	`// Source: http://en.wikipedia.org/wiki/Gray_code`
14	`inline unsigned int grayToBinary(unsigned int num){`	16	`inline unsigned int grayToBinary(unsigned int num){`
15	`unsigned int mask;`	17	`unsigned int mask;`
16	`for(mask = num >> 1; mask != 0; mask = mask >> 1)`	18	`for(mask = num >> 1; mask != 0; mask = mask >> 1)`
17	`num = num ^ mask;`	19	`num = num ^ mask;`
18	`return num;`	20	`return num;`
19	`}`	21	`}`
20		22
21	`// Return the Nth bit of an unsigned integer number`	23	`// Return the Nth bit of an unsigned integer number`
22	`inline bool getBit(int decimal, int N){`	24	`inline bool getBit(int decimal, int N){`
23		25
24	`return decimal & 1 << (N-1);`	26	`return decimal & 1 << (N-1);`
25	`}`	27	`}`
26		28
27	`// Return the position of the least significant bit that is set`	29	`// Return the position of the least significant bit that is set`
28	`inline int leastSignificantBitSet(int x){`	30	`inline int leastSignificantBitSet(int x){`
29	`if(x == 0)`	31	`if(x == 0)`
30	`return 0;`	32	`return 0;`
31		33
32	`int val = 1;`	34	`int val = 1;`
33	`while(x>>=1)`	35	`while(x>>=1)`
34	`val++;`	36	`val++;`
35		37
36	`return val;`	38	`return val;`
37	`}`	39	`}`
38		40
39	`// Compute the power of a number (where the exponent is an integer)`	41	`// Compute the power of a number (where the exponent is an integer)`
40	`inline unsigned int powi(int num, unsigned int exponent){`	42	`inline unsigned int powi(int num, unsigned int exponent){`
41		43
42	`if(exponent == 0)`	44	`if(exponent == 0)`
43	`return 1;`	45	`return 1;`
44		46
45	`float res = num;`	47	`float res = num;`
46	`for(unsigned int i=0; i<exponent-1; i++)`	48	`for(unsigned int i=0; i<exponent-1; i++)`
47	`res *= num;`	49	`res *= num;`
48		50
49	`return res;`	51	`return res;`
50	`}`	52	`}`
51		53
52	`// Compute the power of a 2 (where the exponent is an integer)`	54	`// Compute the power of a 2 (where the exponent is an integer)`
53	`inline unsigned int twopowi(unsigned int exponent){`	55	`inline unsigned int twopowi(unsigned int exponent){`
54	`return 1 << exponent;`	56	`return 1 << exponent;`
55	`}`	57	`}`
56		58
57	`inline cv::Vec3b getColorSubpix(const cv::Mat& img, cv::Point2f pt){`	59	`inline cv::Vec3b getColorSubpix(const cv::Mat& img, cv::Point2f pt){`
58	`assert(!img.empty());`	60	`assert(!img.empty());`
59	`assert(img.channels() == 3);`	61	`assert(img.channels() == 3);`
60		62
61	`int x = (int)pt.x;`	63	`int x = (int)pt.x;`
62	`int y = (int)pt.y;`	64	`int y = (int)pt.y;`
63		65
64	`int x0 = cv::borderInterpolate(x, img.cols, cv::BORDER_REFLECT_101);`	66	`int x0 = cv::borderInterpolate(x, img.cols, cv::BORDER_REFLECT_101);`
65	`int x1 = cv::borderInterpolate(x+1, img.cols, cv::BORDER_REFLECT_101);`	67	`int x1 = cv::borderInterpolate(x+1, img.cols, cv::BORDER_REFLECT_101);`
66	`int y0 = cv::borderInterpolate(y, img.rows, cv::BORDER_REFLECT_101);`	68	`int y0 = cv::borderInterpolate(y, img.rows, cv::BORDER_REFLECT_101);`
67	`int y1 = cv::borderInterpolate(y+1, img.rows, cv::BORDER_REFLECT_101);`	69	`int y1 = cv::borderInterpolate(y+1, img.rows, cv::BORDER_REFLECT_101);`
68		70
69	`float a = pt.x - (float)x;`	71	`float a = pt.x - (float)x;`
70	`float c = pt.y - (float)y;`	72	`float c = pt.y - (float)y;`
71		73
72	`uchar b = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[0] * a) * (1.f - c)`	74	`uchar b = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[0] * a) * (1.f - c)`
73	`+ (img.at<cv::Vec3b>(y1, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[0] * a) * c);`	75	`+ (img.at<cv::Vec3b>(y1, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[0] * a) * c);`
74	`uchar g = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[1] * a) * (1.f - c)`	76	`uchar g = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[1] * a) * (1.f - c)`
75	`+ (img.at<cv::Vec3b>(y1, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[1] * a) * c);`	77	`+ (img.at<cv::Vec3b>(y1, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[1] * a) * c);`
76	`uchar r = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[2] * a) * (1.f - c)`	78	`uchar r = (uchar)cvRound((img.at<cv::Vec3b>(y0, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y0, x1)[2] * a) * (1.f - c)`
77	`+ (img.at<cv::Vec3b>(y1, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[2] * a) * c);`	79	`+ (img.at<cv::Vec3b>(y1, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[2] * a) * c);`
78		80
79	`return cv::Vec3b(b, g, r);`	81	`return cv::Vec3b(b, g, r);`
80	`}`	82	`}`
81		83
82	`inline cv::Mat computePhaseVector(unsigned int length, float phase, float pitch){`	84	`inline cv::Mat computePhaseVector(unsigned int length, float phase, float pitch){`
83		85
84	`cv::Mat phaseVector(length, 1, CV_8UC3);`	86	`cv::Mat phaseVector(length, 1, CV_8UC3);`
85	`//phaseVector.setTo(0);`	87	`//phaseVector.setTo(0);`
86		88
87	`const float pi = M_PI;`	89	`const float pi = M_PI;`
88		90
89	`// Loop through vector`	91	`// Loop through vector`
90	`for(int i=0; i<phaseVector.rows; i++){`	92	`for(int i=0; i<phaseVector.rows; i++){`
91	`// Amplitude of channels`	93	`// Amplitude of channels`
92	`float amp = 0.5(1+cos(2pi*i/pitch - phase));`	94	`float amp = 0.5(1+cos(2pi*i/pitch - phase));`
93	`phaseVector.at<cv::Vec3b>(i, 0) = cv::Vec3b(255.0amp,255.0amp,255.0*amp);`	95	`phaseVector.at<cv::Vec3b>(i, 0) = cv::Vec3b(255.0amp,255.0amp,255.0*amp);`
94	`}`	96	`}`
95		97
96	`return phaseVector;`	98	`return phaseVector;`
97	`}`	99	`}`
98		100
99	`inline cv::Mat getPhase(const cv::Mat I1, const cv::Mat I2, const cv::Mat I3){`	101	`inline cv::Mat getPhase(const cv::Mat I1, const cv::Mat I2, const cv::Mat I3){`
100		102
101	`cv::Mat_<float> I1_(I1);`	103	`cv::Mat_<float> I1_(I1);`
102	`cv::Mat_<float> I2_(I2);`	104	`cv::Mat_<float> I2_(I2);`
103	`cv::Mat_<float> I3_(I3);`	105	`cv::Mat_<float> I3_(I3);`
104		106
105	`cv::Mat phase;`	107	`cv::Mat phase;`
106		108
107	`// One call approach`	109	`// One call approach`
108	`cv::phase(2.0I1_-I3_-I2_, sqrt(3.0)(I2_-I3_), phase);`	110	`cv::phase(2.0I1_-I3_-I2_, sqrt(3.0)(I2_-I3_), phase);`
109	`return phase;`	111	`return phase;`
110		112
111	`}`	113	`}`
112		114
113	`// Phase unwrapping by means of a phase cue`	115	`// Phase unwrapping by means of a phase cue`
114	`inline cv::Mat unwrapWithCue(const cv::Mat &up, const cv::Mat &upCue, float nPhases){`	116	`inline cv::Mat unwrapWithCue(const cv::Mat &up, const cv::Mat &upCue, float nPhases){`
115		117
116	`const float pi = M_PI;`	118	`const float pi = M_PI;`
117		119
118	`// Determine fringe order`	120	`// Determine fringe order`
119	`cv::Mat P = (upCuenPhases-up)/(2.0pi);`	121	`cv::Mat P = (upCuenPhases-up)/(2.0pi);`
120		122
121	`// Round to integers`	123	`// Round to integers`
122	`P.convertTo(P, CV_8U);`	124	`P.convertTo(P, CV_8U);`
123	`P.convertTo(P, CV_32F);`	125	`P.convertTo(P, CV_32F);`
124		126
125	`// Add to phase`	127	`// Add to phase`
126	`cv::Mat upUnwrapped = up + P2pi;`	128	`cv::Mat upUnwrapped = up + P2pi;`
127		129
128	`// Scale to range [0; 2pi]`	130	`// Scale to range [0; 2pi]`
129	`upUnwrapped *= 1.0/nPhases;`	131	`upUnwrapped *= 1.0/nPhases;`
130		132
131	`return upUnwrapped;`	133	`return upUnwrapped;`
132	`}`	134	`}`
133		135
134	`// Absolute phase and magnitude from N frames`	136	`// Absolute phase and magnitude from N frames`
135	`inline std::vector<cv::Mat> getDFTComponents(const std::vector<cv::Mat> &frames){`	137	`inline std::vector<cv::Mat> getDFTComponents(const std::vector<cv::Mat> &frames){`
136		138
137	`unsigned int N = frames.size();`	139	`unsigned int N = frames.size();`
138		140
139	`// std::vector<cv::Mat> framesReverse = frames;`	141	`// std::vector<cv::Mat> framesReverse = frames;`
140	`// std::reverse(framesReverse.begin(), framesReverse.end());`	142	`// std::reverse(framesReverse.begin(), framesReverse.end());`
141		143
142	`// DFT approach`	144	`// DFT approach`
143	`cv::Mat I;`	145	`cv::Mat I;`
144	`cv::merge(frames, I);`	146	`cv::merge(frames, I);`
145	`unsigned int w = I.cols;`	147	`unsigned int w = I.cols;`
146	`unsigned int h = I.rows;`	148	`unsigned int h = I.rows;`
147	`I = I.reshape(1, h*w);`	149	`I = I.reshape(1, h*w);`
148	`I.convertTo(I, CV_32F);`	150	`I.convertTo(I, CV_32F);`
-		151
-		152	`QTime time;`
-		153	`time.start();`
-		154
149	`cv::Mat fI;`	155	`cv::Mat fI;`
150	`cv::dft(I, fI, cv::DFT_ROWS + cv::DFT_COMPLEX_OUTPUT);`	156	`cv::dft(I, fI, cv::DFT_ROWS + cv::DFT_COMPLEX_OUTPUT);`
-		157
-		158	`std::cout << "elapsed: " << time.restart() << std::endl;`
-		159
151	`fI = fI.reshape(N*2, h);`	160	`fI = fI.reshape(N*2, h);`
152		161
153	`std::vector<cv::Mat> fIcomp;`	162	`std::vector<cv::Mat> fIcomp;`
154	`cv::split(fI, fIcomp);`	163	`cv::split(fI, fIcomp);`
155		164
156	`return fIcomp;`	165	`return fIcomp;`
157		166
158	`}`	167	`}`
159		168
160	`#endif // ALGORITHMTOOLS_H`	169	`#endif // ALGORITHMTOOLS_H`
161		170

Subversion Repositories seema-scanner

(root)/src/algorithm/algorithmtools.h – Rev 207 → 236