Subversion Repositories seema-scanner

Rev

Rev 126 | Rev 128 | Go to most recent revision | Only display areas with differences | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 126 Rev 127
1
#include "AlgorithmLineShift.h"
1
#include "AlgorithmLineShift.h"
2
#include <cmath>
2
#include <cmath>
3
#include "cvtools.h"
3
#include "cvtools.h"
4
 
4
 
5
#include <opencv2/imgproc/imgproc.hpp>
5
#include <opencv2/imgproc/imgproc.hpp>
6
 
6
 
7
#ifndef log2f
7
#ifndef log2f
8
#define log2f(x) (log(x)/log(2.0))
8
#define log2f(x) (log(x)/log(2.0))
9
#endif
9
#endif
10
 
10
 
11
static unsigned int nLineShifts = 20; // number of columns over which each line is shifted
11
static unsigned int nLineShifts = 8; // number of columns over which each line is shifted
12
 
12
 
13
/*
13
/*
14
 * The purpose of this function is to convert an unsigned
14
 * The purpose of this function is to convert an unsigned
15
 * binary number to reflected binary Gray code.
15
 * binary number to reflected binary Gray code.
16
 *
16
 *
17
 * The operator >> is shift right. The operator ^ is exclusive or.
17
 * The operator >> is shift right. The operator ^ is exclusive or.
18
 * Source: http://en.wikipedia.org/wiki/Gray_code
18
 * Source: http://en.wikipedia.org/wiki/Gray_code
19
 */
19
 */
20
static unsigned int binaryToGray(unsigned int num) {
20
static unsigned int binaryToGray(unsigned int num) {
21
    return (num >> 1) ^ num;
21
    return (num >> 1) ^ num;
22
}
22
}
23
 
23
 
24
/*
24
/*
25
 * From Wikipedia: http://en.wikipedia.org/wiki/Gray_code
25
 * From Wikipedia: http://en.wikipedia.org/wiki/Gray_code
26
 * The purpose of this function is to convert a reflected binary
26
 * The purpose of this function is to convert a reflected binary
27
 * Gray code number to a binary number.
27
 * Gray code number to a binary number.
28
 */
28
 */
29
static unsigned int grayToBinary(unsigned int num){
29
static unsigned int grayToBinary(unsigned int num){
30
    unsigned int mask;
30
    unsigned int mask;
31
    for(mask = num >> 1; mask != 0; mask = mask >> 1)
31
    for(mask = num >> 1; mask != 0; mask = mask >> 1)
32
        num = num ^ mask;
32
        num = num ^ mask;
33
    return num;
33
    return num;
34
}
34
}
35
 
35
 
36
/*
36
/*
37
 * Return the Nth bit of an unsigned integer number
37
 * Return the Nth bit of an unsigned integer number
38
 */
38
 */
39
static bool getBit(int decimal, int N){
39
static bool getBit(int decimal, int N){
40
 
40
 
41
    return decimal & 1 << (N-1);
41
    return decimal & 1 << (N-1);
42
}
42
}
43
 
43
 
44
/*
44
/*
45
 * Return the number of bits set in an integer
45
 * Return the number of bits set in an integer
46
 */
46
 */
47
static int countBits(int n) {
47
static int countBits(int n) {
48
  unsigned int c; // c accumulates the total bits set in v
48
  unsigned int c; // c accumulates the total bits set in v
49
  for (c = 0; n>0; c++)
49
  for (c = 0; n>0; c++)
50
    n &= n - 1; // clear the least significant bit set
50
    n &= n - 1; // clear the least significant bit set
51
  return c;
51
  return c;
52
}
52
}
53
 
53
 
54
/*
54
/*
55
 * Return the position of the least significant bit that is set
55
 * Return the position of the least significant bit that is set
56
 */
56
 */
57
static int leastSignificantBitSet(int x){
57
static int leastSignificantBitSet(int x){
58
  if(x == 0)
58
  if(x == 0)
59
      return 0;
59
      return 0;
60
 
60
 
61
  int val = 1;
61
  int val = 1;
62
  while(x>>=1)
62
  while(x>>=1)
63
      val++;
63
      val++;
64
 
64
 
65
  return val;
65
  return val;
66
}
66
}
67
 
67
 
68
static inline unsigned int powi(int num, unsigned int exponent){
68
static inline unsigned int powi(int num, unsigned int exponent){
69
 
69
 
70
    if(exponent == 0)
70
    if(exponent == 0)
71
        return 1;
71
        return 1;
72
 
72
 
73
    float res = num;
73
    float res = num;
74
    for(unsigned int i=0; i<exponent-1; i++)
74
    for(unsigned int i=0; i<exponent-1; i++)
75
        res *= num;
75
        res *= num;
76
 
76
 
77
    return res;
77
    return res;
78
}
78
}
79
 
79
 
80
static inline unsigned int twopowi(unsigned int exponent){
80
static inline unsigned int twopowi(unsigned int exponent){
81
 
81
 
82
    return 1 << exponent;
82
    return 1 << exponent;
83
}
83
}
84
 
84
 
85
// Algorithm
85
// Algorithm
86
AlgorithmLineShift::AlgorithmLineShift(unsigned int _screenCols, unsigned int _screenRows) : Algorithm(_screenCols, _screenRows){
86
AlgorithmLineShift::AlgorithmLineShift(unsigned int _screenCols, unsigned int _screenRows) : Algorithm(_screenCols, _screenRows){
87
 
87
 
88
    int nTotalBits = ceilf(log2f((float)screenCols));
88
    int nTotalBits = ceilf(log2f((float)screenCols));
89
 
89
 
-
 
90
    // determine the necessary Gray code bits and add some robustness
90
    nGrayBits = nTotalBits - floorf(log2f((float)nLineShifts));
91
    nGrayBits = nTotalBits - floorf(log2f((float)nLineShifts)) + 2;
91
 
92
 
92
    N = 2 + 2*nGrayBits + nLineShifts;
93
    N = 2 + 2*nGrayBits + nLineShifts;
93
 
94
 
94
    // all on pattern
95
    // all on pattern
95
    cv::Mat allOn(1, screenCols, CV_8UC3, cv::Scalar::all(255));
96
    cv::Mat allOn(1, screenCols, CV_8UC3, cv::Scalar::all(255));
96
    patterns.push_back(allOn);
97
    patterns.push_back(allOn);
97
 
98
 
98
    // all off pattern
99
    // all off pattern
99
    cv::Mat allOff(1, screenCols, CV_8UC3, cv::Scalar::all(0));
100
    cv::Mat allOff(1, screenCols, CV_8UC3, cv::Scalar::all(0));
100
    patterns.push_back(allOff);
101
    patterns.push_back(allOff);
101
 
102
 
102
 
103
 
103
    // Gray code patterns
104
    // Gray code patterns
104
    for(unsigned int p=0; p<nGrayBits; p++){
105
    for(unsigned int p=0; p<nGrayBits; p++){
105
        cv::Mat pattern(1, screenCols, CV_8UC3);
106
        cv::Mat pattern(1, screenCols, CV_8UC3);
106
        cv::Mat patternInv(1, screenCols, CV_8UC3);
107
        cv::Mat patternInv(1, screenCols, CV_8UC3);
107
 
108
 
108
        for(unsigned int j=0; j<screenCols; j++){
109
        for(unsigned int j=0; j<screenCols; j++){
109
 
110
 
110
            unsigned int jGray = binaryToGray(j);
111
            unsigned int jGray = binaryToGray(j);
111
            // Amplitude of channels
112
            // Amplitude of channels
112
            int bit = (int)getBit(jGray, nTotalBits-p);
113
            int bit = (int)getBit(jGray, nTotalBits-p);
113
            pattern.at<cv::Vec3b>(0,j) = cv::Vec3b(255.0*bit,255.0*bit,255.0*bit);
114
            pattern.at<cv::Vec3b>(0,j) = cv::Vec3b(255.0*bit,255.0*bit,255.0*bit);
114
            int invBit = bit^1;
115
            int invBit = bit^1;
115
            patternInv.at<cv::Vec3b>(0,j) = cv::Vec3b(255.0*invBit,255.0*invBit,255.0*invBit);
116
            patternInv.at<cv::Vec3b>(0,j) = cv::Vec3b(255.0*invBit,255.0*invBit,255.0*invBit);
116
        }
117
        }
117
        patterns.push_back(pattern);
118
        patterns.push_back(pattern);
118
        patterns.push_back(patternInv);
119
        patterns.push_back(patternInv);
119
    }
120
    }
120
 
121
 
121
    // line shifts
122
    // line shifts
122
    for(unsigned int p=0; p<nLineShifts; p++){
123
    for(unsigned int p=0; p<nLineShifts; p++){
123
        cv::Mat pattern(1, screenCols, CV_8UC3, cv::Scalar(0));
124
        cv::Mat pattern(1, screenCols, CV_8UC3, cv::Scalar(0));
124
 
125
 
125
        for(unsigned int j=p; j<screenCols; j+= nLineShifts)
126
        for(unsigned int j=p; j<screenCols; j+= nLineShifts)
126
            pattern.at<cv::Vec3b>(0, j) = cv::Vec3b(255, 255, 255);
127
            pattern.at<cv::Vec3b>(0, j) = cv::Vec3b(255, 255, 255);
127
 
128
 
128
        patterns.push_back(pattern);
129
        patterns.push_back(pattern);
129
    }
130
    }
130
 
131
 
131
}
132
}
132
 
133
 
133
cv::Mat AlgorithmLineShift::getEncodingPattern(unsigned int depth){
134
cv::Mat AlgorithmLineShift::getEncodingPattern(unsigned int depth){
134
    return patterns[depth];
135
    return patterns[depth];
135
}
136
}
136
 
137
 
137
static cv::Vec3b getColorSubpix(const cv::Mat& img, cv::Point2f pt){
138
static cv::Vec3b getColorSubpix(const cv::Mat& img, cv::Point2f pt){
138
    assert(!img.empty());
139
    assert(!img.empty());
139
    assert(img.channels() == 3);
140
    assert(img.channels() == 3);
140
 
141
 
141
    int x = (int)pt.x;
142
    int x = (int)pt.x;
142
    int y = (int)pt.y;
143
    int y = (int)pt.y;
143
 
144
 
144
    int x0 = cv::borderInterpolate(x,   img.cols, cv::BORDER_REFLECT_101);
145
    int x0 = cv::borderInterpolate(x,   img.cols, cv::BORDER_REFLECT_101);
145
    int x1 = cv::borderInterpolate(x+1, img.cols, cv::BORDER_REFLECT_101);
146
    int x1 = cv::borderInterpolate(x+1, img.cols, cv::BORDER_REFLECT_101);
146
    int y0 = cv::borderInterpolate(y,   img.rows, cv::BORDER_REFLECT_101);
147
    int y0 = cv::borderInterpolate(y,   img.rows, cv::BORDER_REFLECT_101);
147
    int y1 = cv::borderInterpolate(y+1, img.rows, cv::BORDER_REFLECT_101);
148
    int y1 = cv::borderInterpolate(y+1, img.rows, cv::BORDER_REFLECT_101);
148
 
149
 
149
    float a = pt.x - (float)x;
150
    float a = pt.x - (float)x;
150
    float c = pt.y - (float)y;
151
    float c = pt.y - (float)y;
151
 
152
 
152
    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)
153
    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)
153
                           + (img.at<cv::Vec3b>(y1, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[0] * a) * c);
154
                           + (img.at<cv::Vec3b>(y1, x0)[0] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[0] * a) * c);
154
    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)
155
    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)
155
                           + (img.at<cv::Vec3b>(y1, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[1] * a) * c);
156
                           + (img.at<cv::Vec3b>(y1, x0)[1] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[1] * a) * c);
156
    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)
157
    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)
157
                           + (img.at<cv::Vec3b>(y1, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[2] * a) * c);
158
                           + (img.at<cv::Vec3b>(y1, x0)[2] * (1.f - a) + img.at<cv::Vec3b>(y1, x1)[2] * a) * c);
158
 
159
 
159
    return cv::Vec3b(b, g, r);
160
    return cv::Vec3b(b, g, r);
160
}
161
}
161
 
162
 
162
void getlineCenters(const cv::Mat& linesScanLine, const cv::Mat& codeScanLine, std::vector<cv::Vec2f>& lineCenters){
163
void getlineCenters(const cv::Mat& linesScanLine, const cv::Mat& codeScanLine, std::vector<cv::Vec2f>& lineCenters){
163
 
164
 
164
    int nCols = linesScanLine.cols;
165
    int nCols = linesScanLine.cols;
165
 
166
 
166
    // fourth order derivative filtering
167
    // finite derivative
167
    cv::Mat der(1, nCols, CV_32F);
168
    cv::Mat g(1, nCols, CV_32F);
-
 
169
    cv::Mat h(1, nCols, CV_32F);
168
    for(int i=2; i<nCols-2; i++)
170
    for(int i=2; i<nCols-2; i++){
169
        der.at<float>(0, i) = linesScanLine.at<unsigned char>(0, i-2)+linesScanLine.at<unsigned char>(0, i-1)-
171
//        g.at<float>(0, i) = linesScanLine.at<unsigned char>(0, i+2)+linesScanLine.at<unsigned char>(0, i+1)-
170
                              linesScanLine.at<unsigned char>(0, i+1)-linesScanLine.at<unsigned char>(0, i+2);
172
//                            linesScanLine.at<unsigned char>(0, i-1)-linesScanLine.at<unsigned char>(0, i-2);
-
 
173
        g.at<float>(0, i) = 1.0*linesScanLine.at<unsigned char>(0, i+2)+8.0*linesScanLine.at<unsigned char>(0, i+1)-
-
 
174
                            8.0*linesScanLine.at<unsigned char>(0, i-1)-1.0*linesScanLine.at<unsigned char>(0, i-2);
-
 
175
        h.at<float>(0, i) = -1.0*linesScanLine.at<unsigned char>(0, i+2)+16.0*linesScanLine.at<unsigned char>(0, i+1)-
-
 
176
                            30.0*linesScanLine.at<unsigned char>(0, i)+
-
 
177
                            16.0*linesScanLine.at<unsigned char>(0, i-1)-1.0*linesScanLine.at<unsigned char>(0, i-2);
171
 
178
    }
172
//    cvtools::writeMat(codeScanLine, "codeScanLine.mat", "codeScanLine");
179
//    cvtools::writeMat(codeScanLine, "codeScanLine.mat", "codeScanLine");
173
//    cvtools::writeMat(linesScanLine, "linesScanLine.mat", "linesScanLine");
180
//    cvtools::writeMat(linesScanLine, "linesScanLine.mat", "linesScanLine");
174
//    cvtools::writeMat(der, "der.mat", "der");
181
//    cvtools::writeMat(der, "der.mat", "der");
175
 
182
 
176
    for(int i=0; i<nCols; i++){
183
    for(int i=0; i<nCols-1; i++){
-
 
184
//        float fLeft = linesScanLine.at<unsigned char>(0, i-1);
177
        int code = codeScanLine.at<int>(0, i);
185
        float fI = linesScanLine.at<unsigned char>(0, i);
-
 
186
//        float fRight = linesScanLine.at<unsigned char>(0, i+1);
-
 
187
 
-
 
188
        float gI = g.at<float>(0, i);
-
 
189
        float gRight = g.at<float>(0, i+1);
-
 
190
 
-
 
191
        float hI = h.at<float>(0, i);
178
 
192
 
-
 
193
        int codeI = codeScanLine.at<int>(0, i);
179
        if((code != -1) && (der.at<float>(0, i) < 0.0) && (der.at<float>(0, i+1) > 0.0) && ((der.at<float>(0, i+1) - der.at<float>(0, i)) > 20.0)){
194
        //int codeRight = codeScanLine.at<int>(0, i+1);
180
 
195
 
-
 
196
        if((codeI != -1) && (fI > 1) && (gI >= 0.0) && (gRight <= 0.0) && (gRight < gI) && (hI < -2.0)){
181
            lineCenters.push_back(cv::Vec2f(i, code));
197
            float delta = gI/(gI - gRight);
182
            // TODO: subpixel interpolation, non-max suppression
198
            lineCenters.push_back(cv::Vec2f(i + delta, codeI));
183
        }
199
        }
184
 
200
 
185
    }
201
    }
186
 
202
 
187
}
203
}
188
 
204
 
189
void AlgorithmLineShift::get3DPoints(SMCalibrationParameters calibration, const std::vector<cv::Mat>& frames0, const std::vector<cv::Mat>& frames1, std::vector<cv::Point3f>& Q, std::vector<cv::Vec3b>& color){
205
void AlgorithmLineShift::get3DPoints(SMCalibrationParameters calibration, const std::vector<cv::Mat>& frames0, const std::vector<cv::Mat>& frames1, std::vector<cv::Point3f>& Q, std::vector<cv::Vec3b>& color){
190
 
206
 
191
    assert(frames0.size() == N);
207
    assert(frames0.size() == N);
192
    assert(frames1.size() == N);
208
    assert(frames1.size() == N);
193
 
209
 
194
    int frameRows = frames0[0].rows;
210
    int frameRows = frames0[0].rows;
195
    int frameCols = frames0[0].cols;
211
    int frameCols = frames0[0].cols;
196
 
212
 
197
    // rectifying homographies (rotation+projections)
213
    // rectifying homographies (rotation+projections)
198
    cv::Size frameSize(frameCols, frameRows);
214
    cv::Size frameSize(frameCols, frameRows);
199
    cv::Mat R, T;
215
    cv::Mat R, T;
200
    // stereoRectify segfaults unless R is double precision
216
    // stereoRectify segfaults unless R is double precision
201
    cv::Mat(calibration.R1).convertTo(R, CV_64F);
217
    cv::Mat(calibration.R1).convertTo(R, CV_64F);
202
    cv::Mat(calibration.T1).convertTo(T, CV_64F);
218
    cv::Mat(calibration.T1).convertTo(T, CV_64F);
203
    cv::Mat R0, R1, P0, P1, QRect;
219
    cv::Mat R0, R1, P0, P1, QRect;
204
    cv::stereoRectify(calibration.K0, calibration.k0, calibration.K1, calibration.k1, frameSize, R, T, R0, R1, P0, P1, QRect, 0);
220
    cv::stereoRectify(calibration.K0, calibration.k0, calibration.K1, calibration.k1, frameSize, R, T, R0, R1, P0, P1, QRect, 0);
205
 
221
 
206
    // interpolation maps
222
    // interpolation maps
207
    cv::Mat map0X, map0Y, map1X, map1Y;
223
    cv::Mat map0X, map0Y, map1X, map1Y;
208
    cv::initUndistortRectifyMap(calibration.K0, calibration.k0, R0, P0, frameSize, CV_32F, map0X, map0Y);
224
    cv::initUndistortRectifyMap(calibration.K0, calibration.k0, R0, P0, frameSize, CV_32F, map0X, map0Y);
209
    cv::initUndistortRectifyMap(calibration.K1, calibration.k1, R1, P1, frameSize, CV_32F, map1X, map1Y);
225
    cv::initUndistortRectifyMap(calibration.K1, calibration.k1, R1, P1, frameSize, CV_32F, map1X, map1Y);
210
 
226
 
211
    // gray-scale and remap
227
    // gray-scale and remap
212
    std::vector<cv::Mat> frames0Rect(N);
228
    std::vector<cv::Mat> frames0Rect(N);
213
    std::vector<cv::Mat> frames1Rect(N);
229
    std::vector<cv::Mat> frames1Rect(N);
214
    for(int i=0; i<N; i++){
230
    for(int i=0; i<N; i++){
215
        cv::Mat temp;
231
        cv::Mat temp;
216
        cv::cvtColor(frames0[i], temp, CV_BayerBG2GRAY);
232
        cv::cvtColor(frames0[i], temp, CV_BayerBG2GRAY);
217
        cv::remap(temp, frames0Rect[i], map0X, map0Y, CV_INTER_LINEAR);
233
        cv::remap(temp, frames0Rect[i], map0X, map0Y, CV_INTER_LINEAR);
218
        cv::cvtColor(frames1[i], temp, CV_BayerBG2GRAY);
234
        cv::cvtColor(frames1[i], temp, CV_BayerBG2GRAY);
219
        cv::remap(temp, frames1Rect[i], map1X, map1Y, CV_INTER_LINEAR);
235
        cv::remap(temp, frames1Rect[i], map1X, map1Y, CV_INTER_LINEAR);
220
    }
236
    }
221
 
237
 
222
    //cvtools::writeMat(frames0Rect[0], "frames0Rect_0.mat", "frames0Rect_0");
238
    //cvtools::writeMat(frames0Rect[0], "frames0Rect_0.mat", "frames0Rect_0");
223
    //cvtools::writeMat(frames0Rect[1], "frames0Rect_1.mat", "frames0Rect_1");
239
    //cvtools::writeMat(frames0Rect[1], "frames0Rect_1.mat", "frames0Rect_1");
224
    //cvtools::writeMat(frames0Rect[20], "frames0Rect_20.mat", "frames0Rect_20");
240
    //cvtools::writeMat(frames0Rect[20], "frames0Rect_20.mat", "frames0Rect_20");
225
    //cvtools::writeMat(frames0Rect[21], "frames0Rect_21.mat", "frames0Rect_21");
241
    //cvtools::writeMat(frames0Rect[21], "frames0Rect_21.mat", "frames0Rect_21");
226
 
242
 
227
    // color debayer and remap
243
    // color debayer and remap
228
    cv::Mat color0Rect, color1Rect;
244
    cv::Mat color0Rect, color1Rect;
229
    cv::cvtColor(frames0[0], color0Rect, CV_BayerBG2RGB);
245
    cv::cvtColor(frames0[0], color0Rect, CV_BayerBG2RGB);
230
    cv::remap(color0Rect, color0Rect, map0X, map0Y, CV_INTER_LINEAR);
246
    cv::remap(color0Rect, color0Rect, map0X, map0Y, CV_INTER_LINEAR);
231
 
247
 
232
    cv::cvtColor(frames1[0], color1Rect, CV_BayerBG2RGB);
248
    cv::cvtColor(frames1[0], color1Rect, CV_BayerBG2RGB);
233
    cv::remap(color1Rect, color1Rect, map1X, map1Y, CV_INTER_LINEAR);
249
    cv::remap(color1Rect, color1Rect, map1X, map1Y, CV_INTER_LINEAR);
234
 
250
 
235
    int frameRectRows = frames0Rect[0].rows;
251
    int frameRectRows = frames0Rect[0].rows;
236
    int frameRectCols = frames0Rect[0].cols;
252
    int frameRectCols = frames0Rect[0].cols;
237
 
253
 
238
    // occlusion masks
254
    // occlusion masks
239
    cv::Mat occlusion0Rect, occlusion1Rect;
255
    cv::Mat occlusion0Rect, occlusion1Rect;
240
    cv::subtract(frames0Rect[0], frames0Rect[1], occlusion0Rect);
256
    cv::subtract(frames0Rect[0], frames0Rect[1], occlusion0Rect);
241
    occlusion0Rect = (occlusion0Rect > 20) & (occlusion0Rect < 250);
257
    occlusion0Rect = (occlusion0Rect > 20) & (occlusion0Rect < 250);
242
    cv::subtract(frames1Rect[0], frames1Rect[1], occlusion1Rect);
258
    cv::subtract(frames1Rect[0], frames1Rect[1], occlusion1Rect);
243
    occlusion1Rect = (occlusion1Rect > 20) & (occlusion1Rect < 250);
259
    occlusion1Rect = (occlusion1Rect > 20) & (occlusion1Rect < 250);
244
 
260
 
245
//    cvtools::writeMat(occlusion0Rect, "occlusion0Rect.mat", "occlusion0Rect");
261
//    cvtools::writeMat(occlusion0Rect, "occlusion0Rect.mat", "occlusion0Rect");
246
//    cvtools::writeMat(occlusion1Rect, "occlusion1Rect.mat", "occlusion1Rect");
262
//    cvtools::writeMat(occlusion1Rect, "occlusion1Rect.mat", "occlusion1Rect");
247
 
263
 
248
    // erode occlusion masks
264
    // erode occlusion masks
249
    cv::Mat strel = cv::getStructuringElement(cv::MORPH_ELLIPSE, cv::Size(2,2));
265
    cv::Mat strel = cv::getStructuringElement(cv::MORPH_ELLIPSE, cv::Size(2,2));
250
    cv::erode(occlusion0Rect, occlusion0Rect, strel);
266
    cv::erode(occlusion0Rect, occlusion0Rect, strel);
251
    cv::erode(occlusion1Rect, occlusion1Rect, strel);
267
    cv::erode(occlusion1Rect, occlusion1Rect, strel);
252
 
268
 
253
    // correct for texture modulation and ambient
269
//    // correct for texture modulation and ambient
254
    cv::Mat A0 = frames0Rect[1];
270
//    cv::Mat A0 = frames0Rect[1];
255
    cv::Mat M0 = frames0Rect[0]-frames0Rect[1];
271
//    cv::Mat M0 = frames0Rect[0]-frames0Rect[1];
256
 
272
 
257
    cv::divide(256.0, M0, M0, CV_32F);
273
//    cv::divide(256.0, M0, M0, CV_32F);
258
    cv::Mat A1 = frames1Rect[1];
274
//    cv::Mat A1 = frames1Rect[1];
259
    cv::Mat M1 = frames1Rect[0]-frames1Rect[1];
275
//    cv::Mat M1 = frames1Rect[0]-frames1Rect[1];
260
    cv::divide(256.0, M1, M1, CV_32F);
276
//    cv::divide(256.0, M1, M1, CV_32F);
261
 
277
 
262
    for(int i=2; i<N; i++){
278
//    for(int i=2; i<N; i++){
263
        cv::multiply(frames0Rect[i]-A0, M0, frames0Rect[i], 1.0, CV_8UC1);
279
//        cv::multiply(frames0Rect[i]-A0, M0, frames0Rect[i], 1.0, CV_8UC1);
264
        cv::multiply(frames1Rect[i]-A1, M1, frames1Rect[i], 1.0, CV_8UC1);
280
//        cv::multiply(frames1Rect[i]-A1, M1, frames1Rect[i], 1.0, CV_8UC1);
265
    }
281
//    }
266
 
282
 
267
    //cvtools::writeMat(frames0Rect[22], "frames0Rect_22.mat", "frames0Rect_22");
283
    //cvtools::writeMat(frames0Rect[22], "frames0Rect_22.mat", "frames0Rect_22");
268
    //cvtools::writeMat(frames0Rect[23], "frames0Rect_23.mat", "frames0Rect_23");
284
    //cvtools::writeMat(frames0Rect[23], "frames0Rect_23.mat", "frames0Rect_23");
269
 
285
 
270
    // divide into Gray coding frames and line shift frames
286
    // divide into Gray coding frames and line shift frames
271
    std::vector<cv::Mat> frames0GrayCode(frames0Rect.begin()+2, frames0Rect.begin()+2+2*nGrayBits);
287
    std::vector<cv::Mat> frames0GrayCode(frames0Rect.begin()+2, frames0Rect.begin()+2+2*nGrayBits);
272
    std::vector<cv::Mat> frames0LineShift(frames0Rect.begin()+2+2*nGrayBits, frames0Rect.end());
288
    std::vector<cv::Mat> frames0LineShift(frames0Rect.begin()+2+2*nGrayBits, frames0Rect.end());
273
    std::vector<cv::Mat> frames1GrayCode(frames1Rect.begin()+2, frames1Rect.begin()+2+2*nGrayBits);
289
    std::vector<cv::Mat> frames1GrayCode(frames1Rect.begin()+2, frames1Rect.begin()+2+2*nGrayBits);
274
    std::vector<cv::Mat> frames1LineShift(frames1Rect.begin()+2+2*nGrayBits, frames1Rect.end());
290
    std::vector<cv::Mat> frames1LineShift(frames1Rect.begin()+2+2*nGrayBits, frames1Rect.end());
275
 
291
 
276
    // decode patterns
292
    // decode patterns
277
    cv::Mat code0Gray(frameRectRows, frameRectCols, CV_32S, cv::Scalar(0));
293
    cv::Mat code0Gray(frameRectRows, frameRectCols, CV_32S, cv::Scalar(0));
278
    cv::Mat code1Gray(frameRectRows, frameRectCols, CV_32S, cv::Scalar(0));
294
    cv::Mat code1Gray(frameRectRows, frameRectCols, CV_32S, cv::Scalar(0));
279
 
295
 
280
    // into gray code
296
    // into gray code
281
    for(int i=0; i<nGrayBits; i++){
297
    for(int i=0; i<nGrayBits; i++){
282
        cv::Mat temp, bit0, bit1;
298
        cv::Mat temp, bit0, bit1;
283
 
299
 
284
        cv::compare(frames0GrayCode[i*2], frames0GrayCode[i*2+1], temp, cv::CMP_GT);
300
        cv::compare(frames0GrayCode[i*2], frames0GrayCode[i*2+1], temp, cv::CMP_GT);
285
        temp.convertTo(bit0, CV_32S, 1.0/255.0);
301
        temp.convertTo(bit0, CV_32S, 1.0/255.0);
286
        cv::add(code0Gray, bit0*twopowi(nGrayBits-i-1), code0Gray, cv::noArray(), CV_32S);
302
        cv::add(code0Gray, bit0*twopowi(nGrayBits-i-1), code0Gray, cv::noArray(), CV_32S);
287
 
303
 
288
        cv::compare(frames1GrayCode[i*2], frames1GrayCode[i*2+1], temp, cv::CMP_GT);
304
        cv::compare(frames1GrayCode[i*2], frames1GrayCode[i*2+1], temp, cv::CMP_GT);
289
        temp.convertTo(bit1, CV_32S, 1.0/255.0);
305
        temp.convertTo(bit1, CV_32S, 1.0/255.0);
290
        cv::add(code1Gray, bit1*twopowi(nGrayBits-i-1), code1Gray, cv::noArray(), CV_32S);
306
        cv::add(code1Gray, bit1*twopowi(nGrayBits-i-1), code1Gray, cv::noArray(), CV_32S);
291
    }
307
    }
292
 
308
 
293
    // convert to standard binary
309
    // convert to standard binary
294
    cv::Mat code0Binary(code0Gray.rows, code0Gray.cols, CV_32S, cv::Scalar(-1));
310
    cv::Mat code0Binary(code0Gray.rows, code0Gray.cols, CV_32S, cv::Scalar(-1));
295
    cv::Mat code1Binary(code1Gray.rows, code1Gray.cols, CV_32S, cv::Scalar(-1));
311
    cv::Mat code1Binary(code1Gray.rows, code1Gray.cols, CV_32S, cv::Scalar(-1));
296
    for(int r=0; r<frameRectRows; r++){
312
    for(int r=0; r<frameRectRows; r++){
297
        for(int c=0; c<frameRectCols; c++){
313
        for(int c=0; c<frameRectCols; c++){
298
            code0Binary.at<int>(r,c) = grayToBinary(code0Gray.at<int>(r,c));
314
            code0Binary.at<int>(r,c) = grayToBinary(code0Gray.at<int>(r,c));
299
            code1Binary.at<int>(r,c) = grayToBinary(code1Gray.at<int>(r,c));
315
            code1Binary.at<int>(r,c) = grayToBinary(code1Gray.at<int>(r,c));
300
        }
316
        }
301
    }
317
    }
302
 
318
 
303
    // set occluded pixels to -1
319
    // set occluded pixels to -1
304
    for(int r=0; r<frameRectRows; r++){
320
    for(int r=0; r<frameRectRows; r++){
305
        for(int c=0; c<frameRectCols; c++){
321
        for(int c=0; c<frameRectCols; c++){
306
            if(occlusion0Rect.at<unsigned char>(r,c) == 0)
322
            if(occlusion0Rect.at<unsigned char>(r,c) == 0)
307
                code0Binary.at<int>(r,c) = -1;
323
                code0Binary.at<int>(r,c) = -1;
308
            if(occlusion1Rect.at<unsigned char>(r,c) == 0)
324
            if(occlusion1Rect.at<unsigned char>(r,c) == 0)
309
                code1Binary.at<int>(r,c) = -1;
325
                code1Binary.at<int>(r,c) = -1;
310
        }
326
        }
311
    }
327
    }
312
 
328
 
313
//cvtools::writeMat(code0Gray, "code0Gray.mat", "code0Gray");
329
//cvtools::writeMat(code0Gray, "code0Gray.mat", "code0Gray");
314
//cvtools::writeMat(code1Gray, "code1Gray.mat", "code1Gray");
330
//cvtools::writeMat(code1Gray, "code1Gray.mat", "code1Gray");
315
//cvtools::writeMat(code0Binary, "code0Binary.mat", "code0Binary");
331
//cvtools::writeMat(code0Binary, "code0Binary.mat", "code0Binary");
316
//cvtools::writeMat(code1Binary, "code1Binary.mat", "code1Binary");
332
//cvtools::writeMat(code1Binary, "code1Binary.mat", "code1Binary");
317
 
333
 
318
    // matching
334
    // matching
319
    std::vector<cv::Vec2f> q0Rect, q1Rect;
335
    std::vector<cv::Vec2f> q0Rect, q1Rect;
320
    for(int s=0; s<nLineShifts; s++){
336
    for(int s=0; s<nLineShifts; s++){
321
 
337
 
322
        cv::Mat lines0 = frames0LineShift[s];
338
        cv::Mat lines0 = frames0LineShift[s];
323
        cv::Mat lines1 = frames1LineShift[s];
339
        cv::Mat lines1 = frames1LineShift[s];
324
 
340
 
325
        for(int row=0; row<frameRectRows; row++){
341
        for(int row=0; row<frameRectRows; row++){
326
 
342
 
327
            // line center data structure containing [x-coordinate (sub-px), region-code]
343
            // line center data structure containing [x-coordinate (sub-px), region-code]
328
            std::vector<cv::Vec2f> lineCenters0, lineCenters1;
344
            std::vector<cv::Vec2f> lineCenters0, lineCenters1;
329
 
345
 
330
            // sorted, unique line centers
346
            // sorted, unique line centers
331
            getlineCenters(lines0.row(row), code0Binary.row(row), lineCenters0);
347
            getlineCenters(lines0.row(row), code0Binary.row(row), lineCenters0);
332
            getlineCenters(lines1.row(row), code1Binary.row(row), lineCenters1);
348
            getlineCenters(lines1.row(row), code1Binary.row(row), lineCenters1);
333
 
349
 
334
//         if(s==0 && row==1300){
350
         if(s==0 && row==1300){
335
//            std::cout << cv::Mat(lineCenters0) << std::endl;
351
            std::cout << cv::Mat(lineCenters0) << std::endl;
336
//            std::cout << cv::Mat(lineCenters1) << std::endl;
352
            std::cout << cv::Mat(lineCenters1) << std::endl;
337
 
353
 
338
//            cvtools::writeMat(lines0.row(row), "lines0.mat", "lines0");
354
            cvtools::writeMat(lines0.row(row), "lines0.mat", "lines0");
339
//            cvtools::writeMat(lines1.row(row), "lines1.mat", "lines1");
355
            cvtools::writeMat(lines1.row(row), "lines1.mat", "lines1");
340
//            cvtools::writeMat(code0Binary.row(row), "code0Binary.mat", "code0Binary");
356
            cvtools::writeMat(code0Binary.row(row), "code0Binary.mat", "code0Binary");
341
//            cvtools::writeMat(code1Binary.row(row), "code1Binary.mat", "code1Binary");
357
            cvtools::writeMat(code1Binary.row(row), "code1Binary.mat", "code1Binary");
342
//         }
358
         }
343
 
359
 
344
            // match and store
360
            // match and store
345
            int i=0, j=0;
361
            int i=0, j=0;
346
            while(i<lineCenters0.size() && j<lineCenters1.size()){
362
            while(i<lineCenters0.size() && j<lineCenters1.size()){
347
 
363
 
348
                if(lineCenters0[i][1] == lineCenters1[j][1]){
364
                if(lineCenters0[i][1] == lineCenters1[j][1]){
349
                    q0Rect.push_back(cv::Point2f(lineCenters0[i][0], row));
365
                    q0Rect.push_back(cv::Point2f(lineCenters0[i][0], row));
350
                    q1Rect.push_back(cv::Point2f(lineCenters1[j][0], row));
366
                    q1Rect.push_back(cv::Point2f(lineCenters1[j][0], row));
351
                    i += 1;
367
                    i += 1;
352
                    j += 1;
368
                    j += 1;
353
                } else if(lineCenters0[i][1] < lineCenters1[j][1]){
369
                } else if(lineCenters0[i][1] < lineCenters1[j][1]){
354
                    i += 1;
370
                    i += 1;
355
                } else if(lineCenters0[i][1] > lineCenters1[j][1]){
371
                } else if(lineCenters0[i][1] > lineCenters1[j][1]){
356
                    j += 1;
372
                    j += 1;
357
                }
373
                }
358
            }
374
            }
359
 
375
 
360
        }
376
        }
361
    }
377
    }
362
 
378
 
363
    int nMatches = q0Rect.size();
379
    int nMatches = q0Rect.size();
364
 
380
 
365
    if(nMatches < 1){
381
    if(nMatches < 1){
366
        Q.resize(0);
382
        Q.resize(0);
367
        color.resize(0);
383
        color.resize(0);
368
 
384
 
369
        return;
385
        return;
370
    }
386
    }
371
 
387
 
372
    // retrieve color information (at integer coordinates)
388
    // retrieve color information (at integer coordinates)
373
    color.resize(nMatches);
389
    color.resize(nMatches);
374
    for(int i=0; i<nMatches; i++){
390
    for(int i=0; i<nMatches; i++){
375
 
391
 
376
        cv::Vec3b c0 = color0Rect.at<cv::Vec3b>(q0Rect[i][1], q0Rect[i][0]);
392
        cv::Vec3b c0 = color0Rect.at<cv::Vec3b>(q0Rect[i][1], q0Rect[i][0]);
377
        cv::Vec3b c1 = color1Rect.at<cv::Vec3b>(q1Rect[i][1], q1Rect[i][0]);
393
        cv::Vec3b c1 = color1Rect.at<cv::Vec3b>(q1Rect[i][1], q1Rect[i][0]);
378
//        cv::Vec3b c0 = getColorSubpix(color0Rect, q0Rect[i]);
394
//        cv::Vec3b c0 = getColorSubpix(color0Rect, q0Rect[i]);
379
//        cv::Vec3b c1 = getColorSubpix(color1Rect, q0Rect[i]);
395
//        cv::Vec3b c1 = getColorSubpix(color1Rect, q0Rect[i]);
380
 
396
 
381
        color[i] = 0.5*c0 + 0.5*c1;
397
        color[i] = 0.5*c0 + 0.5*c1;
382
    }
398
    }
383
 
399
 
384
    // triangulate points
400
    // triangulate points
385
    cv::Mat QMatHomogenous, QMat;
401
    cv::Mat QMatHomogenous, QMat;
386
 
402
 
387
    cv::triangulatePoints(P0, P1, q0Rect, q1Rect, QMatHomogenous);
403
    cv::triangulatePoints(P0, P1, q0Rect, q1Rect, QMatHomogenous);
388
    cvtools::convertMatFromHomogeneous(QMatHomogenous, QMat);
404
    cvtools::convertMatFromHomogeneous(QMatHomogenous, QMat);
389
 
405
 
390
    // undo rectifying rotation
406
    // undo rectifying rotation
391
    cv::Mat R0Inv;
407
    cv::Mat R0Inv;
392
    cv::Mat(R0.t()).convertTo(R0Inv, CV_32F);
408
    cv::Mat(R0.t()).convertTo(R0Inv, CV_32F);
393
    QMat = R0Inv*QMat;
409
    QMat = R0Inv*QMat;
394
 
410
 
395
    cvtools::matToPoints3f(QMat, Q);
411
    cvtools::matToPoints3f(QMat, Q);
396
 
412
 
397
}
413
}
398
 
414