1 /*
2  *  Copyright (c) 2014 Dmitry Kazakov <dimula73@gmail.com>
3  *
4  *  This program is free software; you can redistribute it and/or modify
5  *  it under the terms of the GNU General Public License as published by
6  *  the Free Software Foundation; either version 2 of the License, or
7  *  (at your option) any later version.
8  *
9  *  This program is distributed in the hope that it will be useful,
10  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
11  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
12  *  GNU General Public License for more details.
13  *
14  *  You should have received a copy of the GNU General Public License
15  *  along with this program; if not, write to the Free Software
16  *  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
17  */
18 
19 #ifndef __KIS_GRID_INTERPOLATION_TOOLS_H
20 #define __KIS_GRID_INTERPOLATION_TOOLS_H
21 
22 #include <limits>
23 #include <algorithm>
24 
25 #include <QImage>
26 
27 #include "kis_algebra_2d.h"
28 #include "kis_four_point_interpolator_forward.h"
29 #include "kis_four_point_interpolator_backward.h"
30 #include "kis_iterator_ng.h"
31 #include "kis_random_sub_accessor.h"
32 
33 //#define DEBUG_PAINTING_POLYGONS
34 
35 #ifdef DEBUG_PAINTING_POLYGONS
36 #include <QPainter>
37 #endif /* DEBUG_PAINTING_POLYGONS */
38 
39 namespace GridIterationTools {
40 
calcGridDimension(int start,int end,const int pixelPrecision)41 inline int calcGridDimension(int start, int end, const int pixelPrecision)
42 {
43     const int alignmentMask = ~(pixelPrecision - 1);
44 
45     int alignedStart = (start + pixelPrecision - 1) & alignmentMask;
46     int alignedEnd = end & alignmentMask;
47 
48     int size = 0;
49 
50     if (alignedEnd > alignedStart) {
51         size = (alignedEnd - alignedStart) / pixelPrecision + 1;
52         size += alignedStart != start;
53         size += alignedEnd != end;
54     } else {
55         size = 2 + (end - start >= pixelPrecision);
56     }
57 
58     return size;
59 }
60 
calcGridSize(const QRect & srcBounds,const int pixelPrecision)61 inline QSize calcGridSize(const QRect &srcBounds, const int pixelPrecision) {
62     return QSize(calcGridDimension(srcBounds.x(), srcBounds.right(), pixelPrecision),
63                  calcGridDimension(srcBounds.y(), srcBounds.bottom(), pixelPrecision));
64 }
65 
66 template <class ProcessPolygon, class ForwardTransform>
67 struct CellOp
68 {
CellOpCellOp69     CellOp(ProcessPolygon &_polygonOp, ForwardTransform &_transformOp)
70         : polygonOp(_polygonOp),
71           transformOp(_transformOp)
72     {
73     }
74 
processPointCellOp75     inline void processPoint(int col, int row,
76                              int prevCol, int prevRow,
77                              int colIndex, int rowIndex) {
78 
79         QPointF dstPosF = transformOp(QPointF(col, row));
80         currLinePoints << dstPosF;
81 
82         if (rowIndex >= 1 && colIndex >= 1) {
83             QPolygonF srcPolygon;
84 
85             srcPolygon << QPointF(prevCol, prevRow);
86             srcPolygon << QPointF(col, prevRow);
87             srcPolygon << QPointF(col, row);
88             srcPolygon << QPointF(prevCol, row);
89 
90             QPolygonF dstPolygon;
91 
92             dstPolygon << prevLinePoints.at(colIndex - 1);
93             dstPolygon << prevLinePoints.at(colIndex);
94             dstPolygon << currLinePoints.at(colIndex);
95             dstPolygon << currLinePoints.at(colIndex - 1);
96 
97             polygonOp(srcPolygon, dstPolygon);
98         }
99 
100     }
101 
nextLineCellOp102     inline void nextLine() {
103         std::swap(prevLinePoints, currLinePoints);
104 
105         // we are erasing elements for not free'ing the occupied
106         // memory, which is more efficient since we are going to fill
107         // the vector again
108         currLinePoints.erase(currLinePoints.begin(), currLinePoints.end());
109     }
110 
111     QVector<QPointF> prevLinePoints;
112     QVector<QPointF> currLinePoints;
113     ProcessPolygon &polygonOp;
114     ForwardTransform &transformOp;
115 };
116 
117 template <class ProcessCell>
processGrid(ProcessCell & cellOp,const QRect & srcBounds,const int pixelPrecision)118 void processGrid(ProcessCell &cellOp,
119                  const QRect &srcBounds,
120                  const int pixelPrecision)
121 {
122     if (srcBounds.isEmpty()) return;
123 
124     const int alignmentMask = ~(pixelPrecision - 1);
125 
126     int prevRow = std::numeric_limits<int>::max();
127     int prevCol = std::numeric_limits<int>::max();
128 
129     int rowIndex = 0;
130     int colIndex = 0;
131 
132     for (int row = srcBounds.top(); row <= srcBounds.bottom();) {
133         for (int col = srcBounds.left(); col <= srcBounds.right();) {
134 
135             cellOp.processPoint(col, row,
136                                 prevCol, prevRow,
137                                 colIndex, rowIndex);
138 
139             prevCol = col;
140             col += pixelPrecision;
141             colIndex++;
142 
143             if (col > srcBounds.right() &&
144                 col <= srcBounds.right() + pixelPrecision - 1) {
145 
146                 col = srcBounds.right();
147             } else {
148                 col &= alignmentMask;
149             }
150         }
151 
152         cellOp.nextLine();
153         colIndex = 0;
154 
155         prevRow = row;
156         row += pixelPrecision;
157         rowIndex++;
158 
159         if (row > srcBounds.bottom() &&
160             row <= srcBounds.bottom() + pixelPrecision - 1) {
161 
162             row = srcBounds.bottom();
163         } else {
164             row &= alignmentMask;
165         }
166     }
167 }
168 
169 template <class ProcessPolygon, class ForwardTransform>
processGrid(ProcessPolygon & polygonOp,ForwardTransform & transformOp,const QRect & srcBounds,const int pixelPrecision)170 void processGrid(ProcessPolygon &polygonOp, ForwardTransform &transformOp,
171                  const QRect &srcBounds, const int pixelPrecision)
172 {
173     CellOp<ProcessPolygon, ForwardTransform> cellOp(polygonOp, transformOp);
174     processGrid(cellOp, srcBounds, pixelPrecision);
175 }
176 
177 struct PaintDevicePolygonOp
178 {
PaintDevicePolygonOpPaintDevicePolygonOp179     PaintDevicePolygonOp(KisPaintDeviceSP srcDev, KisPaintDeviceSP dstDev)
180         : m_srcDev(srcDev), m_dstDev(dstDev) {}
181 
operatorPaintDevicePolygonOp182     void operator() (const QPolygonF &srcPolygon, const QPolygonF &dstPolygon) {
183         this->operator() (srcPolygon, dstPolygon, dstPolygon);
184     }
185 
operatorPaintDevicePolygonOp186     void operator() (const QPolygonF &srcPolygon, const QPolygonF &dstPolygon, const QPolygonF &clipDstPolygon) {
187         QRect boundRect = clipDstPolygon.boundingRect().toAlignedRect();
188         if (boundRect.isEmpty()) return;
189 
190         KisSequentialIterator dstIt(m_dstDev, boundRect);
191         KisRandomSubAccessorSP srcAcc = m_srcDev->createRandomSubAccessor();
192 
193         KisFourPointInterpolatorBackward interp(srcPolygon, dstPolygon);
194 
195         int y = boundRect.top();
196         interp.setY(y);
197 
198         while (dstIt.nextPixel()) {
199             int newY = dstIt.y();
200 
201             if (y != newY) {
202                 y = newY;
203                 interp.setY(y);
204             }
205 
206             QPointF srcPoint(dstIt.x(), y);
207 
208             if (clipDstPolygon.containsPoint(srcPoint, Qt::OddEvenFill)) {
209 
210                 interp.setX(srcPoint.x());
211                 QPointF dstPoint = interp.getValue();
212 
213                 // brain-blowing part:
214                 //
215                 // since the interpolator does the inverted
216                 // transformation we read data from "dstPoint"
217                 // (which is non-transformed) and write it into
218                 // "srcPoint" (which is transformed position)
219 
220                 srcAcc->moveTo(dstPoint);
221                 srcAcc->sampledOldRawData(dstIt.rawData());
222             }
223 
224         }
225 
226     }
227 
228     KisPaintDeviceSP m_srcDev;
229     KisPaintDeviceSP m_dstDev;
230 };
231 
232 struct QImagePolygonOp
233 {
QImagePolygonOpQImagePolygonOp234     QImagePolygonOp(const QImage &srcImage, QImage &dstImage,
235                     const QPointF &srcImageOffset,
236                     const QPointF &dstImageOffset)
237         : m_srcImage(srcImage), m_dstImage(dstImage),
238           m_srcImageOffset(srcImageOffset),
239           m_dstImageOffset(dstImageOffset),
240           m_srcImageRect(m_srcImage.rect()),
241           m_dstImageRect(m_dstImage.rect())
242     {
243     }
244 
operatorQImagePolygonOp245     void operator() (const QPolygonF &srcPolygon, const QPolygonF &dstPolygon) {
246         this->operator() (srcPolygon, dstPolygon, dstPolygon);
247     }
248 
operatorQImagePolygonOp249     void operator() (const QPolygonF &srcPolygon, const QPolygonF &dstPolygon, const QPolygonF &clipDstPolygon) {
250         QRect boundRect = clipDstPolygon.boundingRect().toAlignedRect();
251         KisFourPointInterpolatorBackward interp(srcPolygon, dstPolygon);
252 
253         for (int y = boundRect.top(); y <= boundRect.bottom(); y++) {
254             interp.setY(y);
255             for (int x = boundRect.left(); x <= boundRect.right(); x++) {
256 
257                 QPointF srcPoint(x, y);
258                 if (clipDstPolygon.containsPoint(srcPoint, Qt::OddEvenFill)) {
259 
260                     interp.setX(srcPoint.x());
261                     QPointF dstPoint = interp.getValue();
262 
263                     // about srcPoint/dstPoint hell please see a
264                     // comment in PaintDevicePolygonOp::operator() ()
265 
266                     srcPoint -= m_dstImageOffset;
267                     dstPoint -= m_srcImageOffset;
268 
269                     QPoint srcPointI = srcPoint.toPoint();
270                     QPoint dstPointI = dstPoint.toPoint();
271 
272                     if (!m_dstImageRect.contains(srcPointI)) continue;
273                     if (!m_srcImageRect.contains(dstPointI)) continue;
274 
275                     m_dstImage.setPixel(srcPointI, m_srcImage.pixel(dstPointI));
276                 }
277             }
278         }
279 
280 #ifdef DEBUG_PAINTING_POLYGONS
281         QPainter gc(&m_dstImage);
282         gc.setPen(Qt::red);
283         gc.setOpacity(0.5);
284 
285         gc.setBrush(Qt::green);
286         gc.drawPolygon(clipDstPolygon.translated(-m_dstImageOffset));
287 
288         gc.setBrush(Qt::blue);
289         //gc.drawPolygon(dstPolygon.translated(-m_dstImageOffset));
290 
291 #endif /* DEBUG_PAINTING_POLYGONS */
292 
293     }
294 
295     const QImage &m_srcImage;
296     QImage &m_dstImage;
297     QPointF m_srcImageOffset;
298     QPointF m_dstImageOffset;
299 
300     QRect m_srcImageRect;
301     QRect m_dstImageRect;
302 };
303 
304 /*************************************************************/
305 /*      Iteration through precalculated grid                 */
306 /*************************************************************/
307 
308 /**
309  *    A-----B         The polygons will be in the following order:
310  *    |     |
311  *    |     |         polygon << A << B << D << C;
312  *    C-----D
313  */
calculateCellIndexes(int col,int row,const QSize & gridSize)314 inline QVector<int> calculateCellIndexes(int col, int row, const QSize &gridSize)
315 {
316     const int tl = col + row * gridSize.width();
317     const int tr = tl + 1;
318     const int bl = tl + gridSize.width();
319     const int br = bl + 1;
320 
321     QVector<int> cellIndexes;
322     cellIndexes << tl;
323     cellIndexes << tr;
324     cellIndexes << br;
325     cellIndexes << bl;
326 
327     return cellIndexes;
328 }
329 
pointToIndex(const QPoint & cellPt,const QSize & gridSize)330 inline int pointToIndex(const QPoint &cellPt, const QSize &gridSize)
331 {
332     return cellPt.x() +
333         cellPt.y() * gridSize.width();
334 }
335 
336 namespace Private {
pointPolygonIndexToColRow(QPoint baseColRow,int index)337     inline QPoint pointPolygonIndexToColRow(QPoint baseColRow, int index)
338     {
339         static QVector<QPoint> pointOffsets;
340         if (pointOffsets.isEmpty()) {
341             pointOffsets << QPoint(0,0);
342             pointOffsets << QPoint(1,0);
343             pointOffsets << QPoint(1,1);
344             pointOffsets << QPoint(0,1);
345         }
346 
347         return baseColRow + pointOffsets[index];
348     }
349 
350     struct PointExtension {
351         int near;
352         int far;
353     };
354 }
355 
356 template <class IndexesOp>
getOrthogonalPointApproximation(const QPoint & cellPt,const QVector<QPointF> & originalPoints,const QVector<QPointF> & transformedPoints,IndexesOp indexesOp,QPointF * srcPoint,QPointF * dstPoint)357 bool getOrthogonalPointApproximation(const QPoint &cellPt,
358                                 const QVector<QPointF> &originalPoints,
359                                 const QVector<QPointF> &transformedPoints,
360                                 IndexesOp indexesOp,
361                                 QPointF *srcPoint,
362                                 QPointF *dstPoint)
363 {
364     QVector<Private::PointExtension> extensionPoints;
365     Private::PointExtension ext;
366 
367     // left
368     if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(-1, 0))) >= 0 &&
369         (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(-2, 0))) >= 0) {
370 
371         extensionPoints << ext;
372     }
373     // top
374     if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(0, -1))) >= 0 &&
375         (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(0, -2))) >= 0) {
376 
377         extensionPoints << ext;
378     }
379     // right
380     if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(1, 0))) >= 0 &&
381         (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(2, 0))) >= 0) {
382 
383         extensionPoints << ext;
384     }
385     // bottom
386     if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(0, 1))) >= 0 &&
387         (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(0, 2))) >= 0) {
388 
389         extensionPoints << ext;
390     }
391 
392     if (extensionPoints.isEmpty()) {
393         // top-left
394         if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(-1, -1))) >= 0 &&
395             (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(-2, -2))) >= 0) {
396 
397             extensionPoints << ext;
398         }
399         // top-right
400         if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(1, -1))) >= 0 &&
401             (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(2, -2))) >= 0) {
402 
403             extensionPoints << ext;
404         }
405         // bottom-right
406         if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(1, 1))) >= 0 &&
407             (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(2, 2))) >= 0) {
408 
409             extensionPoints << ext;
410         }
411         // bottom-left
412         if ((ext.near = indexesOp.tryGetValidIndex(cellPt + QPoint(-1, 1))) >= 0 &&
413             (ext.far = indexesOp.tryGetValidIndex(cellPt + QPoint(-2, 2))) >= 0) {
414 
415             extensionPoints << ext;
416         }
417     }
418 
419     if (extensionPoints.isEmpty()) {
420         return false;
421     }
422 
423     int numResultPoints = 0;
424     *srcPoint = indexesOp.getSrcPointForce(cellPt);
425     *dstPoint = QPointF();
426 
427     Q_FOREACH (const Private::PointExtension &ext, extensionPoints) {
428         QPointF near = transformedPoints[ext.near];
429         QPointF far = transformedPoints[ext.far];
430 
431         QPointF nearSrc = originalPoints[ext.near];
432         QPointF farSrc = originalPoints[ext.far];
433 
434         QPointF base1 = nearSrc - farSrc;
435         QPointF base2 = near - far;
436 
437         QPointF pt = near +
438             KisAlgebra2D::transformAsBase(*srcPoint - nearSrc, base1, base2);
439 
440         *dstPoint += pt;
441         numResultPoints++;
442     }
443 
444     *dstPoint /= numResultPoints;
445 
446     return true;
447 }
448 
449 template <class PolygonOp, class IndexesOp>
450 struct IncompletePolygonPolicy {
451 
tryProcessPolygonIncompletePolygonPolicy452     static inline bool tryProcessPolygon(int col, int row,
453                                          int numExistingPoints,
454                                          PolygonOp &polygonOp,
455                                          IndexesOp &indexesOp,
456                                          const QVector<int> &polygonPoints,
457                                          const QVector<QPointF> &originalPoints,
458                                          const QVector<QPointF> &transformedPoints)
459     {
460         if (numExistingPoints >= 4) return false;
461         if (numExistingPoints == 0) return true;
462 
463         QPolygonF srcPolygon;
464         QPolygonF dstPolygon;
465 
466         for (int i = 0; i < 4; i++) {
467             const int index = polygonPoints[i];
468 
469             if (index >= 0) {
470                 srcPolygon << originalPoints[index];
471                 dstPolygon << transformedPoints[index];
472             } else {
473                 QPoint cellPt = Private::pointPolygonIndexToColRow(QPoint(col, row), i);
474                 QPointF srcPoint;
475                 QPointF dstPoint;
476                 bool result =
477                     getOrthogonalPointApproximation(cellPt,
478                                                     originalPoints,
479                                                     transformedPoints,
480                                                     indexesOp,
481                                                     &srcPoint,
482                                                     &dstPoint);
483 
484                 if (!result) {
485                     //dbgKrita << "*NOT* found any valid point" << allSrcPoints[pointToIndex(cellPt)] << "->" << ppVar(pt);
486                     break;
487                 } else {
488                     srcPolygon << srcPoint;
489                     dstPolygon << dstPoint;
490                 }
491             }
492         }
493 
494         if (dstPolygon.size() == 4) {
495             QPolygonF srcClipPolygon(srcPolygon.intersected(indexesOp.srcCropPolygon()));
496 
497             KisFourPointInterpolatorForward forwardTransform(srcPolygon, dstPolygon);
498             for (int i = 0; i < srcClipPolygon.size(); i++) {
499                 const QPointF newPt = forwardTransform.map(srcClipPolygon[i]);
500                 srcClipPolygon[i] = newPt;
501             }
502 
503             polygonOp(srcPolygon, dstPolygon, srcClipPolygon);
504         }
505 
506         return true;
507     }
508 };
509 
510 template <class PolygonOp, class IndexesOp>
511 struct AlwaysCompletePolygonPolicy {
512 
tryProcessPolygonAlwaysCompletePolygonPolicy513     static inline bool tryProcessPolygon(int col, int row,
514                                          int numExistingPoints,
515                                          PolygonOp &polygonOp,
516                                          IndexesOp &indexesOp,
517                                          const QVector<int> &polygonPoints,
518                                          const QVector<QPointF> &originalPoints,
519                                          const QVector<QPointF> &transformedPoints)
520     {
521         Q_UNUSED(col);
522         Q_UNUSED(row);
523         Q_UNUSED(polygonOp);
524         Q_UNUSED(indexesOp);
525         Q_UNUSED(polygonPoints);
526         Q_UNUSED(originalPoints);
527         Q_UNUSED(transformedPoints);
528 
529         KIS_ASSERT_RECOVER_NOOP(numExistingPoints == 4);
530         return false;
531     }
532 };
533 
534 struct RegularGridIndexesOp {
535 
RegularGridIndexesOpRegularGridIndexesOp536     RegularGridIndexesOp(const QSize &gridSize)
537         : m_gridSize(gridSize)
538     {
539     }
540 
calculateMappedIndexesRegularGridIndexesOp541     inline QVector<int> calculateMappedIndexes(int col, int row,
542                                                int *numExistingPoints) const {
543 
544         *numExistingPoints = 4;
545         QVector<int> cellIndexes =
546             GridIterationTools::calculateCellIndexes(col, row, m_gridSize);
547 
548         return cellIndexes;
549     }
550 
tryGetValidIndexRegularGridIndexesOp551     inline int tryGetValidIndex(const QPoint &cellPt) const {
552         Q_UNUSED(cellPt);
553 
554         KIS_ASSERT_RECOVER_NOOP(0 && "Not applicable");
555         return -1;
556     }
557 
getSrcPointForceRegularGridIndexesOp558     inline QPointF getSrcPointForce(const QPoint &cellPt) const {
559         Q_UNUSED(cellPt);
560 
561         KIS_ASSERT_RECOVER_NOOP(0 && "Not applicable");
562         return QPointF();
563     }
564 
srcCropPolygonRegularGridIndexesOp565     inline const QPolygonF srcCropPolygon() const {
566         KIS_ASSERT_RECOVER_NOOP(0 && "Not applicable");
567         return QPolygonF();
568     }
569 
570     QSize m_gridSize;
571 };
572 
573 /**
574  * There is a weird problem in fetching correct bounds of the polygon.
575  * If the rightmost (bottommost) point of the polygon is integral, then
576  * QRectF() will end exactly on it, but when converting into QRect the last
577  * point will not be taken into account. It happens due to the difference
578  * between center-point/topleft-point point representation. In many cases
579  * the latter is expected, but we don't work with it in Qt/Krita.
580  */
adjustAlignedPolygon(QPolygonF & polygon)581 inline void adjustAlignedPolygon(QPolygonF &polygon)
582 {
583     static const qreal eps = 1e-5;
584     static const  QPointF p1(eps, 0.0);
585     static const  QPointF p2(eps, eps);
586     static const  QPointF p3(0.0, eps);
587 
588     polygon[1] += p1;
589     polygon[2] += p2;
590     polygon[3] += p3;
591 }
592 
593 template <template <class PolygonOp, class IndexesOp> class IncompletePolygonPolicy,
594           class PolygonOp,
595           class IndexesOp>
iterateThroughGrid(PolygonOp & polygonOp,IndexesOp & indexesOp,const QSize & gridSize,const QVector<QPointF> & originalPoints,const QVector<QPointF> & transformedPoints)596 void iterateThroughGrid(PolygonOp &polygonOp,
597                         IndexesOp &indexesOp,
598                         const QSize &gridSize,
599                         const QVector<QPointF> &originalPoints,
600                         const QVector<QPointF> &transformedPoints)
601 {
602     QVector<int> polygonPoints(4);
603 
604     for (int row = 0; row < gridSize.height() - 1; row++) {
605         for (int col = 0; col < gridSize.width() - 1; col++) {
606             int numExistingPoints = 0;
607 
608             polygonPoints = indexesOp.calculateMappedIndexes(col, row, &numExistingPoints);
609 
610             if (!IncompletePolygonPolicy<PolygonOp, IndexesOp>::
611                  tryProcessPolygon(col, row,
612                                    numExistingPoints,
613                                    polygonOp,
614                                    indexesOp,
615                                    polygonPoints,
616                                    originalPoints,
617                                    transformedPoints)) {
618 
619                 QPolygonF srcPolygon;
620                 QPolygonF dstPolygon;
621 
622                 for (int i = 0; i < 4; i++) {
623                     const int index = polygonPoints[i];
624                     srcPolygon << originalPoints[index];
625                     dstPolygon << transformedPoints[index];
626                 }
627 
628                 adjustAlignedPolygon(srcPolygon);
629                 adjustAlignedPolygon(dstPolygon);
630 
631                 polygonOp(srcPolygon, dstPolygon);
632             }
633         }
634     }
635 }
636 
637 }
638 
639 #endif /* __KIS_GRID_INTERPOLATION_TOOLS_H */
640