dust3d/src/regionfiller.cpp

1497 lines
55 KiB
C++
Raw Blame History

This file contains ambiguous Unicode characters!

This file contains ambiguous Unicode characters that may be confused with others in your current locale. If your use case is intentional and legitimate, you can safely ignore this warning. Use the Escape button to highlight these characters.

#include <numeric>
#include <cmath>
#include <QDebug>
#include <set>
#include "regionfiller.h"
//
// This is an implementation of the following paper
//
// <Filling N-Sided Regions by Quad Meshes for Subdivision Surfaces>
// A. Nasri1, M. Sabin2 and Z. Yasseen1
//
#define isEven(n) (0 == (n) % 2)
#define isOdd(n) (!isEven(n))
RegionFiller::RegionFiller(const std::vector<Node> *vertices,
const std::vector<std::vector<size_t>> *polylines) :
m_sourceVertices(vertices),
m_sourcePolylines(new std::vector<std::vector<size_t>>(*polylines))
{
}
RegionFiller::~RegionFiller()
{
delete m_sourcePolylines;
}
bool RegionFiller::resolveEvenSidedEvenSumOfSegmentsAndBothL2L2AreEven()
{
return resolveOddSidedEvenSumOfSegments(m_sideNum / 2);
}
std::vector<size_t> RegionFiller::createPointsToMapBetween(size_t fromIndex,
size_t toIndex,
size_t segments,
std::map<std::pair<size_t, size_t>, std::vector<size_t>> *map)
{
size_t a = fromIndex;
size_t b = toIndex;
bool needReverse = false;
if (a > b) {
std::swap(a, b);
needReverse = true;
}
auto key = std::make_pair(a, b);
auto findPoints = map->find(key);
if (findPoints != map->end()) {
return needReverse ? reversed(findPoints->second) : findPoints->second;
}
std::vector<size_t> newPointIndices;
createPointsBetween(a, b, segments, &newPointIndices);
map->insert({key, newPointIndices});
return needReverse ? reversed(newPointIndices) : newPointIndices;
}
void RegionFiller::createPointsBetween(size_t fromIndex,
size_t toIndex,
size_t segments,
std::vector<size_t> *newPointIndices)
{
if (fromIndex == toIndex || 0 == segments)
return;
auto fromVertex = m_oldAndNewVertices[fromIndex];
auto toVertex = m_oldAndNewVertices[toIndex];
float stepFactor = 1.0 / segments;
newPointIndices->push_back(fromIndex);
for (size_t i = 1; i <= segments - 1; ++i) {
float factor = stepFactor * i;
RegionFiller::Node node;
node.position = fromVertex.position * (1.0 - factor) + toVertex.position * factor;
node.radius = fromVertex.radius * (1.0 - factor) + toVertex.radius * factor;
if (m_centerSources.find(toVertex.source) != m_centerSources.end()) {
node.source = fromVertex.source;
} else if (m_centerSources.find(fromVertex.source) != m_centerSources.end()) {
node.source = toVertex.source;
} else {
if (factor > 0.5) {
node.source = toVertex.source;
} else {
node.source = fromVertex.source;
}
}
size_t newIndex = m_oldAndNewVertices.size();
m_oldAndNewVertices.push_back(node);
newPointIndices->push_back(newIndex);
}
newPointIndices->push_back(toIndex);
}
std::vector<size_t> RegionFiller::createPointsBetween(size_t fromIndex,
size_t toIndex,
size_t segments)
{
std::vector<size_t> newPointIndices;
createPointsBetween(fromIndex, toIndex, segments, &newPointIndices);
return newPointIndices;
}
void RegionFiller::collectEdgePoints(size_t polyline, int startPos, int stopPos,
std::vector<size_t> *pointIndices)
{
const auto &edgeIndices = (*m_sourcePolylines)[polyline];
int step = startPos < stopPos ? 1 : -1;
int totalSteps = std::abs(startPos - stopPos) + 1;
int current = startPos;
for (int i = 0; i < totalSteps; ++i) {
pointIndices->push_back(edgeIndices[current]);
current += step;
}
}
std::vector<size_t> RegionFiller::collectEdgePoints(size_t polyline, int startPos, int stopPos)
{
std::vector<size_t> pointIndices;
collectEdgePoints(polyline, startPos, stopPos, &pointIndices);
return pointIndices;
}
std::vector<size_t> RegionFiller::reversed(const std::vector<size_t> &pointIndices)
{
std::vector<size_t> newPointIndices = pointIndices;
std::reverse(newPointIndices.begin(), newPointIndices.end());
return newPointIndices;
}
float RegionFiller::averageRadius(size_t *maxNodeIndex)
{
float sumOfRadius = 0;
size_t num = 0;
float maxRadius = 0;
for (const auto &polyline: *m_sourcePolylines) {
for (const auto &it: polyline) {
const auto &vertex = (*m_sourceVertices)[it];
if (vertex.radius >= maxRadius) {
maxRadius = vertex.radius;
if (nullptr != maxNodeIndex)
*maxNodeIndex = it;
}
sumOfRadius += vertex.radius;
++num;
}
}
if (0 == num)
return 0;
auto radius = sumOfRadius / num;
return radius;
}
bool RegionFiller::resolveOddSidedEvenSumOfSegments(int siUpperBound)
{
auto si = [&](int i) {
int sum = 0;
for (int j = 1; j <= siUpperBound; ++j) {
sum += std::pow(-1, j + 1) *
m_sideSegmentNums[(i + 2 * j - 1) % m_sideSegmentNums.size()];
}
return sum;
};
auto di = [&](int i) {
return std::ceil(si(i) * 0.5);
};
std::vector<int> gRaySegmentNums(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
gRaySegmentNums[i] = di(i);
if (gRaySegmentNums[i] <= 0) {
qDebug() << "resolveOddSidedEvenSumOfSegments failed, di:" << gRaySegmentNums[i];
return false;
}
}
std::vector<int> eIndices(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
auto length = gRaySegmentNums[(i + 1) % gRaySegmentNums.size()];
eIndices[i] = m_sideSegmentNums[i] - length;
}
std::set<size_t> cornerVertices;
std::vector<std::tuple<int, int, int>> edgeSegments; // tuple<begin, end, edgeIndex>
for (int i = 0; i < m_sideNum; ++i) {
if (gRaySegmentNums[i] <= 0) {
edgeSegments.push_back(std::make_tuple(0, (*m_sourcePolylines)[i].size() - 1, i));
continue;
}
const auto &pointEdgeIndex = eIndices[i];
cornerVertices.insert((*m_sourcePolylines)[i][pointEdgeIndex]);
if (pointEdgeIndex > 0)
edgeSegments.push_back(std::make_tuple(0, pointEdgeIndex, i));
if (pointEdgeIndex < (int)(*m_sourcePolylines)[i].size() - 1)
edgeSegments.push_back(std::make_tuple(pointEdgeIndex, (*m_sourcePolylines)[i].size() - 1, i));
}
std::vector<std::pair<QVector3D, float>> gPositionAndAreas;
for (size_t i = 0; i < edgeSegments.size(); ++i) {
size_t j = (i + 1) % edgeSegments.size();
size_t m = (j + 2) % edgeSegments.size();
const auto &segmentI = edgeSegments[i];
const auto &segmentJ = edgeSegments[j];
const auto &segmentM = edgeSegments[m];
std::vector<size_t> firstSegmentPoints;
std::vector<size_t> secondSegmentPoints;
std::vector<size_t> forthSegmentPoints;
collectEdgePoints(std::get<2>(segmentI), std::get<0>(segmentI), std::get<1>(segmentI),
&firstSegmentPoints);
collectEdgePoints(std::get<2>(segmentJ), std::get<0>(segmentJ), std::get<1>(segmentJ),
&secondSegmentPoints);
collectEdgePoints(std::get<2>(segmentM), std::get<0>(segmentM), std::get<1>(segmentM),
&forthSegmentPoints);
if (cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)]) != cornerVertices.end() &&
cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentI)][std::get<1>(segmentI)]) != cornerVertices.end()) {
const auto &p1 = (*m_sourceVertices)[(*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)]];
const auto &p2 = (*m_sourceVertices)[(*m_sourcePolylines)[std::get<2>(segmentI)][std::get<1>(segmentI)]];
gPositionAndAreas.push_back(std::make_pair((p1.position + p2.position) * 0.5,
(firstSegmentPoints.size() - 1) * (secondSegmentPoints.size() - 1) * 0.5));
continue;
}
if (cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)]) != cornerVertices.end() &&
cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentJ)][std::get<1>(segmentJ)]) != cornerVertices.end()) {
++i;
const auto &p1 = (*m_sourceVertices)[(*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)]];
const auto &p2 = (*m_sourceVertices)[(*m_sourcePolylines)[std::get<2>(segmentI)][std::get<1>(segmentI)]];
const auto &p3 = (*m_sourceVertices)[(*m_sourcePolylines)[std::get<2>(segmentJ)][std::get<1>(segmentJ)]];
gPositionAndAreas.push_back(std::make_pair(p1.position + p2.position - p3.position,
(firstSegmentPoints.size() - 1) * (secondSegmentPoints.size() - 1)));
}
}
QVector3D gTop;
float gBottom = 0;
for (const auto &it: gPositionAndAreas) {
gTop += it.first / it.second;
gBottom += 1.0 / it.second;
}
auto gPosition = gTop / gBottom;
size_t gIndex = m_oldAndNewVertices.size();
Node node;
node.position = gPosition;
node.radius = averageRadius(&node.source);
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
std::map<std::pair<size_t, size_t>, std::vector<size_t>> map;
for (size_t i = 0; i < edgeSegments.size(); ++i) {
size_t j = (i + 1) % edgeSegments.size();
size_t m = (j + 2) % edgeSegments.size();
const auto &segmentI = edgeSegments[i];
const auto &segmentJ = edgeSegments[j];
const auto &segmentM = edgeSegments[m];
std::vector<size_t> firstSegmentPoints;
std::vector<size_t> secondSegmentPoints;
std::vector<size_t> forthSegmentPoints;
collectEdgePoints(std::get<2>(segmentI), std::get<0>(segmentI), std::get<1>(segmentI),
&firstSegmentPoints);
collectEdgePoints(std::get<2>(segmentJ), std::get<0>(segmentJ), std::get<1>(segmentJ),
&secondSegmentPoints);
collectEdgePoints(std::get<2>(segmentM), std::get<0>(segmentM), std::get<1>(segmentM),
&forthSegmentPoints);
if (cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)]) != cornerVertices.end() &&
cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentI)][std::get<1>(segmentI)]) != cornerVertices.end()) {
std::vector<std::vector<size_t>> region = {
firstSegmentPoints,
createPointsToMapBetween((*m_sourcePolylines)[std::get<2>(segmentI)][std::get<1>(segmentI)],
gIndex, forthSegmentPoints.size() - 1, &map),
createPointsToMapBetween(gIndex,
(*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)], secondSegmentPoints.size() - 1, &map)
};
m_newRegions.push_back(region);
continue;
}
if (cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)]) != cornerVertices.end() &&
cornerVertices.find((*m_sourcePolylines)[std::get<2>(segmentJ)][std::get<1>(segmentJ)]) != cornerVertices.end()) {
++i;
std::vector<std::vector<size_t>> region = {
firstSegmentPoints,
secondSegmentPoints,
createPointsToMapBetween((*m_sourcePolylines)[std::get<2>(segmentJ)][std::get<1>(segmentJ)],
gIndex, firstSegmentPoints.size() - 1, &map),
createPointsToMapBetween(gIndex,
(*m_sourcePolylines)[std::get<2>(segmentI)][std::get<0>(segmentI)], secondSegmentPoints.size() - 1, &map)
};
m_newRegions.push_back(region);
}
}
return true;
}
bool RegionFiller::resolveEvenSidedEvenSumOfSegmentsAndBothL2L2AreOdd()
{
return resolveOddSidedOddSumOfSegments(m_sideNum / 2);
}
bool RegionFiller::resolveOddSidedOddSumOfSegments(int siUpperBound)
{
auto si = [&](int i) {
int sum = 0;
for (int j = 1; j <= siUpperBound; ++j) {
sum += std::pow(-1, j + 1) *
(m_sideSegmentNums[(i + 2 * j - 1) % m_sideSegmentNums.size()] - 1);
}
return sum;
};
auto di = [&](int i) {
return std::round(si(i) * 0.5);
};
std::vector<int> gRaySegmentNums(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
gRaySegmentNums[i] = di(i);
if (gRaySegmentNums[i] <= 0) {
qDebug() << "resolveOddSidedOddSumOfSegments failed, di:" << gRaySegmentNums[i];
return false;
}
}
std::vector<int> areas(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
areas[i] = gRaySegmentNums[i] * gRaySegmentNums[(i + 1) % m_sideNum];
}
std::vector<std::vector<int>> eIndices(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
auto length = gRaySegmentNums[(i + 1) % gRaySegmentNums.size()];
eIndices[i] = std::vector<int>{(int)m_sideSegmentNums[i] - length - 1,
(int)m_sideSegmentNums[i] - length
};
if (eIndices[i][0] < 0 || eIndices[i][1] < 0)
return false;
}
std::vector<float> eEdgeLengths(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
eEdgeLengths[i] = ((*m_sourceVertices)[(*m_sourcePolylines)[i][eIndices[i][0]]].position -
(*m_sourceVertices)[(*m_sourcePolylines)[i][eIndices[i][1]]].position).length();
}
float averageDislocationEdgeLength = std::accumulate(eEdgeLengths.begin(), eEdgeLengths.end(), 0.0) / eEdgeLengths.size();
auto gi = [&](int i) {
int j = (i + m_sideNum - 1) % m_sideNum;
return (*m_sourceVertices)[(*m_sourcePolylines)[i][eIndices[i][0]]].position +
(*m_sourceVertices)[(*m_sourcePolylines)[j][eIndices[j][1]]].position -
(*m_sourceVertices)[(*m_sourcePolylines)[i][0]].position;
};
auto gOffsetTargets = [&](int i) {
int h = (i + m_sideNum - 1) % m_sideNum;
return ((*m_sourceVertices)[(*m_sourcePolylines)[i][eIndices[i][0]]].position +
(*m_sourceVertices)[(*m_sourcePolylines)[h][eIndices[h][1]]].position) * 0.5;
};
QVector3D gTop;
float gBottom = 0;
for (int i = 0; i < m_sideNum; ++i) {
gTop += gi(i) / areas[i];
gBottom += 1.0 / areas[i];
}
auto gPosition = gTop / gBottom;
size_t gSource = 0;
auto gRadius = averageRadius(&gSource);
m_centerSources.insert(gSource);
std::vector<std::vector<std::vector<size_t>>> segmentsWithG(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
segmentsWithG[i].resize(2);
}
std::vector<size_t> gFace(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
int j = (i + 1) % m_sideNum;
auto offsetTarget = gOffsetTargets(j);
auto gPositionPerEdge = gPosition + (offsetTarget - gPosition).normalized() * averageDislocationEdgeLength * 0.7;
size_t gIndex = m_oldAndNewVertices.size();
gFace[j] = gIndex;
Node node;
node.position = gPositionPerEdge;
node.radius = gRadius;
node.source = gSource;
m_oldAndNewVertices.push_back(node);
createPointsBetween(gIndex, (*m_sourcePolylines)[i][eIndices[i][1]], gRaySegmentNums[i], &segmentsWithG[i][1]);
createPointsBetween(gIndex, (*m_sourcePolylines)[j][eIndices[j][0]], gRaySegmentNums[j], &segmentsWithG[j][0]);
}
m_newFaces.push_back(gFace);
std::vector<std::vector<std::vector<size_t>>> reversedSegmentsWithG(m_sideNum);
reversedSegmentsWithG = segmentsWithG;
for (int i = 0; i < m_sideNum; ++i) {
for (int x = 0; x < 2; ++x)
std::reverse(reversedSegmentsWithG[i][x].begin(), reversedSegmentsWithG[i][x].end());
}
for (int i = 0; i < m_sideNum; ++i) {
int j = (i + 1) % m_sideNum;
std::vector<size_t> e0c1;
collectEdgePoints(i, eIndices[i][1], (*m_sourcePolylines)[i].size() - 1, &e0c1);
std::vector<size_t> c1e1;
collectEdgePoints(j, 0, eIndices[j][0], &c1e1);
std::vector<std::vector<size_t>> region = {
e0c1,
c1e1,
reversedSegmentsWithG[j][0],
segmentsWithG[i][1]
};
m_newRegions.push_back(region);
}
// rectangle slices
for (int i = 0; i < m_sideNum; ++i) {
int j = (i + 1) % m_sideNum;
std::vector<size_t> e1e2 = {
(*m_sourcePolylines)[i][eIndices[i][0]],
(*m_sourcePolylines)[i][eIndices[i][1]]
};
std::vector<size_t> gCoreEdge = {
gFace[j],
gFace[i]
};
std::vector<std::vector<size_t>> region = {
e1e2,
reversedSegmentsWithG[i][1],
gCoreEdge,
segmentsWithG[i][0],
};
m_newRegions.push_back(region);
}
return true;
}
bool RegionFiller::resolveEvenSideEvenSumOfSegmentsAndL1IsOddL2IsEven()
{
std::rotate(m_sourcePolylines->begin(), m_sourcePolylines->begin() + 1, m_sourcePolylines->end());
prepare();
return resolveEvenSideEvenSumOfSegmentsAndL1IsEvenL2IsOdd();
}
bool RegionFiller::resolveEvenSideEvenSumOfSegmentsAndL1IsEvenL2IsOdd()
{
auto modifiedSideSegmentNums = m_sideSegmentNums;
for (int i = 0; i < m_sideNum; i += 2) {
modifiedSideSegmentNums[i] -= 1;
}
auto si = [&](int i) {
int sum = 0;
int halfSideNum = m_sideNum / 2;
for (int j = 1; j <= halfSideNum; ++j) {
sum += std::pow(-1, j + 1) *
(modifiedSideSegmentNums[(i + 2 * j - 1) % m_sideSegmentNums.size()]);
}
return sum;
};
auto di = [&](int i) {
return std::round(si(i) * 0.5);
};
std::vector<int> gRaySegmentNums(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
gRaySegmentNums[i] = di(i);
if (gRaySegmentNums[i] <= 0) {
qDebug() << "resolveEvenSideEvenSumOfSegmentsAndL1IsEvenL2IsOdd failed, di:" << gRaySegmentNums[i];
return false;
}
}
std::vector<int> areas(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
areas[i] = gRaySegmentNums[i] * gRaySegmentNums[(i + 1) % m_sideNum];
}
std::vector<std::vector<int>> eIndices(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
auto length = gRaySegmentNums[(i + gRaySegmentNums.size() - 1) % gRaySegmentNums.size()];
if (0 == i % 2) {
eIndices[i] = std::vector<int>{(int)length,
(int)length + 1
};
} else {
eIndices[i] = std::vector<int>{(int)length,
(int)length
};
}
if (eIndices[i][0] < 0 || eIndices[i][1] < 0) {
qDebug() << "resolveEvenSideEvenSumOfSegmentsAndL1IsEvenL2IsOdd failed";
return false;
}
}
std::vector<float> eEdgeLengths(m_sideNum, (float)0.0);
int validEdgeNum = 0;
for (int i = 0; i < m_sideNum; ++i) {
const auto &firstIndex = eIndices[i][0];
const auto &secondIndex = eIndices[i][1];
if (firstIndex == secondIndex)
continue;
validEdgeNum++;
eEdgeLengths[i] = ((*m_sourceVertices)[(*m_sourcePolylines)[i][firstIndex]].position -
(*m_sourceVertices)[(*m_sourcePolylines)[i][secondIndex]].position).length();
}
float averageDislocationEdgeLength = std::accumulate(eEdgeLengths.begin(), eEdgeLengths.end(), 0.0) / validEdgeNum;
auto gi = [&](int i) {
int h = (i + m_sideNum - 1) % m_sideNum;
return (*m_sourceVertices)[(*m_sourcePolylines)[i][eIndices[i][0]]].position +
(*m_sourceVertices)[(*m_sourcePolylines)[h][eIndices[h][1]]].position -
(*m_sourceVertices)[(*m_sourcePolylines)[i][0]].position;
};
auto gOffsetTargets = [&](int i) {
int h = (i + m_sideNum - 1) % m_sideNum;
return ((*m_sourceVertices)[(*m_sourcePolylines)[i][eIndices[i][0]]].position +
(*m_sourceVertices)[(*m_sourcePolylines)[h][eIndices[h][1]]].position) * 0.5;
};
QVector3D gTop;
float gBottom = 0;
for (int i = 0; i < m_sideNum; ++i) {
gTop += gi(i) / areas[i];
gBottom += 1.0 / areas[i];
}
auto gPosition = gTop / gBottom;
size_t gSource = 0;
auto gRadius = averageRadius(&gSource);
std::vector<size_t> gFace(m_sideNum / 2);
for (size_t k = 0; k < gFace.size(); ++k) {
int i = k * 2;
int h = (i + m_sideNum - 1) % m_sideNum;
auto offsetTarget = gOffsetTargets(i) * 0.5 + gOffsetTargets(h) * 0.5;
size_t gIndex = m_oldAndNewVertices.size();
auto gPositionPerEdge = gPosition + (offsetTarget - gPosition).normalized() * averageDislocationEdgeLength * 0.7;
Node node;
node.position = gPositionPerEdge;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
gFace[k] = gIndex;
}
m_newFaces.push_back(gFace);
std::vector<std::vector<std::vector<size_t>>> segmentsWithG(m_sideNum);
for (int i = 0; i < m_sideNum; ++i) {
segmentsWithG[i].resize(2);
}
for (int i = 0; i < m_sideNum; ++i) {
int k = i / 2;
const auto &g0Index = gFace[k % gFace.size()];
const auto &g1Index = gFace[(k + 1) % gFace.size()];
if (eIndices[i][0] == eIndices[i][1]) {
createPointsBetween(g1Index, (*m_sourcePolylines)[i][eIndices[i][0]], gRaySegmentNums[i],
&segmentsWithG[i][0]);
segmentsWithG[i][1] = segmentsWithG[i][0];
} else {
createPointsBetween(g0Index, (*m_sourcePolylines)[i][eIndices[i][0]], gRaySegmentNums[i],
&segmentsWithG[i][0]);
createPointsBetween(g1Index, (*m_sourcePolylines)[i][eIndices[i][1]], gRaySegmentNums[i],
&segmentsWithG[i][1]);
}
}
std::vector<std::vector<std::vector<size_t>>> reversedSegmentsWithG(m_sideNum);
reversedSegmentsWithG = segmentsWithG;
for (int i = 0; i < m_sideNum; ++i) {
for (int x = 0; x < 2; ++x)
std::reverse(reversedSegmentsWithG[i][x].begin(), reversedSegmentsWithG[i][x].end());
}
for (int i = 0; i < m_sideNum; ++i) {
int j = (i + 1) % m_sideNum;
std::vector<size_t> e0c1;
collectEdgePoints(i, eIndices[i][1], (*m_sourcePolylines)[i].size() - 1, &e0c1);
std::vector<size_t> c1e1;
collectEdgePoints(j, 0, eIndices[j][0], &c1e1);
std::vector<std::vector<size_t>> region = {
e0c1, c1e1, reversedSegmentsWithG[j][0], segmentsWithG[i][1]
};
m_newRegions.push_back(region);
}
// Rectangle slices
for (size_t k = 0; k < gFace.size(); ++k) {
int i = k * 2;
std::vector<size_t> e1e2 = {
(*m_sourcePolylines)[i][eIndices[i][0]],
(*m_sourcePolylines)[i][eIndices[i][1]]
};
std::vector<size_t> g1g0 = {
gFace[(k + 1) % gFace.size()],
gFace[k % gFace.size()]
};
std::vector<std::vector<size_t>> region = {
g1g0, segmentsWithG[i][0], e1e2, reversedSegmentsWithG[i][1]
};
m_newRegions.push_back(region);
}
return true;
}
bool RegionFiller::resolveQuadrilateralRegionMismatchSimilarCase(int m, int n, int p, int q, bool pqSwapped)
{
return resolveQuadrilateralRegionWithIntegerSolution(m, n, p, q, pqSwapped);
}
bool RegionFiller::resolveQuadrilateralRegionWithNonIntegerSolution(int m, int n, int p, int q, bool pqSwapped)
{
float aPlusB = m_sideSegmentNums[m] - m_sideSegmentNums[n];
float aFloat = ((m_sideSegmentNums[m] - m_sideSegmentNums[n]) + (m_sideSegmentNums[p] - m_sideSegmentNums[q])) / 2.0;
float bFloat = (m_sideSegmentNums[m] - m_sideSegmentNums[n]) - aFloat;
float a = std::ceil(aFloat);
float b = std::ceil(bFloat);
int c0EdgeIndex = pqSwapped ? (*m_sourcePolylines)[m].size() - 1 : 0;
int c1EdgeIndex = pqSwapped ? (*m_sourcePolylines)[p].size() - 1 : 0;
int c2EdgeIndex = pqSwapped ? (*m_sourcePolylines)[n].size() - 1 : 0;
int c3EdgeIndex = pqSwapped ? (*m_sourcePolylines)[q].size() - 1 : 0;
int c0AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[q].size() - 1;
int c1AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[m].size() - 1;
int c2AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[p].size() - 1;
int c3AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[n].size() - 1;
const auto &c0Index = (*m_sourcePolylines)[m][c0EdgeIndex];
const auto &c1Index = (*m_sourcePolylines)[p][c1EdgeIndex];
const auto &c0 = (*m_sourceVertices)[c0Index].position;
const auto &c1 = (*m_sourceVertices)[c1Index].position;
int f0 = m_sideSegmentNums[n] / 2;
int f1 = m_sideSegmentNums[n] - f0;
int f2 = std::ceil((m_sideSegmentNums[p] - a) / 2.0);
int f3 = m_sideSegmentNums[p] - a - f2;
if (f3 < 0) {
qDebug() << "resolveQuadrilateralRegionWithNonIntegerSolution failed, f3:" << f3;
return false;
}
int e0EdgeIndex = f0;
if (pqSwapped)
e0EdgeIndex = (*m_sourcePolylines)[m].size() - 1 - e0EdgeIndex;
const auto &e0Index = (*m_sourcePolylines)[m][e0EdgeIndex];
int e2EdgeIndex = f0 + aPlusB;
if (pqSwapped)
e2EdgeIndex = (*m_sourcePolylines)[m].size() - 1 - e2EdgeIndex;
const auto &e2Index = (*m_sourcePolylines)[m][e2EdgeIndex];
int e3EdgeIndex = f2;
if (pqSwapped)
e3EdgeIndex = (*m_sourcePolylines)[p].size() - 1 - e3EdgeIndex;
const auto &e3Index = (*m_sourcePolylines)[p][e3EdgeIndex];
int e4EdgeIndex = (*m_sourcePolylines)[p].size() - 1 - f3;
if (pqSwapped)
e4EdgeIndex = (*m_sourcePolylines)[p].size() - 1 - e4EdgeIndex;
const auto &e4Index = (*m_sourcePolylines)[p][e4EdgeIndex];
int e5EdgeIndex = f3;
if (pqSwapped)
e5EdgeIndex = (*m_sourcePolylines)[q].size() - 1 - e5EdgeIndex;
const auto &e5Index = (*m_sourcePolylines)[q][e5EdgeIndex];
int e6EdgeIndex = (*m_sourcePolylines)[q].size() - 1 - f2;
if (pqSwapped)
e6EdgeIndex = (*m_sourcePolylines)[q].size() - 1 - e6EdgeIndex;
const auto &e6Index = (*m_sourcePolylines)[q][e6EdgeIndex];
int e7EdgeIndex = f1;
if (pqSwapped)
e7EdgeIndex = (*m_sourcePolylines)[n].size() - 1 - e7EdgeIndex;
const auto &e7Index = (*m_sourcePolylines)[n][e7EdgeIndex];
const auto &e0 = (*m_sourceVertices)[e0Index].position;
const auto &e2 = (*m_sourceVertices)[e2Index].position;
const auto &e3 = (*m_sourceVertices)[e3Index].position;
const auto &e4 = (*m_sourceVertices)[e4Index].position;
const auto &e5 = (*m_sourceVertices)[e5Index].position;
const auto &e6 = (*m_sourceVertices)[e6Index].position;
const auto &e7 = (*m_sourceVertices)[e7Index].position;
auto g1InitialPosition = e0 + e6 - c0;
auto g2InitialPosition = e2 + e3 - c1;
auto gpInitialPosition1 = (e5 + e4) * 0.5;
auto gpInitialPosition2 = gpInitialPosition1;
auto g1a = f0 * f2;
auto g2a = f1 * f2;
auto gp1a = f0 * f3 + f0 * b;
auto gp2a = f1 * f3 + f1 * a;
float weightG1 = 1.0 / g1a;
float weightG2 = 1.0 / g2a;
float weightGp1 = 1.0 / gp1a;
float weightGp2 = 1.0 / gp2a;
size_t gSource = 0;
auto gRadius = averageRadius(&gSource);
auto gPosition = (g1InitialPosition * weightG1 +
g2InitialPosition * weightG2 +
gpInitialPosition1 * weightGp1 +
gpInitialPosition2 * weightGp2) /
(weightG1 + weightG2 + weightGp1 + weightGp2);
auto averageOffsetDistance = e0.distanceToPoint(e2);
auto g1Position = gPosition + (((e0 + e6) * 0.5) - gPosition).normalized() * averageOffsetDistance * 0.7;
auto g2Position = gPosition + (((e2 + e3) * 0.5) - gPosition).normalized() * averageOffsetDistance * 0.7;
auto gpPosition = gPosition + (e7 - gPosition).normalized() * averageOffsetDistance * 0.7;
auto g1Index = m_oldAndNewVertices.size();
{
Node node;
node.position = g1Position;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
}
auto g2Index = m_oldAndNewVertices.size();
{
Node node;
node.position = g2Position;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
}
auto gpIndex = m_oldAndNewVertices.size();
if (0 != f3) {
Node node;
node.position = gpPosition;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
} else {
gpIndex = e7Index;
}
std::map<std::pair<size_t, size_t>, std::vector<size_t>> map;
{
std::vector<std::vector<size_t>> region = {
collectEdgePoints(m, c0EdgeIndex, e0EdgeIndex),
createPointsToMapBetween(e0Index, g1Index, f2, &map),
createPointsToMapBetween(g1Index, e6Index, f0, &map),
collectEdgePoints(q, e6EdgeIndex, c0AlternativeEdgeIndex)
};
m_newRegions.push_back(region);
}
{
std::vector<std::vector<size_t>> region = {
collectEdgePoints(m, e0EdgeIndex, e2EdgeIndex),
createPointsToMapBetween(e2Index, g2Index, f2, &map),
createPointsToMapBetween(g2Index, g1Index, aPlusB, &map),
createPointsToMapBetween(g1Index, e0Index, f2, &map)
};
m_newRegions.push_back(region);
}
{
std::vector<std::vector<size_t>> region = {
collectEdgePoints(m, e2EdgeIndex, c1AlternativeEdgeIndex),
collectEdgePoints(p, c1EdgeIndex, e3EdgeIndex),
createPointsToMapBetween(e3Index, g2Index, f1, &map),
createPointsToMapBetween(g2Index, e2Index, f2, &map)
};
m_newRegions.push_back(region);
}
if (gpIndex != e7Index) {
std::vector<std::vector<size_t>> region = {
createPointsToMapBetween(g2Index, e3Index, f1, &map),
collectEdgePoints(p, e3EdgeIndex, e4EdgeIndex),
createPointsToMapBetween(e4Index, gpIndex, f1, &map),
createPointsToMapBetween(gpIndex, g2Index, a, &map)
};
m_newRegions.push_back(region);
} else {
std::vector<std::vector<size_t>> region = {
createPointsToMapBetween(g2Index, e3Index, f1, &map),
collectEdgePoints(p, e3EdgeIndex, c2AlternativeEdgeIndex),
collectEdgePoints(n, c2EdgeIndex, e7EdgeIndex),
createPointsToMapBetween(gpIndex, g2Index, a, &map)
};
m_newRegions.push_back(region);
}
if (gpIndex != e7Index) {
std::vector<std::vector<size_t>> region = {
collectEdgePoints(p, e4EdgeIndex, c2AlternativeEdgeIndex),
collectEdgePoints(n, c2EdgeIndex, e7EdgeIndex),
createPointsToMapBetween(e7Index, gpIndex, f3, &map),
createPointsToMapBetween(gpIndex, e4Index, f1, &map)
};
m_newRegions.push_back(region);
}
if (gpIndex != e7Index) {
std::vector<std::vector<size_t>> region = {
createPointsToMapBetween(gpIndex, e7Index, f3, &map),
collectEdgePoints(n, e7EdgeIndex, c3AlternativeEdgeIndex),
collectEdgePoints(q, c3EdgeIndex, e5EdgeIndex),
createPointsToMapBetween(e5Index, gpIndex, f0, &map)
};
m_newRegions.push_back(region);
}
if (gpIndex != e7Index) {
std::vector<std::vector<size_t>> region = {
createPointsToMapBetween(g1Index, gpIndex, b, &map),
createPointsToMapBetween(gpIndex, e5Index, f0, &map),
collectEdgePoints(q, e5EdgeIndex, e6EdgeIndex),
createPointsToMapBetween(e6Index, g1Index, f0, &map)
};
m_newRegions.push_back(region);
} else {
std::vector<std::vector<size_t>> region = {
createPointsToMapBetween(g1Index, gpIndex, b, &map),
collectEdgePoints(n, e7EdgeIndex, c3AlternativeEdgeIndex),
collectEdgePoints(q, c3EdgeIndex, e6EdgeIndex),
createPointsToMapBetween(e6Index, g1Index, f0, &map)
};
m_newRegions.push_back(region);
}
{
std::vector<std::vector<size_t>> region = {
createPointsToMapBetween(g1Index, g2Index, aPlusB, &map),
createPointsToMapBetween(g2Index, gpIndex, a, &map),
createPointsToMapBetween(gpIndex, g1Index, b, &map)
};
m_newRegions.push_back(region);
}
return true;
}
bool RegionFiller::resolveQuadrilateralRegionWithIntegerSolution(int m, int n, int p, int q, bool pqSwapped)
{
// a + b = m n
// a b = p q
// (a + b) + (a - b) = (m - n) + (p - q)
// a = ((m - n) + (p - q)) / 2
int a = ((m_sideSegmentNums[m] - m_sideSegmentNums[n]) + (m_sideSegmentNums[p] - m_sideSegmentNums[q])) / 2;
int b = (m_sideSegmentNums[m] - m_sideSegmentNums[n]) - a;
bool isMismatchSimilar = 0 == b;
int f0 = m_sideSegmentNums[n] / 2;
int f1 = m_sideSegmentNums[n] - f0;
int f2 = (m_sideSegmentNums[q] - b) / 2;
int f3 = m_sideSegmentNums[q] - b - f2;
if (f3 < 0) {
qDebug() << "resolveQuadrilateralRegionWithIntegerSolution failed, f3:" << f3;
return false;
}
float a1 = f0 * f2;
float a2 = f0 * b;
float a3 = f0 * f3;
float a4 = f1 * f2;
float a5 = f1 * a;
float a6 = f1 * f3;
float a7 = a * f2;
float a8 = b * f2;
float a9 = a * b;
int c0EdgeIndex = pqSwapped ? (*m_sourcePolylines)[m].size() - 1 : 0;
int c1EdgeIndex = pqSwapped ? (*m_sourcePolylines)[p].size() - 1 : 0;
int c2EdgeIndex = pqSwapped ? (*m_sourcePolylines)[n].size() - 1 : 0;
int c3EdgeIndex = pqSwapped ? (*m_sourcePolylines)[q].size() - 1 : 0;
int c0AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[q].size() - 1;
int c1AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[m].size() - 1;
int c2AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[p].size() - 1;
int c3AlternativeEdgeIndex = pqSwapped ? 0 : (*m_sourcePolylines)[n].size() - 1;
const auto &c0Index = (*m_sourcePolylines)[m][c0EdgeIndex];
const auto &c1Index = (*m_sourcePolylines)[p][c1EdgeIndex];
const auto &c2Index = (*m_sourcePolylines)[n][c2EdgeIndex];
const auto &c3Index = (*m_sourcePolylines)[q][c3EdgeIndex];
const auto &c0 = (*m_sourceVertices)[c0Index].position;
const auto &c1 = (*m_sourceVertices)[c1Index].position;
const auto &c2 = (*m_sourceVertices)[c2Index].position;
const auto &c3 = (*m_sourceVertices)[c3Index].position;
int e0EdgeIndex = f0;
if (pqSwapped)
e0EdgeIndex = (*m_sourcePolylines)[m].size() - 1 - e0EdgeIndex;
const auto &e0Index = (*m_sourcePolylines)[m][e0EdgeIndex];
int e1EdgeIndex = f0 + a;
if (pqSwapped)
e1EdgeIndex = (*m_sourcePolylines)[m].size() - 1 - e1EdgeIndex;
const auto &e1Index = (*m_sourcePolylines)[m][e1EdgeIndex];
int e2EdgeIndex = f0 + a + b;
if (pqSwapped)
e2EdgeIndex = (*m_sourcePolylines)[m].size() - 1 - e2EdgeIndex;
const auto &e2Index = (*m_sourcePolylines)[m][e2EdgeIndex];
int e3EdgeIndex = f2;
if (pqSwapped)
e3EdgeIndex = (*m_sourcePolylines)[p].size() - 1 - e3EdgeIndex;
const auto &e3Index = (*m_sourcePolylines)[p][e3EdgeIndex];
int e4EdgeIndex = f2 + a;
if (pqSwapped)
e4EdgeIndex = (*m_sourcePolylines)[p].size() - 1 - e4EdgeIndex;
const auto &e4Index = (*m_sourcePolylines)[p][e4EdgeIndex];
int e5EdgeIndex = f3;
if (pqSwapped)
e5EdgeIndex = (*m_sourcePolylines)[q].size() - 1 - e5EdgeIndex;
const auto &e5Index = (*m_sourcePolylines)[q][e5EdgeIndex];
int e6EdgeIndex = f3 + b;
if (pqSwapped)
e6EdgeIndex = (*m_sourcePolylines)[q].size() - 1 - e6EdgeIndex;
const auto &e6Index = (*m_sourcePolylines)[q][e6EdgeIndex];
int e7EdgeIndex = f1;
if (pqSwapped)
e7EdgeIndex = (*m_sourcePolylines)[n].size() - 1 - e7EdgeIndex;
const auto &e7Index = (*m_sourcePolylines)[n][e7EdgeIndex];
const auto &e0 = (*m_sourceVertices)[e0Index].position;
const auto &e1 = (*m_sourceVertices)[e1Index].position;
const auto &e2 = (*m_sourceVertices)[e2Index].position;
const auto &e3 = (*m_sourceVertices)[e3Index].position;
const auto &e4 = (*m_sourceVertices)[e4Index].position;
const auto &e5 = (*m_sourceVertices)[e5Index].position;
const auto &e6 = (*m_sourceVertices)[e6Index].position;
const auto &e7 = (*m_sourceVertices)[e7Index].position;
auto g1 = e1 + e5 - c0;
auto gn = e1 + e4 - c1;
auto g2 = e4 + e7 - c2;
auto gp = e5 + e7 - c3;
auto weight1 = 1.0 / (a1 + a2 + a7 + a9);
auto weightn = 1.0 / (a4 + a5 + a8 + a9);
auto weight2 = 1.0 / (a6);
auto weightp = 1.0 / (a3);
auto gPosition = (g1 * weight1 + gn * weightn + g2 * weight2 + gp * weightp) /
(weight1 + weightn + weight2 + weightp);
size_t gSource = 0;
auto gRadius = averageRadius(&gSource);
auto averageOffsetDistance = e0.distanceToPoint(e2);
if (!isMismatchSimilar)
averageOffsetDistance *= 0.5;
auto g1Position = gPosition + (e6 - gPosition).normalized() * averageOffsetDistance * 0.7;
auto gnPosition = gPosition + (e1 - gPosition).normalized() * averageOffsetDistance * 0.7;
auto g2Position = gPosition + (e3 - gPosition).normalized() * averageOffsetDistance * 0.7;
auto gpPosition = gPosition + (e7 - gPosition).normalized() * averageOffsetDistance * 0.7;
if (isMismatchSimilar) {
gpPosition = ((*m_sourceVertices)[e0Index].position +
(*m_sourceVertices)[e5Index].position +
(*m_sourceVertices)[e7Index].position +
(*m_sourceVertices)[e4Index].position) / 4.0;
gnPosition = ((*m_sourceVertices)[e0Index].position +
gpPosition +
(*m_sourceVertices)[e3Index].position +
(*m_sourceVertices)[c1Index].position) / 4.0;
}
size_t g1Index = m_oldAndNewVertices.size();
if (!isMismatchSimilar) {
Node node;
node.position = g1Position;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
}
size_t gnIndex = m_oldAndNewVertices.size();
{
Node node;
node.position = gnPosition;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
}
size_t g2Index = m_oldAndNewVertices.size();
if (!isMismatchSimilar) {
Node node;
node.position = g2Position;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
}
size_t gpIndex = m_oldAndNewVertices.size();
{
Node node;
node.position = gpPosition;
node.radius = gRadius;
node.source = gSource;
m_centerSources.insert(node.source);
m_oldAndNewVertices.push_back(node);
}
if (isMismatchSimilar) {
g1Index = gpIndex;
g2Index = gnIndex;
}
std::map<std::pair<size_t, size_t>, std::vector<size_t>> map;
std::vector<std::vector<size_t>> region1 = {
collectEdgePoints(m, c0EdgeIndex, e0EdgeIndex),
createPointsToMapBetween(e0Index, g1Index, f2, &map),
createPointsToMapBetween(g1Index, e6Index, f0, &map),
collectEdgePoints(q, e6EdgeIndex, c0AlternativeEdgeIndex)
};
m_newRegions.push_back(region1);
if (e5Index != e6Index) {
std::vector<std::vector<size_t>> region2 = {
createPointsToMapBetween(e6Index, g1Index, f0, &map),
createPointsToMapBetween(g1Index, gpIndex, b, &map),
createPointsToMapBetween(gpIndex, e5Index, f0, &map),
collectEdgePoints(q, e5EdgeIndex, e6EdgeIndex)
};
m_newRegions.push_back(region2);
}
std::vector<std::vector<size_t>> region3 = {
createPointsToMapBetween(e5Index, gpIndex, f0, &map),
createPointsToMapBetween(gpIndex, e7Index, f3, &map),
collectEdgePoints(n, e7EdgeIndex, c3AlternativeEdgeIndex),
collectEdgePoints(q, c3EdgeIndex, e5EdgeIndex)
};
m_newRegions.push_back(region3);
std::vector<std::vector<size_t>> region4 = {
collectEdgePoints(m, e2EdgeIndex, c1AlternativeEdgeIndex),
collectEdgePoints(p, c1EdgeIndex, e3EdgeIndex),
createPointsToMapBetween(e3Index, g2Index, f1, &map),
createPointsToMapBetween(g2Index, e2Index, f2, &map)
};
m_newRegions.push_back(region4);
if (e3Index != e4Index) {
std::vector<std::vector<size_t>> region5 = {
createPointsToMapBetween(g2Index, e3Index, f1, &map),
collectEdgePoints(p, e3EdgeIndex, e4EdgeIndex),
createPointsToMapBetween(e4Index, gpIndex, f1, &map),
createPointsToMapBetween(gpIndex, g2Index, a, &map)
};
m_newRegions.push_back(region5);
}
std::vector<std::vector<size_t>> region6 = {
createPointsToMapBetween(gpIndex, e4Index, f1, &map),
collectEdgePoints(p, e4EdgeIndex, c2AlternativeEdgeIndex),
collectEdgePoints(n, c2EdgeIndex, e7EdgeIndex),
createPointsToMapBetween(e7Index, gpIndex, f3, &map)
};
m_newRegions.push_back(region6);
if (g1Index != gnIndex && e0Index != e1Index) {
std::vector<std::vector<size_t>> region7 = {
collectEdgePoints(m, e0EdgeIndex, e1EdgeIndex),
createPointsToMapBetween(e1Index, gnIndex, f2, &map),
createPointsToMapBetween(gnIndex, g1Index, a, &map),
createPointsToMapBetween(g1Index, e0Index, f2, &map)
};
m_newRegions.push_back(region7);
}
if (g2Index != gnIndex && e1Index != e2Index) {
std::vector<std::vector<size_t>> region8 = {
collectEdgePoints(m, e1EdgeIndex, e2EdgeIndex),
createPointsToMapBetween(e2Index, g2Index, f2, &map),
createPointsToMapBetween(g2Index, gnIndex, b, &map),
createPointsToMapBetween(gnIndex, e1Index, f2, &map)
};
m_newRegions.push_back(region8);
}
if (!isMismatchSimilar) {
std::vector<std::vector<size_t>> region9 = {
createPointsToMapBetween(g1Index, gnIndex, a, &map),
createPointsToMapBetween(gnIndex, g2Index, b, &map),
createPointsToMapBetween(g2Index, gpIndex, a, &map),
createPointsToMapBetween(gpIndex, g1Index, b, &map)
};
m_newRegions.push_back(region9);
}
return true;
}
bool RegionFiller::resolveQuadrilateralRegion()
{
int diff02 = std::abs(m_sideSegmentNums[0] - m_sideSegmentNums[2]);
int diff13 = std::abs(m_sideSegmentNums[1] - m_sideSegmentNums[3]);
int m, n, p, q;
if (diff02 >= diff13) {
if (m_sideSegmentNums[0] > m_sideSegmentNums[2]) {
m = 0;
n = 2;
p = 1;
q = 3;
} else {
m = 2;
n = 0;
p = 3;
q = 1;
}
} else {
if (m_sideSegmentNums[1] > m_sideSegmentNums[3]) {
m = 1;
n = 3;
p = 2;
q = 0;
} else {
m = 3;
n = 1;
p = 0;
q = 2;
}
}
bool pqSwapped = false;
if (m_sideSegmentNums[p] < m_sideSegmentNums[q]) {
std::swap(p, q);
pqSwapped = true;
}
if (diff02 != diff13) {
if (isEven(m_sumOfSegments)) {
return resolveQuadrilateralRegionWithIntegerSolution(m, n, p, q, pqSwapped);
} else {
return resolveQuadrilateralRegionWithNonIntegerSolution(m, n, p, q, pqSwapped);
}
} else {
return resolveQuadrilateralRegionMismatchSimilarCase(m, n, p, q, pqSwapped);
}
}
const std::vector<RegionFiller::Node> &RegionFiller::getOldAndNewVertices()
{
return m_oldAndNewVertices;
}
const std::vector<std::vector<size_t>> &RegionFiller::getNewFaces()
{
return m_newFaces;
}
void RegionFiller::prepare()
{
m_sideNum = m_sourcePolylines->size();
m_sideSegmentNums.clear();
for (const auto &it: *m_sourcePolylines) {
auto l = (int)it.size() - 1;
m_sideSegmentNums.push_back(l);
}
m_sumOfSegments = std::accumulate(m_sideSegmentNums.begin(), m_sideSegmentNums.end(), 0);
}
bool RegionFiller::resolveQuadrilateralRegionDirectCase()
{
std::vector<std::vector<size_t>> region = {
collectEdgePoints(0, 0, (*m_sourcePolylines)[0].size() - 1),
collectEdgePoints(1, 0, (*m_sourcePolylines)[1].size() - 1),
collectEdgePoints(2, 0, (*m_sourcePolylines)[2].size() - 1),
collectEdgePoints(3, 0, (*m_sourcePolylines)[3].size() - 1)
};
m_newRegions.push_back(region);
return true;
}
void RegionFiller::fillWithoutPartition()
{
m_oldAndNewVertices = *m_sourceVertices;
m_newFaces.clear();
convertPolylinesToFaces();
}
bool RegionFiller::fill()
{
m_oldAndNewVertices = *m_sourceVertices;
prepare();
if (4 == m_sideNum) {
if (m_sideSegmentNums[0] == m_sideSegmentNums[2] &&
m_sideSegmentNums[1] == m_sideSegmentNums[3]) {
if (!resolveQuadrilateralRegionDirectCase()) {
qDebug() << "resolveQuadrilateralRegionDirectCase failed";
return false;
}
} else {
if (!resolveQuadrilateralRegion()) {
qDebug() << "resolveQuadrilateralRegion failed";
return false;
}
}
} else {
if (isOdd(m_sideNum)) {
if (isEven(m_sumOfSegments)) {
if (!resolveOddSidedEvenSumOfSegments(m_sideNum)) {
qDebug() << "resolveOddSidedEvenSumOfSegments failed, m_sideNum:" << m_sideNum;
return false;
}
} else {
if (!resolveOddSidedOddSumOfSegments(m_sideNum)) {
qDebug() << "resolveOddSidedOddSumOfSegments failed, m_sideNum:" << m_sideNum;
return false;
}
}
} else {
int l1 = 0;
int l2 = 0;
int halfSideNum = m_sideNum / 2;
for (int i = 1; i <= halfSideNum; ++i) {
l1 += m_sideSegmentNums[(i * 2 - 1) % m_sideNum];
l2 += m_sideSegmentNums[(i * 2 - 1 - 1) % m_sideNum];
}
if (isEven(l1) && isEven(l2)) {
if (!resolveEvenSidedEvenSumOfSegmentsAndBothL2L2AreEven()) {
qDebug() << "resolveEvenSidedEvenSumOfSegmentsAndBothL2L2AreEven failed";
return false;
}
} else if (isOdd(l1) && isOdd(l2)) {
if (!resolveEvenSidedEvenSumOfSegmentsAndBothL2L2AreOdd()) {
qDebug() << "resolveEvenSidedEvenSumOfSegmentsAndBothL2L2AreOdd failed";
return false;
}
} else if (isEven(l1) && isOdd(l2)) {
if (!resolveEvenSideEvenSumOfSegmentsAndL1IsEvenL2IsOdd()) {
qDebug() << "resolveEvenSideEvenSumOfSegmentsAndL1IsEvenL2IsOdd failed";
return false;
}
} else {
if (!resolveEvenSideEvenSumOfSegmentsAndL1IsOddL2IsEven()) {
qDebug() << "resolveEvenSideEvenSumOfSegmentsAndL1IsOddL2IsEven failed";
return false;
}
}
}
}
if (!convertRegionsToPatches()) {
qDebug() << "convertRegionsToPatches failed";
return false;
}
return true;
}
void RegionFiller::convertRegionsToFaces()
{
if (m_newRegions.empty()) {
qDebug() << "No region, should not happen";
return;
}
for (const auto &region: m_newRegions) {
std::vector<size_t> face;
for (const auto &line: region) {
for (size_t i = 1; i < line.size(); ++i) {
face.push_back(line[i]);
}
}
if (face.size() < 3) {
qDebug() << "Face length invalid:" << face.size();
continue;
}
m_newFaces.push_back(face);
}
}
void RegionFiller::convertPolylinesToFaces()
{
if (m_sourcePolylines->empty())
return;
std::vector<size_t> face;
for (const auto &line: *m_sourcePolylines) {
for (size_t i = 1; i < line.size(); ++i) {
face.push_back(line[i]);
}
}
if (face.size() < 3) {
return;
}
m_newFaces.push_back(face);
}
bool RegionFiller::createCoonsPatchFrom(const std::vector<size_t> &c0,
const std::vector<size_t> &c1,
const std::vector<size_t> &d0,
const std::vector<size_t> &d1,
bool fillLastTriangle)
{
auto Lc_Position = [&](int s, int t) {
float factor = (float)t / d0.size();
return (1.0 - factor) * m_oldAndNewVertices[c0[s]].position + factor * m_oldAndNewVertices[c1[s]].position;
};
auto Ld_Position = [&](int s, int t) {
float factor = (float)s / c0.size();
return (1.0 - factor) * m_oldAndNewVertices[d0[t]].position + factor * m_oldAndNewVertices[d1[t]].position;
};
auto B_Position = [&](int s, int t) {
float tFactor = (float)t / d0.size();
float sFactor = (float)s / c0.size();
return m_oldAndNewVertices[c0[0]].position * (1.0 - sFactor) * (1.0 - tFactor) +
m_oldAndNewVertices[c0[c0.size() - 1]].position * sFactor * (1.0 - tFactor) +
m_oldAndNewVertices[c1[0]].position * (1.0 - sFactor) * tFactor +
m_oldAndNewVertices[c1[c1.size() - 1]].position * sFactor * tFactor;
};
auto C_Position = [&](int s, int t) {
return Lc_Position(s, t) + Ld_Position(s, t) - B_Position(s, t);
};
auto Lc_Radius = [&](int s, int t) {
float factor = (float)t / d0.size();
return (1.0 - factor) * m_oldAndNewVertices[c0[s]].radius + factor * m_oldAndNewVertices[c1[s]].radius;
};
auto Ld_Radius = [&](int s, int t) {
float factor = (float)s / c0.size();
return (1.0 - factor) * m_oldAndNewVertices[d0[t]].radius + factor * m_oldAndNewVertices[d1[t]].radius;
};
auto B_Radius = [&](int s, int t) {
float tFactor = (float)t / d0.size();
float sFactor = (float)s / c0.size();
return m_oldAndNewVertices[c0[0]].radius * (1.0 - sFactor) * (1.0 - tFactor) +
m_oldAndNewVertices[c0[c0.size() - 1]].radius * sFactor * (1.0 - tFactor) +
m_oldAndNewVertices[c1[0]].radius * (1.0 - sFactor) * tFactor +
m_oldAndNewVertices[c1[c1.size() - 1]].radius * sFactor * tFactor;
};
auto C_Radius = [&](int s, int t) {
return Lc_Radius(s, t) + Ld_Radius(s, t) - B_Radius(s, t);
};
auto getSource = [&](size_t i, size_t j, float factor) {
if (m_centerSources.find(m_oldAndNewVertices[j].source) != m_centerSources.end()) {
return m_oldAndNewVertices[i].source;
} else if (m_centerSources.find(m_oldAndNewVertices[i].source) != m_centerSources.end()) {
return m_oldAndNewVertices[j].source;
} else {
if (factor > 0.5) {
return m_oldAndNewVertices[j].source;
}
return m_oldAndNewVertices[i].source;
}
};
auto Lc_Source = [&](int s, int t) {
float factor = (float)t / d0.size();
return getSource(c0[s], c1[s], factor);
};
auto C_Source = [&](int s, int t) {
return Lc_Source(s, t);
};
std::vector<std::vector<size_t>> grid(c0.size());
for (int s = 1; s < (int)c0.size() - 1; ++s) {
grid[s].resize(d0.size());
for (int t = 1; t < (int)d0.size() - 1; ++t) {
Node node;
node.position = C_Position(s, t);
node.radius = C_Radius(s, t);
node.source = C_Source(s, t);
grid[s][t] = m_oldAndNewVertices.size();
m_oldAndNewVertices.push_back(node);
}
}
grid[0].resize(d0.size());
grid[c0.size() - 1].resize(d0.size());
for (size_t i = 0; i < c0.size(); ++i) {
grid[i][0] = c0[i];
grid[i][d0.size() - 1] = c1[i];
}
for (size_t i = 0; i < d0.size(); ++i) {
grid[0][i] = d0[i];
grid[c0.size() - 1][i] = d1[i];
}
for (int s = 1; s < (int)c0.size(); ++s) {
for (int t = 1; t < (int)d0.size(); ++t) {
std::vector<size_t> face = {
grid[s - 1][t - 1],
grid[s - 1][t],
grid[s][t],
grid[s][t - 1]
};
m_newFaces.push_back(face);
}
}
if (fillLastTriangle) {
std::vector<size_t> face = {
grid[c0.size() - 1][d0.size() - 2],
grid[c0.size() - 2][d0.size() - 2],
grid[c0.size() - 2][d0.size() - 1]
};
m_newFaces.push_back(face);
}
return true;
}
bool RegionFiller::createCoonsPatch(const std::vector<std::vector<size_t>> &region)
{
// https://en.wikipedia.org/wiki/Coons_patch
if (region.size() != 4) {
if (region.size() == 3) {
createCoonsPatchThreeSidedRegion(region);
return true;
}
qDebug() << "Invalid region edges:" << region.size();
return false;
}
const auto &c0 = region[0];
auto c1 = region[2];
std::reverse(c1.begin(), c1.end());
const auto &d1 = region[1];
auto d0 = region[3];
std::reverse(d0.begin(), d0.end());
if (c0.empty() ||
c0.size() != c1.size() ||
d0.empty() ||
d0.size() != d1.size()) {
qDebug() << "Invalid region size:" << c0.size() << c1.size() << d0.size() << d1.size();
return false;
}
return createCoonsPatchFrom(c0, c1, d0, d1, false);
}
bool RegionFiller::createCoonsPatchThreeSidedRegion(const std::vector<std::vector<size_t>> &region)
{
if (region.size() != 3) {
qDebug() << "Not three sided region";
return false;
}
int longest = 0;
size_t longestLength = region[longest].size();
for (int i = 1; i < 3; ++i) {
if (region[i].size() > longestLength) {
longest = i;
longestLength = region[i].size();
}
}
int c = longest;
int a = (c + 1) % 3;
int b = (a + 1) % 3;
auto c0 = region[a];
std::reverse(c0.begin(), c0.end());
auto c1 = std::vector<size_t>(region[c].begin(), region[c].begin() + c0.size());
auto d0 = region[b];
auto d1 = std::vector<size_t>(region[c].begin() + region[c].size() - d0.size(), region[c].begin() + region[c].size());
std::reverse(d1.begin(), d1.end());
bool fillTriangle = region[c].size() != region[a].size() &&
region[a].size() != region[b].size();
if (c0.size() < d0.size()) {
if (!createCoonsPatchFrom(d0, d1, c0, c1, fillTriangle))
return false;
} else {
if (!createCoonsPatchFrom(c0, c1, d0, d1, fillTriangle))
return false;
}
return true;
}
bool RegionFiller::convertRegionsToPatches()
{
if (m_newRegions.empty()) {
qDebug() << "No region, should not happen";
return false;
}
for (const auto &region: m_newRegions) {
if (!createCoonsPatch(region))
return false;
}
return true;
}