Walsh-Hadamard transform and clearing of photomontage part

pull/49/head
Bellaktris 11 years ago
parent 1aad484ccf
commit b6f379350d
  1. 63
      modules/xphoto/src/algo.hpp
  2. 125
      modules/xphoto/src/annf.hpp
  3. 64
      modules/xphoto/src/dct_image_denoising.cpp
  4. 385
      modules/xphoto/src/gcgraph.hpp
  5. 254
      modules/xphoto/src/inpainting.cpp
  6. 52
      modules/xphoto/src/norm2.hpp
  7. 259
      modules/xphoto/src/photomontage.hpp

@ -0,0 +1,63 @@
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
// Copyright (C) 2009-2011, Willow Garage Inc., all rights reserved.
// Third party copyrights are property of their respective owners.
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of Intel Corporation may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __ALGO_HPP__
#define __ALGO_HPP__
#include <vector>
#include <algorithm>
#include <cmath>
template <typename Tp> static inline int min_idx(std::vector <Tp> vec)
{
return std::min_element(vec.begin(), vec.end()) - vec.begin();
}
static inline int hamming_length(int x)
{
int res = 0;
while (x)
{
res += x&1;
x >>= 1;
}
return res;
}
#endif /* __ALGO_HPP__ */

@ -0,0 +1,125 @@
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
// Copyright (C) 2009-2011, Willow Garage Inc., all rights reserved.
// Third party copyrights are property of their respective owners.
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of Intel Corporation may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __ANNF_HPP__
#define __ANNF_HPP__
#include "algo.hpp"
static void plusToMinusUpdate(const cv::Mat &current, cv::Mat &next, const int dx, const int dy)
{
for (int i = 0; i < next.rows; ++i)
for (int j = 0; j < next.cols; ++j)
{
int y = cv::borderInterpolate(i - dy, next.rows, cv::BORDER_CONSTANT);
int x = cv::borderInterpolate(j - dx, next.cols, cv::BORDER_CONSTANT);
next.at<float>(i, j) = -next.at<float>(y, x)
+ current.at<float>(i, j) - current.at<float>(y, x);
}
}
static void minusToPlusUpdate(const cv::Mat &current, cv::Mat &next, const int dx, const int dy)
{
for (int i = 0; i < next.rows; ++i)
for (int j = 0; j < next.cols; ++j)
{
int y = cv::borderInterpolate(i - dy, next.rows, cv::BORDER_CONSTANT);
int x = cv::borderInterpolate(j - dx, next.cols, cv::BORDER_CONSTANT);
next.at<float>(i, j) = next.at<float>(y, x)
- current.at<float>(i, j) + current.at<float>(y, x);
}
}
static void getWHSeries(const cv::Mat &src, cv::Mat &dst, const int nProjections, const int psize = 8)
{
CV_Assert(nProjections <= psize*psize && src.type() == CV_32FC3);
CV_Assert( hamming_length(psize) == 1 );
std::vector <cv::Mat> projections;
cv::Mat proj;
cv::boxFilter(proj, proj, CV_32F, cv::Size(psize, psize),
cv::Point(-1,-1), false, cv::BORDER_REFLECT);
projections.push_back(proj);
std::vector <int> snake_idx(1, 0);
std::vector <int> snake_idy(1, 0);
for (int k = 1, num = 1; k < psize && num <= nProjections; ++k)
{
int dx[] = { (k % 2 == 0) ? +1 : 0, (k % 2 == 0) ? 0 : -1};
int dy[] = { (k % 2 == 0) ? 0 : +1, (k % 2 == 0) ? -1 : 0};
snake_idx.push_back(snake_idx[num - 1] - dx[1]);
snake_idy.push_back(snake_idy[num++ - 1] - dy[1]);
for (int i = 0; i < k && num < nProjections; ++i, ++num)
{
snake_idx.push_back(snake_idx[num - 1] + dx[0]);
snake_idy.push_back(snake_idy[num - 1] + dy[0]);
}
for (int i = 0; i < k && num < nProjections; ++i, ++num)
{
snake_idx.push_back(snake_idx[num - 1] + dx[1]);
snake_idy.push_back(snake_idy[num - 1] + dy[1]);
}
}
for (int i = 1; i < nProjections; ++i)
{
int dx = (snake_idx[i] - snake_idx[i - 1]);
int dy = (snake_idy[i] - snake_idy[i - 1]);
dx <<= hamming_length(psize - 1) - hamming_length(snake_idx[i - 1] ^ snake_idx[i]);
dy <<= hamming_length(psize - 1) - hamming_length(snake_idy[i - 1] ^ snake_idy[i]);
if (i % 2 == 0)
plusToMinusUpdate(proj, proj, dx, dy);
else
minusToPlusUpdate(proj, proj, dx, dy);
}
cv::merge(projections, dst);
}
#endif /* __ANNF_HPP__ */

@ -57,21 +57,21 @@ namespace cv
struct grayDctDenoisingInvoker : public ParallelLoopBody
{
public:
grayDctDenoisingInvoker(const Mat_<float> &src, std::vector < Mat_<float> > &patches, const double sigma, const int psize);
grayDctDenoisingInvoker(const Mat &src, std::vector <Mat> &patches, const double sigma, const int psize);
~grayDctDenoisingInvoker();
void operator() (const Range &range) const;
private:
const Mat_<float> &src;
std::vector < Mat_<float> > &patches; // image decomposition into sliding patches
const Mat &src;
std::vector <Mat> &patches; // image decomposition into sliding patches
const int psize; // size of block to compute dct
const double sigma; // expected noise standard deviation
const double thresh; // thresholding estimate
};
grayDctDenoisingInvoker::grayDctDenoisingInvoker(const Mat_<float> &src, std::vector < Mat_<float> > &patches,
grayDctDenoisingInvoker::grayDctDenoisingInvoker(const Mat &src, std::vector <Mat> &patches,
const double sigma, const int psize)
: src(src), patches(patches), sigma(sigma), thresh(3*sigma), psize(psize) {}
grayDctDenoisingInvoker::~grayDctDenoisingInvoker(){}
@ -85,7 +85,7 @@ namespace cv
Rect patchNum( x, y, psize, psize );
Mat_<float> patch(psize, psize);
Mat patch(psize, psize, CV_32FC1);
src(patchNum).copyTo( patch );
dct(patch, patch);
@ -96,44 +96,20 @@ namespace cv
}
}
void grayDctDenoising(const Mat_<float> &src, Mat_<float> &dst, const double sigma, const int psize)
void grayDctDenoising(const Mat &src, Mat &dst, const double sigma, const int psize)
{
CV_Assert( src.channels() == 1 );
//Mat_<float> res( src.size(), 0.0f ),
// num( src.size(), 0.0f );
//
//double threshold = 3*sigma;
//
//for (int i = 0; i <= src.rows - psize; ++i)
// for (int j = 0; j <= src.cols - psize; ++j)
// {
// Mat_<float> patch = src( Rect(j, i, psize, psize) ).clone();
//
// dct(patch, patch);
// float * ptr = (float *) patch.data;
// for (int k = 0; k < psize*psize; ++k)
// if (fabs(ptr[k]) < threshold)
// ptr[k] = 0.0f;
// idct(patch, patch);
//
// res( Rect(j, i, psize, psize) ) += patch;
// num( Rect(j, i, psize, psize) ) += Mat_<float>::ones(psize, psize);
// }
//res /= num;
//
//res.convertTo( dst, src.type() );
CV_Assert( src.type() == CV_MAKE_TYPE(CV_32F, 1) );
int npixels = (src.rows - psize)*(src.cols - psize);
std::vector < Mat_<float> > patches;
std::vector <Mat> patches;
for (int i = 0; i < npixels; ++i)
patches.push_back( Mat_<float>(psize, psize) );
patches.push_back( Mat(psize, psize, CV_32FC1) );
parallel_for_( cv::Range(0, npixels),
grayDctDenoisingInvoker(src, patches, sigma, psize) );
Mat_<float> res( src.size(), 0.0f ),
num( src.size(), 0.0f );
Mat res( src.size(), CV_32FC1, 0.0f ),
num( src.size(), CV_32FC1, 0.0f );
for (int k = 0; k < npixels; ++k)
{
@ -141,24 +117,24 @@ namespace cv
int j = k % (src.cols - psize);
res( Rect(j, i, psize, psize) ) += patches[k];
num( Rect(j, i, psize, psize) ) += Mat_<float>::ones(psize, psize);
num( Rect(j, i, psize, psize) ) += Mat::ones(psize, psize, CV_32FC1);
}
res /= num;
res.convertTo( dst, src.type() );
}
void rgbDctDenoising(const Mat_<Vec3f> &src, Mat_<Vec3f> &dst, const double sigma, const int psize)
void rgbDctDenoising(const Mat &src, Mat &dst, const double sigma, const int psize)
{
CV_Assert( src.channels() == 3 );
CV_Assert( src.type() == CV_MAKE_TYPE(CV_32F, 3) );
float M[] = {cvInvSqrt(3), cvInvSqrt(3), cvInvSqrt(3),
cvInvSqrt(2), 0.0f, -cvInvSqrt(2),
cvInvSqrt(6), -2.0f*cvInvSqrt(6), cvInvSqrt(6)};
Mat_<Vec3f>::iterator outIt = dst.begin();
Mat_<Vec3f>::iterator outIt = dst.begin<Vec3f>();
for (Mat_<Vec3f>::const_iterator it = src.begin(); it != src.end(); ++it, ++outIt)
for (Mat_<Vec3f>::const_iterator it = src.begin<Vec3f>(); it != src.end<Vec3f>(); ++it, ++outIt)
{
Vec3f rgb = *it;
*outIt = Vec3f(M[0]*rgb[0] + M[1]*rgb[1] + M[2]*rgb[2],
@ -167,7 +143,7 @@ namespace cv
}
/*************************************/
std::vector < Mat_<float> > mv;
std::vector <Mat> mv;
split(dst, mv);
for (int i = 0; i < mv.size(); ++i)
@ -176,7 +152,7 @@ namespace cv
merge(mv, dst);
/*************************************/
for (Mat_<Vec3f>::iterator it = dst.begin(); it != dst.end(); ++it)
for (Mat_<Vec3f>::iterator it = dst.begin<Vec3f>(); it != dst.end<Vec3f>(); ++it)
{
Vec3f rgb = *it;
*it = Vec3f(M[0]*rgb[0] + M[3]*rgb[1] + M[6]*rgb[2],
@ -202,9 +178,9 @@ namespace cv
src.convertTo(img, xtype);
if ( img.type() == CV_32FC3 )
rgbDctDenoising( Mat_<Vec3f>(img), Mat_<Vec3f>(img), sigma, psize );
rgbDctDenoising( img, img, sigma, psize );
else if ( img.type() == CV_32FC1 )
grayDctDenoising( Mat_<float>(img), Mat_<float>(img), sigma, psize );
grayDctDenoising( img, img, sigma, psize );
else
CV_Assert( false );

@ -0,0 +1,385 @@
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// Intel License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2000, Intel Corporation, all rights reserved.
// Third party copyrights are property of their respective owners.
//
// Redistribution and use in source and binary forms, with or without modification,
// are permitted provided that the following conditions are met:
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of Intel Corporation may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef _CV_GCGRAPH_H_
#define _CV_GCGRAPH_H_
template <class TWeight> class GCGraph
{
public:
GCGraph();
GCGraph( unsigned int vtxCount, unsigned int edgeCount );
~GCGraph();
void create( unsigned int vtxCount, unsigned int edgeCount );
int addVtx();
void addEdges( int i, int j, TWeight w, TWeight revw );
void addTermWeights( int i, TWeight sourceW, TWeight sinkW );
TWeight maxFlow();
bool inSourceSegment( int i );
private:
class Vtx
{
public:
Vtx *next; // initialized and used in maxFlow() only
int parent;
int first;
int ts;
int dist;
TWeight weight;
unsigned char t;
};
class Edge
{
public:
int dst;
int next;
TWeight weight;
};
std::vector<Vtx> vtcs;
std::vector<Edge> edges;
TWeight flow;
};
template <class TWeight>
GCGraph<TWeight>::GCGraph()
{
flow = 0;
}
template <class TWeight>
GCGraph<TWeight>::GCGraph( unsigned int vtxCount, unsigned int edgeCount )
{
create( vtxCount, edgeCount );
}
template <class TWeight>
GCGraph<TWeight>::~GCGraph()
{
}
template <class TWeight>
void GCGraph<TWeight>::create( unsigned int vtxCount, unsigned int edgeCount )
{
vtcs.reserve( vtxCount );
edges.reserve( edgeCount + 2 );
flow = 0;
}
template <class TWeight>
int GCGraph<TWeight>::addVtx()
{
Vtx v;
memset( &v, 0, sizeof(Vtx));
vtcs.push_back(v);
return (int)vtcs.size() - 1;
}
template <class TWeight>
void GCGraph<TWeight>::addEdges( int i, int j, TWeight w, TWeight revw )
{
CV_Assert( i>=0 && i<(int)vtcs.size() );
CV_Assert( j>=0 && j<(int)vtcs.size() );
CV_Assert( w>=0 && revw>=0 );
CV_Assert( i != j );
if( !edges.size() )
edges.resize( 2 );
Edge fromI, toI;
fromI.dst = j;
fromI.next = vtcs[i].first;
fromI.weight = w;
vtcs[i].first = (int)edges.size();
edges.push_back( fromI );
toI.dst = i;
toI.next = vtcs[j].first;
toI.weight = revw;
vtcs[j].first = (int)edges.size();
edges.push_back( toI );
}
template <class TWeight>
void GCGraph<TWeight>::addTermWeights( int i, TWeight sourceW, TWeight sinkW )
{
CV_Assert( i>=0 && i<(int)vtcs.size() );
TWeight dw = vtcs[i].weight;
if( dw > 0 )
sourceW += dw;
else
sinkW -= dw;
flow += (sourceW < sinkW) ? sourceW : sinkW;
vtcs[i].weight = sourceW - sinkW;
}
template <class TWeight>
TWeight GCGraph<TWeight>::maxFlow()
{
const int TERMINAL = -1, ORPHAN = -2;
Vtx stub, *nilNode = &stub, *first = nilNode, *last = nilNode;
int curr_ts = 0;
stub.next = nilNode;
Vtx *vtxPtr = &vtcs[0];
Edge *edgePtr = &edges[0];
std::vector<Vtx*> orphans;
// initialize the active queue and the graph vertices
for( int i = 0; i < (int)vtcs.size(); i++ )
{
Vtx* v = vtxPtr + i;
v->ts = 0;
if( v->weight != 0 )
{
last = last->next = v;
v->dist = 1;
v->parent = TERMINAL;
v->t = v->weight < 0;
}
else
v->parent = 0;
}
first = first->next;
last->next = nilNode;
nilNode->next = 0;
// run the search-path -> augment-graph -> restore-trees loop
for(;;)
{
Vtx* v, *u;
int e0 = -1, ei = 0, ej = 0;
TWeight minWeight, weight;
uchar vt;
// grow S & T search trees, find an edge connecting them
while( first != nilNode )
{
v = first;
if( v->parent )
{
vt = v->t;
for( ei = v->first; ei != 0; ei = edgePtr[ei].next )
{
if( edgePtr[ei^vt].weight == 0 )
continue;
u = vtxPtr+edgePtr[ei].dst;
if( !u->parent )
{
u->t = vt;
u->parent = ei ^ 1;
u->ts = v->ts;
u->dist = v->dist + 1;
if( !u->next )
{
u->next = nilNode;
last = last->next = u;
}
continue;
}
if( u->t != vt )
{
e0 = ei ^ vt;
break;
}
if( u->dist > v->dist+1 && u->ts <= v->ts )
{
// reassign the parent
u->parent = ei ^ 1;
u->ts = v->ts;
u->dist = v->dist + 1;
}
}
if( e0 > 0 )
break;
}
// exclude the vertex from the active list
first = first->next;
v->next = 0;
}
if( e0 <= 0 )
break;
// find the minimum edge weight along the path
minWeight = edgePtr[e0].weight;
CV_Assert( minWeight > 0 );
// k = 1: source tree, k = 0: destination tree
for( int k = 1; k >= 0; k-- )
{
for( v = vtxPtr+edgePtr[e0^k].dst;; v = vtxPtr+edgePtr[ei].dst )
{
if( (ei = v->parent) < 0 )
break;
weight = edgePtr[ei^k].weight;
minWeight = MIN(minWeight, weight);
CV_Assert( minWeight > 0 );
}
weight = fabs(v->weight);
minWeight = MIN(minWeight, weight);
CV_Assert( minWeight > 0 );
}
// modify weights of the edges along the path and collect orphans
edgePtr[e0].weight -= minWeight;
edgePtr[e0^1].weight += minWeight;
flow += minWeight;
// k = 1: source tree, k = 0: destination tree
for( int k = 1; k >= 0; k-- )
{
for( v = vtxPtr+edgePtr[e0^k].dst;; v = vtxPtr+edgePtr[ei].dst )
{
if( (ei = v->parent) < 0 )
break;
edgePtr[ei^(k^1)].weight += minWeight;
if( (edgePtr[ei^k].weight -= minWeight) == 0 )
{
orphans.push_back(v);
v->parent = ORPHAN;
}
}
v->weight = v->weight + minWeight*(1-k*2);
if( v->weight == 0 )
{
orphans.push_back(v);
v->parent = ORPHAN;
}
}
// restore the search trees by finding new parents for the orphans
curr_ts++;
while( !orphans.empty() )
{
Vtx* v2 = orphans.back();
orphans.pop_back();
int d, minDist = INT_MAX;
e0 = 0;
vt = v2->t;
for( ei = v2->first; ei != 0; ei = edgePtr[ei].next )
{
if( edgePtr[ei^(vt^1)].weight == 0 )
continue;
u = vtxPtr+edgePtr[ei].dst;
if( u->t != vt || u->parent == 0 )
continue;
// compute the distance to the tree root
for( d = 0;; )
{
if( u->ts == curr_ts )
{
d += u->dist;
break;
}
ej = u->parent;
d++;
if( ej < 0 )
{
if( ej == ORPHAN )
d = INT_MAX-1;
else
{
u->ts = curr_ts;
u->dist = 1;
}
break;
}
u = vtxPtr+edgePtr[ej].dst;
}
// update the distance
if( ++d < INT_MAX )
{
if( d < minDist )
{
minDist = d;
e0 = ei;
}
for( u = vtxPtr+edgePtr[ei].dst; u->ts != curr_ts; u = vtxPtr+edgePtr[u->parent].dst )
{
u->ts = curr_ts;
u->dist = --d;
}
}
}
if( (v2->parent = e0) > 0 )
{
v2->ts = curr_ts;
v2->dist = minDist;
continue;
}
/* no parent is found */
v2->ts = 0;
for( ei = v2->first; ei != 0; ei = edgePtr[ei].next )
{
u = vtxPtr+edgePtr[ei].dst;
ej = u->parent;
if( u->t != vt || !ej )
continue;
if( edgePtr[ei^(vt^1)].weight && !u->next )
{
u->next = nilNode;
last = last->next = u;
}
if( ej > 0 && vtxPtr+edgePtr[ej].dst == v2 )
{
orphans.push_back(u);
u->parent = ORPHAN;
}
}
}
}
return flow;
}
template <class TWeight>
bool GCGraph<TWeight>::inSourceSegment( int i )
{
CV_Assert( i>=0 && i<(int)vtcs.size() );
return vtcs[i].t == 0;
}
#endif

@ -44,10 +44,7 @@
#include <time.h>
#include "opencv2/xphoto.hpp"
#include "../../../../opencv_main/modules/imgproc/src/gcgraph.hpp"
#include "opencv2/imgproc.hpp"
#include "opencv2/stitching.hpp"
#include "opencv2/core.hpp"
#include "opencv2/core/core_c.h"
@ -55,217 +52,69 @@
#include "opencv2/core/types.hpp"
#include "opencv2/core/types_c.h"
namespace cv
namespace xphotoInternal
{
template <typename Tp1, typename Tp2> static inline Tp2 sqr(Tp1 arg) { return arg * arg; }
template <typename Tp> static inline Tp sqr(Tp arg) { return arg * arg; }
static inline float norm2(const float &a, const float &b) { return sqr(a - b); }
static inline float norm2(const Vec2f &a, const Vec2f &b) { return (a - b).dot(a - b); }
# include "photomontage.hpp"
# include "annf.hpp"
}
static inline float norm2(const Vec3f &a, const Vec3f &b) { return (a - b).dot(a - b); }
static inline float norm2(const Vec4f &a, const Vec4f &b) { return (a - b).dot(a - b); }
#define dist(imgs, l1, l2, idx1, idx2) ( norm2(imgs[l1](idx1), imgs[l2](idx1)) \
+ norm2(imgs[l1](idx2), imgs[l2](idx2)) )
template <typename Tp>
static inline void setWeights(GCGraph <float> &graph,
const std::vector < Mat_<Tp> > &imgs, const int A, const int B,
const int labelA, const int labelB, const int alpha,
const Point &pointA, const Point &pointB)
{
//************************************************************//
//************************************************************//
if (labelA == labelB)
{
double weightAB = dist( imgs, labelA, alpha, pointA, pointB );
graph.addEdges(A, B, weightAB, weightAB);
}
else
{
double weightAX = dist( imgs, labelA, alpha, pointA, pointB );
double weightXB = dist( imgs, alpha, labelB, pointA, pointB );
double weightXSink = dist( imgs, labelA, labelB, pointA, pointB );
int X = graph.addVtx();
graph.addEdges(A, X, weightAX, weightAX);
graph.addEdges(X, B, weightXB, weightXB);
graph.addTermWeights(X, 0, weightXSink);
}
}
template <typename Tp>
static double alphaExpansion(const std::vector < Mat_<Tp> > &imgs,
const std::vector < Mat_<uchar> > &masks,
const Mat_<int> &labeling, const int alpha, Mat_<int> &nlabeling)
{
//************************************************************//
//************************************************************//
const int height = imgs[0].rows;
const int width = imgs[0].cols;
const double infinity = 10000000000;
const int actualEdges = height*(width - 1) + width*(height - 1);
GCGraph <float> graph(height*width + actualEdges, 2*actualEdges);
// terminal links
for (int i = 0; i < height; ++i)
{
const uchar *maskAlphaRow = masks[alpha].ptr(i);
const int *labelRow = (const int *) labeling.ptr(i);
for (int j = 0; j < width; ++j)
graph.addTermWeights( graph.addVtx(),
maskAlphaRow[j] ? 0 : infinity,
masks[ labelRow[j] ](i, j) ? 0 : infinity );
}
// neighbor links
for (int i = 0; i < height - 1; ++i)
{
const int *currentRow = (const int *) labeling.ptr(i);
const int *nextRow = (const int *) labeling.ptr(i + 1);
for (int j = 0; j < width - 1; ++j)
{
setWeights( graph, imgs, i*width + j, i*width + (j + 1),
currentRow[j], currentRow[j + 1], alpha,
Point(i, j), Point(i, j + 1) );
setWeights( graph, imgs, i*width + j, (i + 1)*width + j,
currentRow[j], nextRow[j], alpha,
Point(i, j), Point(i + 1, j) );
}
}
double result = graph.maxFlow();
nlabeling.create( labeling.size() );
for (int i = 0; i < height; ++i)
{
const int *inRow = (const int *) labeling.ptr(i);
int *outRow = (int *) nlabeling.ptr(i);
for (int j = 0; j < width; ++j)
{
bool gPart = graph.inSourceSegment(i*width + j);
outRow[j] = gPart ? inRow[j] : alpha;
}
}
return result;
}
template <typename Tp>
static void shiftMapInpaint(const Mat_<Tp> &src, const Mat_<uchar> &mask, Mat_<Tp> &dst)
namespace cv
{
template <typename Tp, unsigned int cn>
static void shiftMapInpaint(const Mat &src, const Mat &mask, Mat &dst)
{
//************************************************************//
//************************************************************//
const int nTransform = 60; // number of dominant transforms for stitching
const int psize = 8; // single ANNF patch size
const int width = src.cols;
const int height = src.rows;
Mat_<uchar> invMask = 255 - mask;
dilate(invMask, invMask, Mat(), Point(-1,-1), 2);
/** Downsample **/
//...
/** ANNF computation **/
int init = time(NULL);
srand( 1406297336 );
std::cout << init << std::endl;
srand( time(NULL) );
std::vector <Matx33f> transforms; // dominant transforms
for (int i = 0; i < nTransform; ++i)
{
float dx = rand()%width - width/2;
float dy = rand()%height - height/2;
float dx = rand()%src.cols - src.cols/2;
float dy = rand()%src.rows - src.rows/2;
transforms.push_back( Matx33f( 1, 0, dx,
0, 1, dy,
0, 0, 1) );
}
/** Warping **/
std::vector < Mat_<Tp> > imgs( nTransform + 1 ); // source image transformed with transforms[i]
std::vector < Mat_<uchar> > masks( nTransform + 1 ); // validity mask for current shift
std::vector <Mat> images( nTransform + 1 ); // source image transformed with transforms
std::vector <Mat> masks( nTransform + 1 ); // definition domain for current shift
Mat_<uchar> invMask = 255 - mask;
dilate(invMask, invMask, Mat(), Point(-1,-1), 2);
src.copyTo( imgs[0] );
src.convertTo( images[0], CV_32F );
mask.copyTo( masks[0] );
for (int i = 0; i < nTransform; ++i)
{
Mat_<Tp> nsrc( src.size() );
warpPerspective( src, nsrc, transforms[i], src.size(),
INTER_LINEAR, BORDER_CONSTANT, 0 );
warpPerspective( images[0], images[i + 1], transforms[i],
images[0].size(), INTER_LINEAR );
Mat_<uchar> nmask( mask.size(), mask.type() );
warpPerspective( mask, nmask, transforms[i], mask.size(),
INTER_NEAREST, BORDER_CONSTANT, 0 );
nmask &= invMask;
nsrc.copyTo( imgs[i + 1] );
nmask.copyTo( masks[i + 1] );
warpPerspective( masks[0], masks[i + 1], transforms[i],
masks[0].size(), INTER_NEAREST);
masks[i + 1] &= invMask;
}
/** Stitching **/
std::vector <double> costs( nTransform + 1 );
std::vector < Mat_<int> > labelings( nTransform + 1 );
Mat_<int> labeling( height, width, 0 );
double cost = std::numeric_limits<double>::max();
for (int success = false, num = 0; ; success = false)
{
for (int i = 0; i < nTransform + 1; ++i)
costs[i] = alphaExpansion(imgs, masks, labeling, i, labelings[i]);
for (int i = 0; i < nTransform + 1; ++i)
if (costs[i] < 0.98*cost)
{
success = true;
cost = costs[num = i];
}
if (success == false)
break;
labelings[num].copyTo(labeling);
}
for (int k = 0; k < height*width; ++k)
{
int i = k / width;
int j = k % width;
dst(i, j) = imgs[labeling(i, j)](i, j);
}
/** Upsample and refinement **/
//...
Mat photomontageResult;
xphotoInternal::Photomontage < cv::Vec <float, cn> >( images, masks )
.assignResImage(photomontageResult);
photomontageResult.convertTo( dst, dst.type() );
}
template <typename Tp>
void inpaint(const Mat_<Tp> src, const Mat_<uchar> mask, Mat_<Tp> dst, const int algorithmType)
template <typename Tp, unsigned int cn>
void inpaint(const Mat &src, const Mat &mask, Mat &dst, const int algorithmType)
{
//************************************************************//
//************************************************************//
dst.create( src.size() );
dst.create( src.size(), src.type() );
switch ( algorithmType )
{
case INPAINT_SHIFTMAP:
shiftMapInpaint(src, mask, dst);
shiftMapInpaint <Tp, cn>(src, mask, dst);
break;
default:
CV_Assert( false );
@ -281,73 +130,70 @@ namespace cv
*/
void inpaint(const Mat &src, const Mat &mask, Mat &dst, const int algorithmType)
{
//************************************************************//
//************************************************************//
CV_Assert( mask.channels() == 1 && mask.depth() == CV_8U );
CV_Assert( src.rows == mask.rows && src.cols == mask.cols );
switch ( src.type() )
{
case CV_8UC1:
inpaint( Mat_<uchar>(src), Mat_<uchar>(mask), Mat_<uchar>(dst), algorithmType );
inpaint <uchar, 1>( src, mask, dst, algorithmType );
break;
case CV_8UC2:
inpaint( Mat_<Vec2b>(src), Mat_<uchar>(mask), Mat_<Vec2b>(dst), algorithmType );
inpaint <uchar, 2>( src, mask, dst, algorithmType );
break;
case CV_8UC3:
inpaint( Mat_<Vec3b>(src), Mat_<uchar>(mask), Mat_<Vec3b>(dst), algorithmType );
inpaint <uchar, 3>( src, mask, dst, algorithmType );
break;
case CV_8UC4:
inpaint( Mat_<Vec4b>(src), Mat_<uchar>(mask), Mat_<Vec4b>(dst), algorithmType );
inpaint <uchar, 4>( src, mask, dst, algorithmType );
break;
case CV_16SC1:
inpaint( Mat_<short>(src), Mat_<uchar>(mask), Mat_<short>(dst), algorithmType );
inpaint <short, 1>( src, mask, dst, algorithmType );
break;
case CV_16SC2:
inpaint( Mat_<Vec2s>(src), Mat_<uchar>(mask), Mat_<Vec2s>(dst), algorithmType );
inpaint <short, 2>( src, mask, dst, algorithmType );
break;
case CV_16SC3:
inpaint( Mat_<Vec3s>(src), Mat_<uchar>(mask), Mat_<Vec3s>(dst), algorithmType );
inpaint <short, 3>( src, mask, dst, algorithmType );
break;
case CV_16SC4:
inpaint( Mat_<Vec4s>(src), Mat_<uchar>(mask), Mat_<Vec4s>(dst), algorithmType );
inpaint <short, 4>( src, mask, dst, algorithmType );
break;
case CV_32SC1:
inpaint( Mat_<int>(src), Mat_<uchar>(mask), Mat_<int>(dst), algorithmType );
inpaint <int, 1>( src, mask, dst, algorithmType );
break;
case CV_32SC2:
inpaint( Mat_<Vec2i>(src), Mat_<uchar>(mask), Mat_<Vec2i>(dst), algorithmType );
inpaint <int, 2>( src, mask, dst, algorithmType );
break;
case CV_32SC3:
inpaint( Mat_<Vec3i>(src), Mat_<uchar>(mask), Mat_<Vec3i>(dst), algorithmType );
inpaint <int, 3>( src, mask, dst, algorithmType );
break;
case CV_32SC4:
inpaint( Mat_<Vec4i>(src), Mat_<uchar>(mask), Mat_<Vec4i>(dst), algorithmType );
inpaint <int, 4>( src, mask, dst, algorithmType );
break;
case CV_32FC1:
inpaint( Mat_<float>(src), Mat_<uchar>(mask), Mat_<float>(dst), algorithmType);
inpaint <float, 1>( src, mask, dst, algorithmType );
break;
case CV_32FC2:
inpaint( Mat_<Vec2f>(src), Mat_<uchar>(mask), Mat_<Vec2f>(dst), algorithmType );
inpaint <float, 2>( src, mask, dst, algorithmType );
break;
case CV_32FC3:
inpaint( Mat_<Vec3f>(src), Mat_<uchar>(mask), Mat_<Vec3f>(dst), algorithmType );
inpaint <float, 3>( src, mask, dst, algorithmType );
break;
case CV_32FC4:
inpaint( Mat_<Vec4f>(src), Mat_<uchar>(mask), Mat_<Vec4f>(dst), algorithmType );
inpaint <float, 4>( src, mask, dst, algorithmType );
break;
case CV_64FC1:
inpaint( Mat_<double>(src), Mat_<uchar>(mask), Mat_<double>(dst), algorithmType );
inpaint <double, 1>( src, mask, dst, algorithmType );
break;
case CV_64FC2:
inpaint( Mat_<Vec2d>(src), Mat_<uchar>(mask), Mat_<Vec2d>(dst), algorithmType );
inpaint <double, 2>( src, mask, dst, algorithmType );
break;
case CV_64FC3:
inpaint( Mat_<Vec3d>(src), Mat_<uchar>(mask), Mat_<Vec3d>(dst), algorithmType );
inpaint <double, 3>( src, mask, dst, algorithmType );
break;
case CV_64FC4:
inpaint( Mat_<Vec4d>(src), Mat_<uchar>(mask), Mat_<Vec4d>(dst), algorithmType );
inpaint <double, 4>( src, mask, dst, algorithmType );
break;
default:
CV_Assert( false );

@ -0,0 +1,52 @@
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
// Copyright (C) 2009-2011, Willow Garage Inc., all rights reserved.
// Third party copyrights are property of their respective owners.
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of Intel Corporation may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __NORM2_HPP__
#define __NORM2_HPP__
template <typename Tp> static inline Tp sqr(Tp x) { return x*x; }
template <typename Tp, int cn> static inline Tp sqr( cv::Vec<Tp, cn> x) { return x.dot(x); }
template <typename Tp> static inline Tp norm2(const Tp &a, const Tp &b) { return sqr(a - b); }
template <typename Tp, int cn> static inline
Tp norm2(const cv::Vec <Tp, cn> &a, const cv::Vec<Tp, cn> &b) { return sqr(a - b); }
#endif /* __NORM2_HPP__ */

@ -0,0 +1,259 @@
/*M///////////////////////////////////////////////////////////////////////////////////////
//
// IMPORTANT: READ BEFORE DOWNLOADING, COPYING, INSTALLING OR USING.
//
// By downloading, copying, installing or using the software you agree to this license.
// If you do not agree to this license, do not download, install,
// copy or use the software.
//
//
// License Agreement
// For Open Source Computer Vision Library
//
// Copyright (C) 2000-2008, Intel Corporation, all rights reserved.
// Copyright (C) 2009-2011, Willow Garage Inc., all rights reserved.
// Third party copyrights are property of their respective owners.
//
// * Redistribution's of source code must retain the above copyright notice,
// this list of conditions and the following disclaimer.
//
// * Redistribution's in binary form must reproduce the above copyright notice,
// this list of conditions and the following disclaimer in the documentation
// and/or other materials provided with the distribution.
//
// * The name of Intel Corporation may not be used to endorse or promote products
// derived from this software without specific prior written permission.
//
// This software is provided by the copyright holders and contributors "as is" and
// any express or implied warranties, including, but not limited to, the implied
// warranties of merchantability and fitness for a particular purpose are disclaimed.
// In no event shall the Intel Corporation or contributors be liable for any direct,
// indirect, incidental, special, exemplary, or consequential damages
// (including, but not limited to, procurement of substitute goods or services;
// loss of use, data, or profits; or business interruption) however caused
// and on any theory of liability, whether in contract, strict liability,
// or tort (including negligence or otherwise) arising in any way out of
// the use of this software, even if advised of the possibility of such damage.
//
//M*/
#ifndef __PHOTOMONTAGE_HPP__
#define __PHOTOMONTAGE_HPP__
#include <opencv2/core.hpp>
#include "norm2.hpp"
#include "algo.hpp"
#include "annf.hpp"
#include "gcgraph.hpp"
#define GCInfinity 10*1000*1000*1000.0
#define EFFECTIVE_HEIGHT 600
#define EFFECTIVE_WIDTH 800
template <typename Tp> class Photomontage
{
private:
const std::vector <cv::Mat> &images; // vector of images for different labels
const std::vector <cv::Mat> &masks; // vector of definition domains for each image
std::vector <cv::Mat> labelings; // vector of labelings for different expansions
std::vector <double> distances; // vector of max-flow costs for different labeling
const int height;
const int width;
const int type;
const int channels;
const int lsize;
bool multiscale; // if true, Photomontage use coarse-to-fine scheme
double singleExpansion(const cv::Mat &x_i, const int alpha); // single neighbor computing
void gradientDescent(const cv::Mat &x_0, cv::Mat x_n); // gradient descent in alpha-expansion topology
protected:
virtual double dist(const Tp &l1p1, const Tp &l1p2, const Tp &l2p1, const Tp &l2p2);
virtual void setWeights(GCGraph <double> &graph, const cv::Point &pA, const cv::Point &pB, const int lA, const int lB, const int lX);
public:
Photomontage(const std::vector <cv::Mat> &images, const std::vector <cv::Mat> &masks, bool multiscale = true);
virtual ~Photomontage(){};
void assignLabeling(cv::Mat &img);
void assignResImage(cv::Mat &img);
};
template <typename Tp> double Photomontage <Tp>::
dist(const Tp &l1p1, const Tp &l1p2, const Tp &l2p1, const Tp &l2p2)
{
return norm2(l1p1, l2p1) + norm2(l1p2, l2p2);
}
template <typename Tp> void Photomontage <Tp>::
setWeights(GCGraph <double> &graph, const cv::Point &pA, const cv::Point &pB, const int lA, const int lB, const int lX)
{
if (lA == lB)
{
/** Link from A to B **/
double weightAB = dist( images[lA].at<Tp>(pA),
images[lA].at<Tp>(pB),
images[lX].at<Tp>(pA),
images[lX].at<Tp>(pB) );
graph.addEdges(pA.y*width + pA.x, pB.y*width + pB.x, weightAB, weightAB);
}
else
{
int X = graph.addVtx();
/** Link from X to sink **/
double weightXS = dist( images[lA].at<Tp>(pA),
images[lA].at<Tp>(pB),
images[lB].at<Tp>(pA),
images[lB].at<Tp>(pB) );
graph.addTermWeights(X, 0, weightXS);
/** Link from A to X **/
double weightAX = dist( images[lA].at<Tp>(pA),
images[lA].at<Tp>(pB),
images[lX].at<Tp>(pA),
images[lX].at<Tp>(pB) );
graph.addEdges(pA.y*width + pA.x, X, weightAX, weightAX);
/** Link from X to B **/
double weightXB = dist( images[lX].at<Tp>(pA),
images[lX].at<Tp>(pB),
images[lB].at<Tp>(pA),
images[lB].at<Tp>(pB) );
graph.addEdges(X, pB.y*width + pB.x, weightXB, weightXB);
}
}
template <typename Tp> double Photomontage <Tp>::
singleExpansion(const cv::Mat &x_i, const int alpha)
{
int actualEdges = (height - 1)*width + height*(width - 1);
GCGraph <double> graph(actualEdges + height*width, 2*actualEdges);
/** Terminal links **/
for (int i = 0; i < height; ++i)
{
const uchar *maskAlphaRow = masks[alpha].ptr <uchar>(i);
const int *labelRow = (const int *) x_i.ptr <int>(i);
for (int j = 0; j < width; ++j)
graph.addTermWeights( graph.addVtx(),
maskAlphaRow[j] ? 0 : GCInfinity,
masks[ labelRow[j] ].at<uchar>(i, j) ? 0 : GCInfinity );
}
/** Neighbor links **/
for (int i = 0; i < height - 1; ++i)
{
const int *currentRow = (const int *) x_i.ptr <int>(i);
const int *nextRow = (const int *) x_i.ptr <int>(i + 1);
for (int j = 0; j < width - 1; ++j)
{
setWeights( graph, cv::Point(i, j), cv::Point(i, j + 1), currentRow[j], currentRow[j + 1], alpha );
setWeights( graph, cv::Point(i, j), cv::Point(i + 1, j), currentRow[j], nextRow[j], alpha );
}
setWeights( graph, cv::Point(i, width - 1), cv::Point(i + 1, width - 1),
currentRow[width - 1], nextRow[width - 1], alpha );
}
const int *currentRow = (const int *) x_i.ptr <int>(height - 1);
for (int i = 0; i < width - 1; ++i)
setWeights( graph, cv::Point(height - 1, i), cv::Point(height - 1, i + 1),
currentRow[i], currentRow[i + 1], alpha );
/** Max-flow computation **/
double result = graph.maxFlow();
/** Writing results **/
labelings[alpha].create( height, width, CV_32SC1 );
for (int i = 0; i < height; ++i)
{
const int *inRow = (const int *) x_i.ptr <int>(i);
int *outRow = (int *) labelings[alpha].ptr <int>(i);
for (int j = 0; j < width; ++j)
outRow[j] = graph.inSourceSegment(i*width + j) ? inRow[j] : alpha;
}
return result;
}
template <typename Tp> void Photomontage <Tp>::
gradientDescent(const cv::Mat &x_0, cv::Mat x_n)
{
double optValue = std::numeric_limits<double>::max();
x_0.copyTo(x_n);
for (int num = -1; /**/; num = -1)
{
for (int i = 0; i < lsize; ++i)
distances[i] = singleExpansion(x_n, i);
int minIndex = min_idx(distances);
double minValue = distances[minIndex];
if (minValue < 0.98*optValue)
optValue = distances[num = minIndex];
if (num == -1)
break;
labelings[num].copyTo(x_n);
}
}
template <typename Tp> void Photomontage <Tp>::
assignLabeling(cv::Mat &img)
{
if (multiscale == 0 || (height < EFFECTIVE_HEIGHT && width < EFFECTIVE_WIDTH))
{
img.create(height, width, CV_32SC1);
img.setTo(0);
gradientDescent(img, img);
}
else
{
int l = std::min( cvRound(height/600.0), cvRound(width/800.0) );
img.create( cv::Size(width/l, height/l), CV_32SC1 );
...
img.setTo(0);
gradientDescent(img, img);
resize(img, img, cv::Size(height, width), 0.0, 0.0, cv::INTER_NEAREST);
...
}
}
template <typename Tp> void Photomontage <Tp>::
assignResImage(cv::Mat &img)
{
cv::Mat optimalLabeling;
assignLabeling(optimalLabeling);
img.create( height, width, type );
for (int i = 0; i < height; ++i)
for (int j = 0; j < width; ++j)
img.at<Tp>(i, j) = images[ optimalLabeling.at<int>(i, j) ].at<Tp>(i, j);
}
template <typename Tp> Photomontage <Tp>::
Photomontage(const std::vector <cv::Mat> &images, const std::vector <cv::Mat> &masks, const bool multiscale)
:
images(images), masks(masks), multiscale(multiscale), height(images[0].rows), width(images[0].cols), type(images[0].type()),
channels(images[0].channels()), lsize(images.size()), labelings(images.size()), distances(images.size())
{
CV_Assert(images[0].depth() != CV_8U && masks[0].depth() == CV_8U);
}
#endif /* __PHOTOMONTAGE_HPP__ */
Loading…
Cancel
Save