// Copyright 2020 Rich Altmaier  richalt2@yahoo.com
//
//
//   This program is free software: you can redistribute it and/or modify
//   it under the terms of the GNU General Public License as published by
//   the Free Software Foundation, either version 3 of the License, or
//   (at your option) any later version.
//
//   This program is distributed in the hope that it will be useful,
//   but WITHOUT ANY WARRANTY; without even the implied warranty of
//   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//   GNU General Public License for more details.
//
//   You should have received a copy of the GNU General Public License
//   along with this program.  If not, see <http://www.gnu.org/licenses/>.
//


// show an image bounded by an arbitrary polygon

//nSplines to sweep out image
//nSplines openscad library from Parkinbot on thingiverse,
//  https://www.thingiverse.com/thing:1208001/files 
//Install nSplines in a peer directory and point to these two files:
use <../../curves/nSplines-2019/files/splines.scad>
use <../../curves/nSplines-2019/files/Naca_sweep.scad>
//NOTE: Parkinbot says on Feb 20, 2019 that convex polygons are fine with:
//      sweep(dat, planar_caps = true);





   
//interpolate x, given y between point P and Q
function interpolateX2(y, px, py, qx, qy) =
     px +  (qx - px) * (y - py) / (qy - py);


    
function quicksortI(arr, elementI) =
   //given arr as a list of lists, we want to sort each element on a selected key
   //sort arr using [elementI] as key
   //if arr is [], return []
   //if arr is single element, return that element
  !(len(arr)>1) ? arr :
      let(  pivot   = arr[floor(len(arr)/2)][elementI],
            lesser  = [ for (y = arr) if (y[elementI]  < pivot) y ],
            equal   = [ for (y = arr) if (y[elementI] == pivot) y ],
            greater = [ for (y = arr) if (y[elementI]  > pivot) y ]
      )
      concat( quicksortI(lesser, elementI), equal, quicksortI(greater, elementI) );

//determine inside/outside of an irregular polygon
// polygon defined by a trace, a list of XY points, where we assume first point
//  is joined to last trace entry to form closure. 
// NOTE: this is not an openscad polygon(), which can have holes.
// NOTE: given typical style of trace is this design, points likely must be added
//   to traces to extend the polygon to base of likely images. 
//   This must be done before calling these functions, as application may need. 
// Algorithm: usage of this data is by an image line traversing constant-Y, so
//   1.  at y=Y, search list for pairs of points spanning Y,
//   2.  interpolate the X value at exact y crossing
//   3.  sort the resulting list of [interpolated X, Y] in order of increasing x

function LT_wrapNextPt(trace, i) =
     i >= len(trace)-1 ? 0 : i+1 ; 
  
function LT_findPolyEdges(trace, atY, traceI=0, returnList=[]) =
    //Note: make sure atY is in the same coordinates as the trace list passed! 
    traceI >= len(trace) ? returnList : //have searched up to len-1 : 0.
    let (  thisPt = trace[traceI],
           nextPt = trace[LT_wrapNextPt(trace,traceI)],
           thisPairSpansatY = (atY >= thisPt[1] && atY < nextPt[1] )
                           || (atY < thisPt[1] && atY >= nextPt[1] ),
                     //note: tangents are not crossings! 
           interpX = interpolateX2(atY, thisPt[0], thisPt[1], nextPt[0], nextPt[1]),
              //valid only if thisPairSpansatY == true
           possibleSpan = [interpX, atY],
           newReturnList = thisPairSpansatY ? concat (returnList, [possibleSpan]) :
                           returnList
          //, xx =  echo("LT_findPolyEdges found pair: atY, spanBool, thisPt, nextPt, possibleSpan", atY, thisPairSpansatY, thisPt, nextPt, possibleSpan) 
     )
     LT_findPolyEdges(trace, atY, traceI+1,  newReturnList)
     ;

function LT_findPolyEdgesSortedXasc(trace, atY) =
    let ( Crossings = LT_findPolyEdges(trace, atY),
        CrossingsSorted = quicksortI(Crossings, 0) //sort using index 0
    )
    CrossingsSorted ; 
    
    
        
function LT_polyEdgeCrossingCount(X, CrossingsSorted, CrossingsI=0 ) =
     //return count of CrossingsSorted entries for which  X < 
     len(CrossingsSorted) < 1 ?  0 : //no crossing possible
     CrossingsI >= len(CrossingsSorted) ? CrossingsI : // X crossed all entries
     X < CrossingsSorted[CrossingsI][0] ? CrossingsI : // X did not pass I
     LT_polyEdgeCrossingCount(X, CrossingsSorted, CrossingsI+1 ) ;

function LT_PtInsidePolywCrossings(X, CrossingsForYSorted) =
     let ( crossCount = LT_polyEdgeCrossingCount(X, CrossingsForYSorted),
           InBool = (crossCount % 2) ==1
     )
     InBool; 
     
function LT_PtInsidePoly(trace, XYpt) =
    let ( cs = LT_findPolyEdgesSortedXasc(trace, XYpt[1]),
       crossCount = LT_polyEdgeCrossingCount(XYpt[0], cs),
       InBool = (crossCount % 2) == 1
     )
     InBool ;
function LT_minmaxXatEdgewCrossings(CrossingsForYSorted) =
    //return [min, max] values of X for first crossing and last crossing.
    //   if no crossings, return [-1, -1]
    len(CrossingsForYSorted) < 2 ?  [-1, -1] :
    [ CrossingsForYSorted[0][0],  CrossingsForYSorted[len(CrossingsForYSorted)-1][0] ];


     


//show image inside of a polygon (defined by a trace).
module LT_ImageInsidePoly(ImageArrayInput,  trace_osu, minX_osu, maxX_osu, smallYextMaxval_osu, maxYextMinval_osu, showInPlusYbool, imgMag, backThick_osu, flipImageTopToBottomBool=false, flipImageLeftToRightBool=false) {
      //ImageArray is mapped to the rectangle around the trace_osu:
      //    left-right image ->  minX:maxX
      //    top-bot image -> maxYext (min val) : smallYext (max value)
      //    This rectangle permits changing the aspect ratio of the image.
      //    HOWEVER the trace_osu points must be within this rectangle. 
      //    Typically this rectangle is matched to the file which contained the
      //    trace_osu data.

      //NOTE: commonly this trace has been extended to provide a base beyond BG, for
      //   which the corresponding image is aligned. 
      //Return:
      //  generate an openscad polyhedron in the x-z positive plane,
      //  x-z location will be as defined by trace_osu values. 
      //  base of imageArray is placed at largest z, z=smallYextMaxval_osu, 
      //  top of imageArray is placed at smallest z, z=maxYextMinval_osu, 
      //  and bright pixels become variation in the 
      //      showInPlusYbool ? +y direction : -y direction,
      //  with back of polyhedron at showInPlusYBool? y=-backThick : y=backThick.
      //Polyhedron with image can be rotated around X to align with trace extrusion model. 
      //NOTE: showInPlusYbool placing image on opposite sides means the end view
      //      photos must be left-right reversals of each other. 
      //Make the polyhedron points match the number of pixels.
      // Note the polyhyedron requires a constant number of x traversing points, so
      // where the poly trace cuts off, we compress those edge points into a
      // tiny border.
      // Note: trace may have multiple lobs, but this polyhedron must be single lob,
      //   so finding spaces between lobs is treated by dropping that portion of
      //   the front of the polyhedron nearly to the back, and then differenceing it away. 
    
    ImageArrayTB = flipImageTopToBottomBool ? LT_flipImageTopToBottom(ImageArrayInput) : ImageArrayInput; 
    ImageArray = flipImageLeftToRightBool ? LT_flipImageLeftToRight(ImageArrayTB) : ImageArrayTB; 
    NumXpix = len(ImageArray[0]);  
    NumYpix = len(ImageArray); 
    xIncr_osu = (maxX_osu - minX_osu)/ NumXpix;
    yIncr_osu = (smallYextMaxval_osu - maxYextMinval_osu) / NumYpix;

    echo("LT_ImageInsidePoly: pix x y", NumXpix, NumYpix, "X min, max osu, Xincr", minX_osu, maxX_osu, xIncr_osu, "Y min, max osu, Yincr", maxYextMinval_osu,smallYextMaxval_osu, yIncr_osu,  "showInPlusYbool", showInPlusYbool, "trace len", len(trace_osu) , "flip TB, LR", flipImageTopToBottomBool, flipImageLeftToRightBool);
    assert(len(trace_osu) > 3);
    assert(len(ImageArray) > 2 && len(ImageArray[0]) > 2, "expecting 2 dimensional array of photo image data");
    assert(backThick_osu >= 0);
    epsilonBack = 1.0; //0.05; 
    backThickExtended_osu = backThick_osu + 2 * epsilonBack; 

    
   function PixRowCWposY(imgOneRow, atZ, PolyCrossings, Xincr_osu, Xmax_osu, minMaxX, showInPlusYbool, magnify, backThickE_osu, epsilonB) =
      //image is on +y
      //since polyhedron requires clockwise around z, +y sweeps 0 to large X,
     assert(len(imgOneRow) > 2) 
     concat (
        [ for (i=[0:len(imgOneRow)-1]) 
           let ( NxPixinRow = len(imgOneRow),
                 xAtThisPix = i * Xincr_osu,
                 xInsidePolyBool = LT_PtInsidePolywCrossings(xAtThisPix, PolyCrossings),
                 
                 xGen = xAtThisPix, 
        
                 //note: bounding polygon may have multiple lobs, but this generated
                 // polyhedron can have only one.  So for points between lobs, the
                 // image is black AND descends to within epsilonB of backing, so it
                 // can be difference() removed.
        
                 //In both outside-poly, and in lobs,
                 // push poly surface down to back-epsilon.
                 pixRawVal = xInsidePolyBool ? imgOneRow[i] : 0,
                     //portions outside of poly are black
                 OutsidePolyDepress = xInsidePolyBool ? +0.01 : -backThickE_osu + epsilonB
                 //, xx = atZ > 20.0 && atZ < 22 ? echo("PixRowCWposY: xAtThisPix, xInsidePolyBool", xAtThisPix, xInsidePolyBool, "PolyCrossings", PolyCrossings):0 
               )

           [ xGen, (pixRawVal /255.0 * magnify) +OutsidePolyDepress, atZ ]
    
        ],
        [ [(len(imgOneRow)-1)* Xincr_osu, -backThickE_osu, atZ],  [0, -backThickE_osu, atZ]   ] 
     );
        
   function PixRowCWnegY(imgOneRow, atZ, PolyCrossings, Xincr_osu, Xmax_osu, minMaxX,  showInPlusYbool, magnify, backThickE_osu, epsilonB) =
      //image is on -y
      //since polyhedron requires clockwise around z, -y sweeps large X to 0,
     assert(len(imgOneRow) > 2)
     concat (
        [ for (i=[len(imgOneRow)-1 : -1: 0]) 
           let ( NxPixinRow = len(imgOneRow),
                 xAtThisPix = i * Xincr_osu,
                 xInsidePolyBool = LT_PtInsidePolywCrossings(xAtThisPix, PolyCrossings),
                 
                 xGen = xAtThisPix,
        
                 //note: bounding polygon may have multiple lobs, but this generated
                 // polyhedron can have only one.  So for points between lobs, the
                 // image is black AND descends to within epsilonB of backing, so it
                 // can be difference() removed.
                 //In both outside-poly, and in lobs,
                 // push poly surface down to back-epsilon.
                 pixRawVal = xInsidePolyBool ? imgOneRow[len(imgOneRow)-1-i] : 0,
                 OutsidePolyDepress = xInsidePolyBool ? -0.01 : backThickE_osu-epsilonB
                   //when inside poly, shift a small amount away from 0 so black
                   // pixels don't have 0 thickness!
                   //when in a poly lob, shift poly front face to epsilon away from back.
                 //, xx = atZ ==150 ? echo("PixRowCWnegY: xAtThisPix, xInsidePolyBool", xAtThisPix, xInsidePolyBool, "PolyCrossings", PolyCrossings):0 
               )

           [ xGen, (-pixRawVal /255.0 * magnify) +OutsidePolyDepress, atZ ] 
    
        ],
        [ [ 0, backThickE_osu, atZ], [ (len(imgOneRow) -1) * Xincr_osu, backThickE_osu, atZ]   ]

     );


   function RIgrey_gendat(imgarray, tracePoly, startZ_osu, Zincr_osu, Xmax_osu, Xincr_osu,  showInPlusYbool, imgMag, backThickE_osu, epsilonB) =
      //sweep from top of image, placed at z=startZ_osu (aligned to imgarray y =0) to bottom of image.   
      //each poly outline sweeps as a wraparound, clockwize around Z. 
      //the image is in the x-z axis, with bright spots extending in showInPlusYbool ? +y : -y. 
      //back of image is at y=0. 
     [ for(y = [0:len(imgarray)-1])
         let ( NumX = len(imgarray[0]), 
               NumY = len(imgarray), 
               thisZrow = startZ_osu + Zincr_osu * y,
               sortPolyCrossings = LT_findPolyEdgesSortedXasc(tracePoly, thisZrow),
               minmaxX = LT_minmaxXatEdgewCrossings(sortPolyCrossings),
               atLeastOnePixinPoly = minmaxX[0] != -1,
               OneWrapright_to_left = atLeastOnePixinPoly ? 
                 (showInPlusYbool ?
                     PixRowCWposY(imgarray[y], thisZrow, sortPolyCrossings, Xincr_osu, Xmax_osu, minmaxX, showInPlusYbool, imgMag, backThickE_osu, epsilonB) 
                   : PixRowCWnegY(imgarray[y], thisZrow, sortPolyCrossings, Xincr_osu, Xmax_osu, minmaxX, showInPlusYbool, imgMag, backThickE_osu, epsilonB) 
                 )
               : []
              //, xx = thisZrow > 20 && thisZrow < 22 ? echo("LT_ImageInsidePoly: RIgrey_gendat: thisZrow, minmixatedge", thisZrow, minmaxX) :0
              )
         OneWrapright_to_left
     ]
     ;
     
   function stripNullRows(polySet, i=0, accumNonNull=[]) =
     i >= len(polySet) ? accumNonNull :
     let ( P = polySet[i],
        newAccum = len(P) == 0 ? accumNonNull :
          concat( accumNonNull,  [P])
     )
     stripNullRows(polySet, i+1, newAccum) ;
     
    //testP = [   [], [], [1, 2], [4,5], [] ]; 
    // echo("stripNullRows:", stripNullRows(testP)); 
     
   polys = RIgrey_gendat(ImageArray, trace_osu, maxYextMinval_osu, yIncr_osu, maxX_osu, xIncr_osu,  showInPlusYbool, imgMag, backThickExtended_osu, epsilonBack); 
   polysN = stripNullRows(polys);
   //echo("LT_ImageInsidePoly polygon number", len(polysN), "each poly len", len(polysN[0]) ); 
   //echo("LT_ImageInsidePoly row Z 0", polysN[0]);
   
   ReferenceLen = len(polysN[0]); 
   assert(len(polysN) > 4 && ReferenceLen > 4 ); 
   checkMatching = [ for(i=[0:len(polysN)-1]) if (len(polysN[i]) != ReferenceLen) i ] ;
   assert (len(checkMatching) == 0, str("one or more polys for polyhedron construction not of same length", checkMatching)); 
   
   
   cutThick = epsilonBack * 3.0; 
   cutBacking = showInPlusYbool ?  -cutThick - backThick_osu   : backThick_osu; 
   difference() {
      sweep(polysN, planar_caps = true);  
      translate([minX_osu-xIncr_osu, cutBacking, maxYextMinval_osu-yIncr_osu])
        cube([ xIncr_osu * (NumXpix+2), cutThick, yIncr_osu * (NumYpix+2)]);
   } 
} 

function LT_cropImage(img, skiplinesTop, skiplinesBot, skiprowsLeft, skiprowsRight) =
     assert( len(img[0]) > (skiprowsLeft+skiprowsRight) && len(img) > (skiplinesTop + skiplinesBot)) 
     [ for (y=[skiplinesTop: len(img)-1-skiplinesBot]) 
        [ for (x=[skiprowsLeft: len(img[0])-1-skiprowsRight])
            img[y][x]
        ]
     ] ;
        
function LT_flipImageTopToBottom(img) = 
     assert(len(img) > 2 && len(img[0]) > 2) 
     [ for (y=[len(img)-1: -1: 0])
        img[y] 
     ];
function LT_flipImageLeftToRight(img) =
     assert(len(img) > 2 && len(img[0]) > 2) 
     [ for (y=[0:len(img)-1])
        [ for (x = [len(img[0])-1 : -1 : 0])
           img[y][x]
        ]
     ];
 

        
        

//show image on rectangular plate, from greyscale text array
//  image stands in x-z plane.
module RectImage_fromgreyarray(ImageArray, xrectwidth, imgMag, stdBrightPix) {

    NumXpix = len(ImageArray[0]);  
    NumYpix = len(ImageArray);  
    //echo("RectImage_fromgreyarray:  x y", NumXpix, NumYpix);
    
    scaledBrightPixel = stdBrightPix * xrectwidth;   
    MagnifiedPixel =  scaledBrightPixel * imgMag; 
        
   function PixRowCW(imgarray, Y, numX, magnify) =
      //image is on negative Y side
     concat (
        [ for (x=[len(imgarray)-1 : -1: 0]) 
           let ( 
               )
           [ x, (-imgarray[x] /255.0 * magnify) -0.1, Y ] //shift a bit in -y so 
                //black pixels don't have 0 thickness ! 
    
        ],
        [ [0, 0, Y],  [numX, 0, Y]   ] //backside is at y=0. 
     );


   function RIgrey_gendat(imgarray, imgMag) =
      //sweep from top of image to bottom.   
      //each poly outline sweeps as a wraparound, clockwize around Z. 
      //the image is in the x-z axis, with bright spots extending in -y. 
      //back of image is at y=0. 
     [ for(y = [0:len(imgarray)-1])
         let ( NumX = len(imgarray[0]), 
               NumY = len(imgarray), 
               OneWrapright_to_left = PixRowCW(imgarray[NumY-y-1], y, NumX, imgMag)
              )
         OneWrapright_to_left
     ]
     ;
     
   //echo ("ImageArray 0 ", ImageArray[0]);
   capPixrow = [for (x = [0: len(ImageArray[0])-1])  128] ; 
   //echo ("capPixrow", capPixrow);  
   capImageArray = concat ( [capPixrow], ImageArray,[capPixrow] );  
   //echo ("capImageArray 0", capImageArray[0]);
   //echo ("capImageArray 1", capImageArray[1]);
   
   polys = RIgrey_gendat(capImageArray, MagnifiedPixel);
   //echo ("polys generated", polys[0]);
   //echo ("polys, last Y -1", polys[len(polys)-1-1] ); 
  // echo ("polys, last Y", polys[len(polys)-1] );  

    Xscaling = xrectwidth/NumXpix;
    scale([Xscaling, 1.0, Xscaling])
        sweep(polys, planar_caps=true);   
}