tahoma2d/toonz/sources/toonzlib/tcenterlinepolygonizer.cpp
Campbell Barton b3bd842e04 Make functions static, ensure declarations match headers (#610)
This patch used -Wmissing-declarations warning
to show functions and symbols that had no declarations, and either:

- Make static
- Add to header

This helps avoid possability that declarations and functions get out of sync.
And ensures all source files reference headers correctly.

It also makes sure functions defined with extern "C",
have this defined in the header. An error found in calligraph.h while writing this patch.

This has been applied to toonzlib, to avoid making very large global changes.
If accepted, -Wmissing-declarations warning could be added to CMake.
2016-07-13 21:05:06 +09:00

877 lines
28 KiB
C++

#include "tcenterlinevectP.h"
//==========================================================================
//************************
//* Polygonization *
//************************
//--------------------------------------------------------------------------
//===============================
// Raw Borders Extraction
//===============================
// Raw contour class definition
class RawBorderPoint {
TPoint m_position;
int m_ambiguousTurn; // used to remember cases of multiple turning directions
// in a RawBorder extraction.
public:
RawBorderPoint() : m_ambiguousTurn(0) {}
RawBorderPoint(int i, int j) : m_position(i, j), m_ambiguousTurn(0) {}
inline TPoint pos() const { return m_position; }
inline int x() const { return m_position.x; }
inline int y() const { return m_position.y; }
enum { left = 1, right = 2 }; // Direction taken at ambiguous turning point
inline int getAmbiguous() const { return m_ambiguousTurn; }
inline void setAmbiguous(int direction) { m_ambiguousTurn = direction; }
};
//--------------------------------------------------------------------------
class RawBorder final : public std::vector<RawBorderPoint> {
int m_xExternal; // x coordinate of a specific vertex in the outer
// RawBorder which contains this inner one.
TPointD *m_coordinateSums;
TPointD *m_coordinateSquareSums;
double *m_coordinateMixedSums;
public:
RawBorder() {}
~RawBorder() {}
void setXExternalPixel(int a) { m_xExternal = a; }
int xExternalPixel() { return m_xExternal; }
TPointD *&sums() { return m_coordinateSums; }
TPointD *&sums2() { return m_coordinateSquareSums; }
double *&sumsMix() { return m_coordinateMixedSums; }
};
//--------------------------------------------------------------------------
// Of course we don't want RawBorders to be entirely copied whenever STL
// requires to resize a BorderFamily...
typedef std::vector<RawBorder *> BorderFamily;
typedef std::vector<BorderFamily> BorderList;
//--------------------------------------------------------------------------
//==========================================================================
//============================
// Polygonizer Locals
//============================
namespace {
// Const names
enum { white = 0, black = 1 };
enum { inner = 0, outer = 1, none = 2, invalid = 3 };
}
//=======================================================================================
//-------------------------------
// Raster Data Functions
//-------------------------------
// NOTA: Il tono di un TPixelCM32 rappresenta la transizione tra colore ink e
// colore paint.
// di solito, se il tono e' basso abbiamo un colore ink - che puo' anche
// essere bianco,
// in teoria...
// Sarebbe opportuno che il vettorizzatore riconoscesse - per colormap - il
// colore
// delle strokes.
// Approcci: a) la Signaturemap diventa *piu'* di una bitmap e si seguono
// le outline
// dei singoli colori.
// => Sconnessioni tra i colori adiacenti. Bisogna introdurre la
// distanza tra
// colori per seguire l'outline (ossia, i pixel tenuti a dx
// della outline
// devono essere simili).
//
// b) Una volta vettorizzato tutto, si sceglie il colore della
// stroke.
// E' possibile controllare il colore sui vertici delle
// sequenze semplificate
// e fare una media.
// NOTE: Transparency makes colors fade to white. Full transparent black pixels
// are considered white.
//--------------------------------------------------------------------------
template <typename T>
class PixelEvaluator {
TRasterPT<T> m_ras;
int m_threshold;
public:
PixelEvaluator(const TRasterPT<T> &ras, int threshold)
: m_ras(ras), m_threshold(threshold) {}
inline unsigned char getBlackOrWhite(int x, int y);
};
//--------------------------------------------------------------------------
template <>
inline unsigned char PixelEvaluator<TPixel32>::getBlackOrWhite(int x, int y) {
// return ras->pixels(y)[x].r + 2 * ras->pixels(y)[x].g + ras->pixels(y)[x].b
// <
// threshold * (ras->pixels(y)[x].m / 255.0);
// NOTE: Green is considered twice brighter than red or blue channel.
// Using Value of HSV color model
return std::max(m_ras->pixels(y)[x].r,
std::max(m_ras->pixels(y)[x].g, m_ras->pixels(y)[x].b)) <
m_threshold * (m_ras->pixels(y)[x].m / 255.0);
// Using Lightness of HSV color model
// return (max(ras->pixels(y)[x].r, max(ras->pixels(y)[x].g,
// ras->pixels(y)[x].b)) +
// min(ras->pixels(y)[x].r, min(ras->pixels(y)[x].g,
// ras->pixels(y)[x].b))) / 2.0 <
// threshold * (ras->pixels(y)[x].m / 255.0);
// Using (relative) Luminance
// return 0.2126 * ras->pixels(y)[x].r + 0.7152 * ras->pixels(y)[x].g + 0.0722
// * ras->pixels(y)[x].b <
// threshold * (ras->pixels(y)[x].m / 255.0);
}
template <>
inline unsigned char PixelEvaluator<TPixelGR8>::getBlackOrWhite(int x, int y) {
return m_ras->pixels(y)[x].value < m_threshold;
}
template <>
inline unsigned char PixelEvaluator<TPixelCM32>::getBlackOrWhite(int x, int y) {
return m_ras->pixels(y)[x].getTone() < m_threshold;
}
//--------------------------------------------------------------------------
// Signaturemap format:
// stores a map of bytes, whose first bit represents the color (black/white) of
// corresponding pixel, and
// the rest its 'signature', used as an int to store information.
// NOTE: given a TRaster32, the corresponding Signaturemap constructed is
// intended 0(white)-padded
class Signaturemap {
std::unique_ptr<unsigned char[]> m_array;
int m_rowSize;
int m_colSize;
public:
Signaturemap(const TRasterP &ras, int threshold);
template <typename T>
void readRasterData(const TRasterPT<T> &ras, int threshold);
inline int getRowSize() const { return m_rowSize; }
inline int getColSize() const { return m_colSize; }
unsigned char *pixelByte(int x, int y) {
return &m_array[(y + 1) * m_rowSize + x + 1];
}
bool getBitmapColor(int x, int y) const {
return m_array[(y + 1) * m_rowSize + x + 1] & 1;
}
inline unsigned char getSignature(int x, int y) const {
return m_array[(y + 1) * m_rowSize + x + 1] >> 1;
}
void setSignature(int x, int y, int val) {
unsigned char *pixel = pixelByte(x, y);
*pixel &= 1;
*pixel |= (val << 1); // Si puo' fare meglio??
}
};
//--------------------------------------------------------------------------
Signaturemap::Signaturemap(const TRasterP &ras, int threshold) {
// Extrapolate raster type
TRaster32P rr = (TRaster32P)ras;
TRasterGR8P rgr = (TRasterGR8P)ras;
TRasterCM32P rt = (TRasterCM32P)ras;
assert(rr || rgr || rt);
// Read raster data
if (rr) {
rr->lock();
readRasterData(rr, threshold);
rr->unlock();
} else if (rgr) {
rgr->lock();
readRasterData(rgr, threshold);
rgr->unlock();
} else {
rt->lock();
readRasterData(rt, threshold);
rt->unlock();
}
}
//--------------------------------------------------------------------------
template <typename T>
void Signaturemap::readRasterData(const TRasterPT<T> &ras, int threshold) {
unsigned char *currByte;
int x, y;
PixelEvaluator<T> evaluator(ras, threshold);
m_rowSize = ras->getLx() + 2;
m_colSize = ras->getLy() + 2;
m_array.reset(new unsigned char[m_rowSize * m_colSize]);
memset(m_array.get(), none << 1, m_rowSize);
currByte = m_array.get() + m_rowSize;
for (y = 0; y < ras->getLy(); ++y) {
*currByte = none << 1;
currByte++;
for (x = 0; x < ras->getLx(); ++x, ++currByte)
*currByte = evaluator.getBlackOrWhite(x, y) | (none << 1);
*currByte = none << 1;
currByte++;
}
memset(currByte, none << 1, m_rowSize);
}
//--------------------------------------------------------------------------
// Minority check for amiguous turning directions
inline bool getMinorityCheck(const Signaturemap &ras, int x, int y) {
// Assumes (x,y) is ambiguous case: 2 immediate surrounding pixels are white
// and 2 black
return (ras.getBitmapColor(x + 1, y) + ras.getBitmapColor(x + 1, y - 1) +
ras.getBitmapColor(x - 2, y) + ras.getBitmapColor(x - 2, y - 1) +
ras.getBitmapColor(x - 1, y + 1) + ras.getBitmapColor(x - 1, y - 2) +
ras.getBitmapColor(x, y + 1) + ras.getBitmapColor(x, y - 2)) > 4;
}
//--------------------------------------------------------------------------
// Sets signature of a given border
inline void setSignature(Signaturemap &ras, const RawBorder &border, int val) {
unsigned int j;
int yOld;
// Set border's alpha channel
yOld = border.back().y();
for (j = 0; j < border.size(); ++j) {
if (border[j].y() < yOld) {
ras.setSignature(border[j].x(), border[j].y(), val);
} else if (border[j].y() > yOld) {
ras.setSignature(border[j].x(), yOld, val);
}
yOld = border[j].y();
}
}
//==========================================================================
//-------------------------------
// Raw Borders Extraction
//-------------------------------
// RawBorderPoints correspond to lower-left pixel corners.
// EXAMPLE: (0,0) is the lower-left *corner* of the image, whereas (0,0) also
// represents coordinates of the lower-left *pixel*.
// NOTE: 'Ambiguous turning' vertices are those of kind:
//
// B|W W|B
// -x- -or- -x-
// W|B B|W
//
// Keeping B on the right of our path-seeking direction, we may either turn
// left or right at these points.
static RawBorder *extractPath(Signaturemap &ras, int x0, int y0, int pathType,
int xOuterPixel, int despeckling) {
RawBorder *path = new RawBorder;
int x, y;
short dirX, dirY;
long int area = 0;
bool nextLeftPixel, nextRightPixel;
if (pathType == outer) {
dirX = 0;
dirY = 1;
} else {
dirX = 1;
dirY = 0;
area += y0;
path->setXExternalPixel(xOuterPixel);
}
path->push_back(RawBorderPoint(x0, y0));
// Check here if (x0, y0) is an ambiguous-direction point
nextLeftPixel = ras.getBitmapColor(x0 + (dirY - dirX - 1) / 2,
y0 + (-dirY - dirX - 1) / 2);
nextRightPixel = ras.getBitmapColor(x0 + (-dirX - dirY - 1) / 2,
y0 + (dirX - dirY - 1) / 2);
if ((nextRightPixel == black) && (nextLeftPixel == white))
path->back().setAmbiguous(dirX ? RawBorderPoint::left
: RawBorderPoint::right);
// Begin path extraction
for (x = x0 + dirX, y = y0 + dirY; !(x == x0 && y == y0);
x += dirX, y += dirY) {
path->push_back(RawBorderPoint(x, y));
// Calculate next direction
nextLeftPixel = ras.getBitmapColor(x + (dirX - dirY - 1) / 2,
y + (dirY + dirX - 1) / 2);
nextRightPixel = ras.getBitmapColor(x + (dirX + dirY - 1) / 2,
y + (dirY - dirX - 1) / 2);
if ((nextRightPixel == black) && (nextLeftPixel == black)) {
// Left Turn
std::swap(dirY, dirX);
dirX = -dirX;
} else if ((nextRightPixel == white) && (nextLeftPixel == white)) {
// Right Turn
std::swap(dirY, dirX);
dirY = -dirY;
} else if ((nextRightPixel == white) && (nextLeftPixel == black)) {
// path->back().setAmbiguous();
// Do a surrounding check and connect minority color
if (getMinorityCheck(ras, x, y) == black) {
std::swap(dirY, dirX);
dirY = -dirY;
path->back().setAmbiguous(RawBorderPoint::right);
} // right turn
else {
std::swap(dirY, dirX);
dirX = -dirX;
path->back().setAmbiguous(RawBorderPoint::left);
} // left turn
}
// Also calculate border area
area += y * dirX;
// And sign treated pixel
if (dirY != 0) ras.setSignature(x, y + (dirY - 1) / 2, pathType);
}
// If the inner region's overall area is under a given threshold,
// then erase it (intended as image noise).
if (abs(area) < despeckling) {
setSignature(ras, *path, invalid);
delete path;
path = 0;
}
return path;
}
//--------------------------------------------------------------------------
static BorderList *extractBorders(const TRasterP &ras, int threshold,
int despeckling) {
Signaturemap byteImage(ras, threshold);
BorderList *borderHierarchy = new BorderList;
std::vector<RawBorder *> outerBorders;
std::list<RawBorder *> innerBorders;
RawBorder *foundPath;
int x, y;
bool Color, oldColor;
int xOuterPixel = 0;
bool enteredRegionType;
unsigned char signature;
// Traverse image to extract raw borders
for (y = 0; y < ras->getLy(); ++y) {
oldColor = white;
enteredRegionType = outer;
for (x = 0; x < ras->getLx(); ++x) {
if (oldColor ^ (Color = byteImage.getBitmapColor(x, y))) {
// Region type changes
enteredRegionType = !enteredRegionType;
if ((signature = byteImage.getSignature(x, y)) == none) {
// We've found a border
if ((foundPath = extractPath(byteImage, x, y, !enteredRegionType,
xOuterPixel, despeckling)))
if (enteredRegionType == outer)
innerBorders.push_back(foundPath);
else
outerBorders.push_back(foundPath);
}
// If leaving a white region, remember it - in order to establish
// border hierarchy in the future
if (enteredRegionType == inner && signature != invalid) xOuterPixel = x;
// Invalid pixels got signed by a cut-out path, due to insufficient area
if (signature == invalid)
byteImage.setSignature(x, y, none); // Restore them now
oldColor = Color;
}
}
}
// Now, we have all borders found, but no hierarchy between them.
unsigned int i;
std::list<RawBorder *>::iterator l;
// Build hierarchy
innerBorders.push_front(0); // Just to keep a fixed list head
for (i = 0; i < outerBorders.size(); ++i) {
// Initialize a border family
borderHierarchy->push_back(BorderFamily());
borderHierarchy->back().push_back(outerBorders[i]);
// Reset outerBorders[i]'s signature
setSignature(byteImage, *outerBorders[i], none);
// Now check inner borders for insideness - check if the outerPixel
// remembered in path extraction has been cleared
for (l = innerBorders.begin(), ++l; l != innerBorders.end(); ++l) {
if (byteImage.getSignature((*l)->xExternalPixel(), (**l)[0].y()) ==
none) {
borderHierarchy->back().push_back(*l);
setSignature(byteImage, **l, none);
l = innerBorders.erase(l);
--l;
}
}
}
return borderHierarchy;
}
//--------------------------------------------------------------------------
//==================================
// Calculate optimal polygons
//==================================
// The optimal polygon for a given original border is found like:
// 1) Find couples (i,k(i)), so that k(i) be the largest k:
// d(j,ik) <= 1; for *all* i<j<k. (d is infinite-norm distance)
// It can be shown that such a condition is equivalent to:
// exists line l : d(l,j)<=1/2, for all i<=j<=k(i).
// 2) Clean the above couples - find couples (i,l(i)):
// l(i)=min{k(j)}, j=i..n.
// 3) Calculate clipped couples (i',l'); where i'=i+1, l'=l(i)-1.
// 4) Calculate sums for path penalties.
// 5) Apply optimality algorithm.
// NOTE: Weak simpleness reads like: a set of polygons is weak-simple if no edge
// *crosses* another edge. Superposition and collision of edges with
// vertices
// are still admitted.
// => It can be shown that due to 1) and special conditions on ambiguous
// turnings applied in both 1) and 3), weak simpleness is insured in
// our polygonization.
//--------------------------------------------------------------------------
// Helper functions/classes: circular-indexed vectors
// returns 1 whenever the triple (a,b,c) is 'circular' mod n.
// NOTE: We'll find useful taking (i,i,j) as 1 and (i,j,j) as 0.
inline bool isCircular(int a, int b, int c) {
return a <= c ? a <= b && b < c : c > b || b >= a;
}
//--------------------------------------------------------------------------
// Extracts a 'next corner' array - helps improving overall speed
inline std::unique_ptr<int[]> findNextCorners(RawBorder &path) {
std::unique_ptr<int[]> corners(new int[path.size()]);
// NOTE: 0 is a corner, due to the path extraction procedure.
int currentCorner = 0;
for (int i = path.size() - 1; i >= 0; --i) {
if (path[currentCorner].x() != path[i].x() &&
path[currentCorner].y() != path[i].y())
currentCorner = i + 1;
corners[i] = currentCorner;
}
return corners;
}
//--------------------------------------------------------------------------
// Calculate furthest k satisfying 1) for all fixed i.
inline std::unique_ptr<int[]> furthestKs(RawBorder &path,
std::unique_ptr<int[]> &nextCorners) {
int n = path.size();
std::unique_ptr<int[]> kVector(new int[n]);
enum { left, up, right, down };
int directionsOccurred[4];
nextCorners = findNextCorners(path);
int i, j, k;
TPoint shift;
TPoint leftConstraint, rightConstraint, violatedConstraint;
TPoint newLeftConstraint, newRightConstraint;
TPoint jPoint, jNextPoint, iPoint, direction;
int directionSignature;
for (i = 0; i < n; ++i) {
// Initialize search
leftConstraint = rightConstraint = TPoint();
directionsOccurred[0] = directionsOccurred[1] = directionsOccurred[2] =
directionsOccurred[3] = 0;
j = i;
jNextPoint = iPoint = path[i].pos();
// Search for k(i)
while (1) {
// NOTE: Here using TPoint::operator= is less effective than setting
// its x and y components directly...
jPoint = jNextPoint;
jNextPoint = path[nextCorners[j]].pos();
// Update directions count
directionSignature = jNextPoint.x > jPoint.x
? right
: jNextPoint.x < jPoint.x
? left
: jNextPoint.y > jPoint.y ? up : down;
directionsOccurred[directionSignature] = 1;
// If all 4 axis directions occurred, quit
if (directionsOccurred[left] && directionsOccurred[right] &&
directionsOccurred[up] && directionsOccurred[down]) {
k = j;
goto foundK;
}
// Update displacement from i
shift = jNextPoint - iPoint;
// Test j against constraints
// if(cross(shift, leftConstraint)<0 || cross(shift, rightConstraint)>0)
if (cross(shift, leftConstraint) < 0) {
violatedConstraint = leftConstraint;
break;
}
if (cross(shift, rightConstraint) > 0) {
violatedConstraint = rightConstraint;
break;
}
// Update constraints
if (abs(shift.x) > 1 || abs(shift.y) > 1) {
newLeftConstraint.x =
shift.x + (shift.y < 0 || (shift.y == 0 && shift.x < 0) ? 1 : -1);
newLeftConstraint.y =
shift.y + (shift.x > 0 || (shift.x == 0 && shift.y < 0) ? 1 : -1);
if (cross(newLeftConstraint, leftConstraint) >= 0)
leftConstraint = newLeftConstraint;
newRightConstraint.x =
shift.x + (shift.y > 0 || (shift.y == 0 && shift.x < 0) ? 1 : -1);
newRightConstraint.y =
shift.y + (shift.x < 0 || (shift.x == 0 && shift.y < 0) ? 1 : -1);
if (cross(newRightConstraint, rightConstraint) <= 0)
rightConstraint = newRightConstraint;
}
// Imposing strict constraint for ambiguous turnings, to ensure polygons'
// weak simpleness.
// Has to be defined *outside* abs checks.
if (path[nextCorners[j]].getAmbiguous()) {
if (path[nextCorners[j]].getAmbiguous() == RawBorderPoint::left)
rightConstraint = shift;
else
leftConstraint = shift;
}
j = nextCorners[j];
}
// At this point, constraints are violated by the next corner.
// Then, search for the last k between j and corners[j] not violating them.
direction = convert(normalize(convert(jNextPoint - jPoint)));
k = (j +
cross(jPoint - iPoint, violatedConstraint) /
cross(violatedConstraint, direction)) %
n;
foundK:
kVector[i] = k;
}
return kVector;
}
//--------------------------------------------------------------------------
// Now find the effective intervals inside which we can define possible
// arcs approximating the given raw border:
// for every a in [i,res[i]], the arc connecting border[i] and
// border[a] will be a possible one.
inline std::unique_ptr<int[]> calculateForwardArcs(RawBorder &border,
bool ambiguitiesCheck) {
int const n = (int)border.size();
std::unique_ptr<int[]> nextCorners;
std::unique_ptr<int[]> k = furthestKs(border, nextCorners);
std::unique_ptr<int[]> K(new int[n]);
std::unique_ptr<int[]> res(new int[n]);
// find K[i]= min {k[j]}, j=i..n-1.
for (int i = 0; i < n; ++i) {
int j;
for (j = i, K[i] = k[i]; isCircular(i, j, K[i]); j = (j + 1) % n)
if (isCircular(j, k[j], K[i])) K[i] = k[j];
}
// Finally, we perform the following clean-up operations:
// first, extremities of [i,K[i]] are clipped away, to obtain a
// smoother optimal polygon (and deal with cases like the unitary
// square);
// second, arcs of the kind [i,j] with j<i, become [i,n].
for (int i = n - 1, j = 0; j < n; i = j, ++j) {
res[j] = K[i] < j ? (K[i] == 0 ? n - 1 : n) : K[i] - 1;
}
// Amibiguities check for vertex and edge superpositions. Prevent problems in
// the forecoming
// straight-skeleton thinning process.
if (ambiguitiesCheck) {
for (int i = 1; nextCorners[i] > 0; i = nextCorners[i]) {
if (border[i].getAmbiguous() == RawBorderPoint::right) {
// Check vertices from i (excluded) to res[res[i]]; if in it there
// exists vertex k so that pos(k)==pos(i)...
// This prevents the existence of 0 degree angles in the optimal
// polygon.
int rrPlus1 = (res[res[i] % n] + 1) % n;
for (int j = nextCorners[i];
isCircular(i, j, rrPlus1) &&
j != i; // remember that isCircular(a,a,b) == 1 ...
j = nextCorners[j]) {
if (border[j].getAmbiguous() &&
(border[j].pos() == border[i].pos())) {
res[res[i] % n] = j - 1;
assert((res[i] % n) != j - 1);
// Further, ensure res is increasing
for (int k = res[i] % n; res[k] >= j - 1 && k >= 0; --k) {
res[k] = j - 1;
assert(k != j - 1);
}
break;
}
}
}
}
}
return res;
}
//--------------------------------------------------------------------------
// Let sum[i] and sum2[i] be respectively the sums of vertex coordinates
// from 0 to i, and the sums of their squares; sumsMix contain sums of
// xy terms.
inline void calculateSums(RawBorder &path) {
unsigned int i, n = path.size();
TPointD currentRelativePos;
path.sums() = new TPointD[n + 1];
path.sums2() = new TPointD[n + 1];
path.sumsMix() = new double[n + 1];
path.sums()[0].x = path.sums()[0].y = path.sums2()[0].x = path.sums2()[0].y =
0;
for (i = 1; i < path.size(); ++i) {
currentRelativePos = convert(path[i].pos() - path[0].pos());
path.sums()[i] = path.sums()[i - 1] + currentRelativePos;
path.sums2()[i].x =
path.sums2()[i - 1].x + currentRelativePos.x * currentRelativePos.x;
path.sums2()[i].y =
path.sums2()[i - 1].y + currentRelativePos.y * currentRelativePos.y;
path.sumsMix()[i] =
path.sumsMix()[i - 1] + currentRelativePos.x * currentRelativePos.y;
}
// path[n] is virtually intended as path[0], but we prefer to introduce
// it in the optimality algorithm's count
path.sums()[n].x = path.sums()[n].y = path.sums2()[n].x = path.sums2()[n].y =
0;
}
//--------------------------------------------------------------------------
// Let a,b the index-extremities of an arc of this path.
// Then return its penalty.
inline double penalty(RawBorder &path, int a, int b) {
int n = b - a + 1;
TPointD v =
convert(rotate90(path[b == path.size() ? 0 : b].pos() - path[a].pos()));
TPointD sum = path.sums()[b] - path.sums()[a];
TPointD sum2 = path.sums2()[b] - path.sums2()[a];
double sumMix = path.sumsMix()[b] - path.sumsMix()[a];
double F1 = sum2.x - 2 * sum.x * path[a].x() + n * path[a].x() * path[a].x();
double F2 = sum2.y - 2 * sum.y * path[a].y() + n * path[a].y() * path[a].y();
double F3 = sumMix - sum.x * path[a].y() - sum.y * path[a].x() +
n * path[a].x() * path[a].y();
return sqrt((v.y * v.y * F1 + v.x * v.x * F2 - 2 * v.x * v.y * F3) / n);
}
//--------------------------------------------------------------------------
// NOTA: Il seguente algoritmo di riduzione assicura la semplicita' (debole) dei
// poligoni prodotti.
//
inline void reduceBorder(RawBorder &border, Contour &res,
bool ambiguitiesCheck) {
int n = border.size();
int minPenaltyNext;
std::unique_ptr<int[]> minPenaltyNextArray(new int[n]);
// Calculate preliminary infos
std::unique_ptr<int[]> longestArcFrom =
calculateForwardArcs(border, ambiguitiesCheck);
calculateSums(border);
std::unique_ptr<double[]> penaltyToEnd(new double[n + 1]);
// EXPLANATION:
// The fastest way to extract the optimal reduced border is based on the
// weakly monotonic property of longestArc[].
// The minimal number of its vertices 'm' is easily found by
// traversing the path with the longest step allowed. Let b[] be that
// succession; then, given res[i], it has to be reached by a vertex in
// the interval: {a[i-1], .. , b[i-1]}, where longestArc[a[i-1]]=a[i],
// longestArc[a[i-1]-1]<a[i], and a[m]=n.
// Calculate m
int m = 0;
for (int i = 0; i < n; i = longestArcFrom[i]) ++m;
// Calculate b[]
std::unique_ptr<int[]> b(new int[m + 1]);
b[m] = n;
for (int i = 0, j = 0; j < m; i = longestArcFrom[i], ++j) b[j] = i;
// NOTE: a[] need not be completely found - we just remember the
// a=a[j+1] currently needed.
// Now, build the optimal polygon
for (int j = m - 1, a = n; j >= 0; --j) {
int k;
for (k = b[j]; k >= 0 && longestArcFrom[k] >= a; --k) {
penaltyToEnd[k] = infinity;
for (int i = a; i <= longestArcFrom[k]; ++i) {
double newPenalty = penaltyToEnd[i] + penalty(border, k, i);
if (newPenalty < penaltyToEnd[k]) penaltyToEnd[k] = newPenalty;
minPenaltyNext = i;
}
minPenaltyNextArray[k] = minPenaltyNext;
}
a = k + 1;
}
// Finally, read off the optimal polygon
res.resize(m);
for (int i = 0, j = 0; i < n; i = minPenaltyNextArray[i], ++j) {
res[j] = ContourNode(border[i].x(), border[i].y());
// Ambiguities are still remembered in the output polygon.
if (border[i].getAmbiguous() == RawBorderPoint::left)
res[j].setAttribute(ContourNode::AMBIGUOUS_LEFT);
if (border[i].getAmbiguous() == RawBorderPoint::right)
res[j].setAttribute(ContourNode::AMBIGUOUS_RIGHT);
}
delete[] border.sums();
delete[] border.sums2();
delete[] border.sumsMix();
}
//--------------------------------------------------------------------------
// Reduction caller and list copier.
inline void reduceBorders(BorderList &borders, Contours &result,
bool ambiguitiesCheck) {
unsigned int i, j;
// Initialize output container
result.resize(borders.size());
// Copy results
for (i = 0; i < borders.size(); ++i) {
result[i].resize(borders[i].size());
for (j = 0; j < borders[i].size(); ++j) {
reduceBorder(*borders[i][j], result[i][j], ambiguitiesCheck);
delete borders[i][j];
}
}
}
//--------------------------------------------------------------------------
//===========================
// Polygonization Main
//===========================
// Extracts a polygonal, minimal yet faithful representation of image contours
// Contours* polygonize(const TRasterP &ras){
void polygonize(const TRasterP &ras, Contours &polygons,
VectorizerCoreGlobals &g) {
BorderList *borders;
borders = extractBorders(ras, g.currConfig->m_threshold,
g.currConfig->m_despeckling);
reduceBorders(*borders, polygons, g.currConfig->m_maxThickness > 0.0);
}