/********************************************************************************

   Fotoxx      edit photos and manage collections

   Copyright 2007-2016 Michael Cornelison
   Source URL: http://kornelix.net
   Contact: kornelix@posteo.de

   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/.

*********************************************************************************

   Fotoxx image edit - warp menu functions

   m_unbend                straighten curvature added by pano or improve perspective
   m_perspective           select a tetragon and convert into a rectangle
   m_warp_area             select an area and warp interior with mouse drags
   m_unwarp_CU             select a face in a close-up photo, remove distortion
   m_warp_curved           warp an image or area using a curved transform
   m_warp_linear           warp an image using a linear transform
   m_warp_affine           warp an image using an affine transform
   m_flatbook              flatten a photographed book page
   m_sphere                spherical image projection, variable radius
   m_selective_rescale     rescale while leaving selected areas unchanged
   m_waves                 distort an image with a wave pattern

*********************************************************************************/

#define EX extern                                                                //  enable extern declarations
#include "fotoxx.h"                                                              //  (variables in fotoxx.h are refs)

/********************************************************************************/


//  unbend an image
//  straighten curvature added by pano or improve perspective

float    unbend_lin_horz, unbend_lin_vert;                                       //  unbend values from dialog
float    unbend_cur_horz, unbend_cur_vert;
float    unbend_x1, unbend_x2, unbend_y1, unbend_y2;                             //  unbend axes scaled 0 to 1
int      unbend_hx1, unbend_hy1, unbend_hx2, unbend_hy2;
int      unbend_vx1, unbend_vy1, unbend_vx2, unbend_vy2;

editfunc    EFunbend;


//  menu function

void m_unbend(GtkWidget *, cchar *)                                              //  overhauled
{
   int    unbend_dialog_event(zdialog* zd, cchar *event);
   void * unbend_thread(void *);
   void   unbend_mousefunc();

   F1_help_topic = "unbend_image";

   EFunbend.funcname = "unbend";
   EFunbend.FprevReq = 1;                                                        //  use preview
   EFunbend.Frestart = 1;                                                        //  restart allowed
   EFunbend.threadfunc = unbend_thread;                                          //  thread function
   EFunbend.mousefunc = unbend_mousefunc;                                        //  mouse function
   if (! edit_setup(EFunbend)) return;                                           //  setup edit
   
   PXM_addalpha(E0pxm);                                                          //  15.09
   PXM_addalpha(E1pxm);
   PXM_addalpha(E3pxm);

/***                                                                             //  16.06
          _________________________________________
         |                Unbend                   |
         |   ____                ____              |
         |  |    | [____|-+]    |    | [____|-+]   |     splinvert   spcurvert
         |  |____|              |____|             |
         |                                         |
         |   ____                ____              |
         |  |    | [____|-+]    |    | [____|-+]   |     splinhorz   spcurhorz
         |  |____|              |____|             |
         |                                         |
         |   [ grid ]                              |
         |                         [done] [cancel] |
         |_________________________________________|

***/

   zdialog *zd = zdialog_new(ZTX("Unbend"),Mwin,Bdone,Bcancel,null);
   EFunbend.zd = zd;
   
   zdialog_add_widget(zd,"hbox","hb1","dialog",0,"space=3");
   zdialog_add_widget(zd,"icon","VL","hb1","unbend vert linear.png","size=64|space=5");
   zdialog_add_widget(zd,"spin","splinvert","hb1","-99|99|1|0");
   zdialog_add_widget(zd,"label","space","hb1",0,"space=10");
   zdialog_add_widget(zd,"icon","VC","hb1","unbend vert curved.png","size=64|space=5");
   zdialog_add_widget(zd,"spin","spcurvert","hb1","-99|99|1|0");

   zdialog_add_widget(zd,"hbox","hb2","dialog",0,"space=3");
   zdialog_add_widget(zd,"icon","HL","hb2","unbend horz linear.png","size=64|space=5");
   zdialog_add_widget(zd,"spin","splinhorz","hb2","-99|99|1|0");
   zdialog_add_widget(zd,"label","space","hb2",0,"space=10");
   zdialog_add_widget(zd,"icon","HC","hb2","unbend horz curved.png","size=64|space=5");
   zdialog_add_widget(zd,"spin","spcurhorz","hb2","-99|99|1|0");

   zdialog_add_widget(zd,"hbox","hb3","dialog",0,"space=5");
   zdialog_add_widget(zd,"button","grid","hb3",Bgrid,"space=10");


   unbend_x1 = unbend_x2 = unbend_y1 = unbend_y2 = 0.5;                          //  initial axes thru image middle
   unbend_lin_horz = unbend_lin_vert = 0;
   unbend_cur_horz = unbend_cur_vert = 0;

   currgrid = 2;                                                                 //  use unbend grid
   takeMouse(unbend_mousefunc,dragcursor);                                       //  connect mouse function
   signal_thread();
   zdialog_run(zd,unbend_dialog_event,"save");                                   //  run dialog, parallel
   return;
}


//  dialog event and completion callback function

int unbend_dialog_event(zdialog *zd, cchar *event)
{
   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (zd->zstat)                                                                //  dialog complete
   {
      currgrid = 0;                                                              //  restore normal grid settings       15.09

      if (zd->zstat == 1)
      {                                                                          //  done
         CEF->Fmods = 0;
         if (unbend_cur_vert || unbend_cur_horz ||                               //  image3 modified
             unbend_lin_vert || unbend_lin_horz) {
            CEF->Fmods = 1;
            CEF->Fsaved = 0;                                                     //  done
         }
         edit_fullsize();                                                        //  get full size image
         signal_thread();
         edit_done(0);                                                           //  commit edit
      }
      else edit_cancel(0);                                                       //  discard edit
      draw_toplines(2);                                                          //  erase axes-lines
      return 1;
   }

   if (strstr(event,"splinvert")) {                                              //  get new unbend value
      zdialog_fetch(zd,"splinvert",unbend_lin_vert);
      signal_thread();                                                           //  trigger thread
   }

   if (strstr(event,"splinhorz")) {
      zdialog_fetch(zd,"splinhorz",unbend_lin_horz);
      signal_thread();
   }

   if (strstr(event,"spcurvert")) {
      zdialog_fetch(zd,"spcurvert",unbend_cur_vert);
      signal_thread();
   }

   if (strstr(event,"spcurhorz")) {
      zdialog_fetch(zd,"spcurhorz",unbend_cur_horz);
      signal_thread();
   }

   if (strmatch(event,"grid")) m_gridlines(0,"grid 2");                          //  grid settings dialog

   return 1;
}


//  unbend mouse function                                                        //  adjustable axes

void unbend_mousefunc()
{
   cchar       *close;
   float       dist1, dist2;
   float       mpx = 0, mpy = 0;

   if (LMclick) {                                                                //  left mouse click
      mpx = Mxclick;
      mpy = Myclick;
   }

   if (Mxdrag || Mydrag) {                                                       //  mouse dragged
      mpx = Mxdrag;
      mpy = Mydrag;
   }

   if (! mpx && ! mpy) return;

   mpx = 1.0 * mpx / E3pxm->ww;                                                  //  scale mouse position 0 to 1
   mpy = 1.0 * mpy / E3pxm->hh;

   if (mpx < 0.2 || mpx > 0.8 ) {                                                //  check reasonable position
      if (mpy < 0.1 || mpy > 0.9) return;
   }
   else if (mpy < 0.2 || mpy > 0.8) {
      if (mpx < 0.1 || mpx > 0.9) return;
   }
   else return;

   close = "?";                                                                  //  find closest axis end-point
   dist1 = 2;

   dist2 = mpx * mpx + (mpy-unbend_y1) * (mpy-unbend_y1);
   if (dist2 < dist1) {
      dist1 = dist2;
      close = "left";
   }

   dist2 = (1-mpx) * (1-mpx) + (mpy-unbend_y2) * (mpy-unbend_y2);
   if (dist2 < dist1) {
      dist1 = dist2;
      close = "right";
   }

   dist2 = (mpx-unbend_x1) * (mpx-unbend_x1) + mpy * mpy;
   if (dist2 < dist1) {
      dist1 = dist2;
      close = "top";
   }

   dist2 = (mpx-unbend_x2) * (mpx-unbend_x2) + (1-mpy) * (1-mpy);
   if (dist2 < dist1) {
      dist1 = dist2;
      close = "bottom";
   }

   if (strmatch(close,"left")) unbend_y1 = mpy;                                  //  set new axis end-point
   if (strmatch(close,"right")) unbend_y2 = mpy;
   if (strmatch(close,"top")) unbend_x1 = mpx;
   if (strmatch(close,"bottom")) unbend_x2 = mpx;

   signal_thread();                                                              //  trigger thread

   LMclick = Mxdrag = Mydrag = 0;
   return;
}


//  unbend thread function

void * unbend_thread(void *arg)
{
   void * unbend_wthread(void *);

   while (true)
   {
      thread_idle_loop();                                                        //  wait for work or exit request

      unbend_hx1 = 0;                                                            //  scale axes to E3ww/hh
      unbend_hy1 = unbend_y1 * E3pxm->hh;
      unbend_hx2 = E3pxm->ww;
      unbend_hy2 = unbend_y2 * E3pxm->hh;

      unbend_vx1 = unbend_x1 * E3pxm->ww;
      unbend_vy1 = 0;
      unbend_vx2 = unbend_x2 * E3pxm->ww;
      unbend_vy2 = E3pxm->hh;

      Ntoplines = 2;
      toplines[0].x1 = unbend_hx1;                                               //  lines on window
      toplines[0].y1 = unbend_hy1;
      toplines[0].x2 = unbend_hx2;
      toplines[0].y2 = unbend_hy2;
      toplines[0].type = 2;
      toplines[1].x1 = unbend_vx1;
      toplines[1].y1 = unbend_vy1;
      toplines[1].x2 = unbend_vx2;
      toplines[1].y2 = unbend_vy2;
      toplines[1].type = 2;

      for (int ii = 0; ii < NWT; ii++)                                           //  start worker threads
         start_wthread(unbend_wthread,&Nval[ii]);
      wait_wthreads();                                                           //  wait for completion

      Fpaint2();                                                                 //  update window
   }

   return 0;                                                                     //  not executed, stop g++ warning
}


void * unbend_wthread(void *arg)                                                 //  worker thread function
{
   int         index = *((int *) arg);
   int         vstat, px3, py3, cx3, cy3;
   float       dispx, dispy, dispx2, dispy2;
   float       px1, py1, vx1, vx2, hy1, hy2;
   float       curvert, curhorz, linvert, linhorz;
   float       vpix[4], *pix3;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);

   curvert = int(unbend_cur_vert * 0.01 * E3pxm->hh);                            //  -0.99 to +0.99
   curhorz = int(unbend_cur_horz * 0.01 * E3pxm->ww);
   linvert = int(unbend_lin_vert * 0.0013 * E3pxm->hh);                          //  -0.13 to +0.13
   linhorz = int(unbend_lin_horz * 0.0013 * E3pxm->ww);

   vx1 = unbend_vx1;
   vx2 = unbend_vx2;
   hy1 = unbend_hy1;
   hy2 = unbend_hy2;
   
   for (py3 = index; py3 < E3pxm->hh; py3 += NWT)                                //  step through F3 pixels
   for (px3 = 0; px3 < E3pxm->ww; px3++)
   {
      cx3 = vx1 + (vx2 - vx1) * py3 / E3pxm->hh;                                 //  center of unbend
      cy3 = hy1 + (hy2 - hy1) * px3 / E3pxm->ww;
      dispx = 2.0 * (px3 - cx3) / E3pxm->ww;                                     //  -1.0 ..  0.0 .. +1.0  (roughly)
      dispy = 2.0 * (py3 - cy3) / E3pxm->hh;                                     //  -1.0 ..  0.0 .. +1.0
      dispx2 = -cosf(0.8 * dispx) + 1;                                           //   curved
      dispy2 = -cosf(0.8 * dispy) + 1;

      pix3 = PXMpix(E3pxm,px3,py3);                                              //  output pixel
      px1 = px3;                                                                 //  input pixel = output
      py1 = py3;

      px1 += dispx * dispy * linhorz;                                            //  move input pixel
      py1 += dispy * dispx * linvert;
      px1 -= dispx * dispy2 * curhorz;                                           //  change sign                        16.06
      py1 -= dispy * dispx2 * curvert;

      vstat = vpixel(E1pxm,px1,py1,vpix);
      if (vstat) memcpy(pix3,vpix,pcc);
      else memset(pix3,0,pcc);
   }

   exit_wthread();
   return 0;                                                                     //  not executed, avoid gcc warning
}


/********************************************************************************/

//  Convert a selected tetragon area into a rectangle, converting the
//  rest of the image to match and keeping straight lines straight.

namespace perspective
{
   int         PSP_pixel[4][2];                                                  //  last 0-4 pixels clicked
   int         PSP_npix;                                                         //  count of pixels
   char        PSP_pixlab[4][4] = { " A ", " B ", " C ", " D " };
   int         PSP_corner = 0;

   editfunc    EFperspective;

   int   dialog_event(zdialog *zd, cchar *event);
   void  mousefunc(void);
   void  warpfunc(void);
   void  KBfunc(int key);
}


//  menu function

void m_perspective(GtkWidget *, cchar *)
{
   using namespace perspective;

   cchar  *PSP_message = ZTX(
          " Click the four corners of a tetragon area. Press [apply]. \n"
          " The image is warped to make the tetragon into a rectangle.");

   F1_help_topic = "fix_perspective";

   EFperspective.menufunc = m_perspective;
   EFperspective.funcname = "perspective";
   EFperspective.Frestart = 1;                                                   //  restart allowed
   EFperspective.mousefunc = mousefunc;                                          //  mouse function
   if (! edit_setup(EFperspective)) return;                                      //  setup edit

   PXM_addalpha(E0pxm);                                                          //  15.09
   PXM_addalpha(E1pxm);
   PXM_addalpha(E3pxm);

   PSP_npix = 0;                                                                 //  no pixels yet

   zdialog *zd = zdialog_new(ZTX("Perspective Correction"),Mwin,Bapply,Breset,Btrim,Bdone,null);
   zdialog_add_widget(zd,"label","lab1","dialog",PSP_message,"space=3");

   EFperspective.zd = zd;
   zdialog_run(zd,dialog_event,"save");                                          //  run dialog, parallel

   takeMouse(mousefunc,dragcursor);                                              //  connect mouse function
   return;
}


//  dialog completion callback function

int perspective::dialog_event(zdialog *zd, cchar *event)
{
   using namespace perspective;

   int      ii, px, py;
   int      minx, maxx, miny, maxy;

   if (strmatch(event,"done")) zd->zstat = 3;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 3;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 4;                                  //  from f_open()

   if (strmatch(event,"focus"))                                                  //  toggle mouse capture   12.01
      takeMouse(mousefunc,dragcursor);

   if (! zd->zstat) return 1;                                                    //  wait for completion

   if (zd->zstat == 1)                                                           //  apply
   {
      erase_toptext(102);                                                        //  erase points
      warpfunc();                                                                //  do the warp
      zd->zstat = 0;                                                             //  keep dialog active
   }

   else if (zd->zstat == 2)                                                      //  reset
   {
      edit_reset();
      zd->zstat = 0;

      for (ii = 0; ii < PSP_npix; ii++)                                          //  show pixel labels on image
      {
         px = PSP_pixel[ii][0];
         py = PSP_pixel[ii][1];
         add_toptext(102,px,py,PSP_pixlab[ii],"Sans 8");
      }
      Fpaint2();
   }
   
   else if (zd->zstat == 3)                                                      //  trim                         15.12
   {
      if (PSP_npix < 4) {                                                        //  wait for 4 corners
         zd->zstat = 0;
         return 1;
      }
      
      erase_toptext(102);                                                        //  erase labels
      warpfunc();                                                                //  do the warp
      edit_done(0);

      minx = miny = 99999;                                                       //  find trim limits
      maxx = maxy = 0;
      
      for (ii = 0; ii < 4; ii++) {
         if (PSP_pixel[ii][0] < minx) minx = PSP_pixel[ii][0];
         if (PSP_pixel[ii][0] > maxx) maxx = PSP_pixel[ii][0];
         if (PSP_pixel[ii][1] < miny) miny = PSP_pixel[ii][1];
         if (PSP_pixel[ii][1] > maxy) maxy = PSP_pixel[ii][1];
      }

      trimrect[0] = minx;                                                        //  set parameters for trim function
      trimrect[2] = maxx;
      trimrect[1] = miny;
      trimrect[3] = maxy;
      m_trimrotate(0,"keep");
   }

   else if (zd->zstat == 4) {                                                    //  done
      erase_toptext(102);                                                        //  erase labels
      edit_done(0);                                                              //  commit edit
   }

   else {                                                                        //  cancel
      erase_toptext(102);                                                        //  erase labels
      edit_cancel(0);                                                            //  discard edit
   }

   return 1;
}


//  mouse function - click on 4 corners of tetragon

void perspective::mousefunc(void)
{
   using namespace perspective;

   int         ii, minii, jj, px, py;
   float       dist, distx, disty, mindist;

   if (LMclick)                                                                  //  left click
   {
      for (ii = 0; ii < PSP_npix; ii++)                                          //  check if very near a previous corner
      {
         if (abs(PSP_pixel[ii][0] - Mxclick) < 0.07 * E3pxm->ww
          && abs(PSP_pixel[ii][1] - Myclick) < 0.07 * E3pxm->hh)
         {
            PSP_pixel[ii][0] = Mxclick;                                          //  yes, set new corner position
            PSP_pixel[ii][1] = Myclick;
            PSP_corner = ii;
            goto showcorners;
         }
      }

      if (PSP_npix < 4)                                                          //  if < 4 corners, add a new one
      {
         ii = PSP_npix;                                                          //  next corner to fill
         PSP_pixel[ii][0] = Mxclick;                                             //  save newest corner position
         PSP_pixel[ii][1] = Myclick;
         PSP_corner = ii;
         PSP_npix++;
         goto showcorners;
      }

      mindist = 99999;                                                           //  all 4 corners already specified
      minii = -1;

      for (ii = 0; ii < 4; ii++)                                                 //  find closest corner to click position
      {
         distx = (Mxclick - PSP_pixel[ii][0]);
         disty = (Myclick - PSP_pixel[ii][1]);
         dist = sqrt(distx*distx + disty*disty);
         if (dist < mindist) {
            mindist = dist;
            minii = ii;
         }
      }

      if (minii >= 0) {                                                          //  set new corner position
         ii = minii;
         PSP_pixel[ii][0] = Mxclick;
         PSP_pixel[ii][1] = Myclick;
         PSP_corner = ii;
         goto showcorners;
      }
   }

   else if (RMclick)                                                             //  right click
   {
      mindist = 99999;
      minii = -1;

      for (ii = 0; ii < PSP_npix; ii++)                                          //  find closest corner to click position
      {
         distx = (Mxclick - PSP_pixel[ii][0]);
         disty = (Myclick - PSP_pixel[ii][1]);
         dist = sqrt(distx*distx + disty*disty);
         if (dist < mindist) {
            mindist = dist;
            minii = ii;
         }
      }

      if (minii >= 0) {                                                          //  replace deleted corner with
         ii = minii;                                                             //    last corner
         jj = PSP_npix - 1;
         PSP_pixel[ii][0] = PSP_pixel[jj][0];
         PSP_pixel[ii][1] = PSP_pixel[jj][1];
         PSP_corner = ii;
         --PSP_npix;                                                             //  reduce corner count
         goto showcorners;
      }
   }

showcorners:                                                                     //  show corner labels on image

   erase_toptext(102);

   for (ii = 0; ii < PSP_npix; ii++)
   {
      px = PSP_pixel[ii][0];
      py = PSP_pixel[ii][1];
      add_toptext(102,px,py,PSP_pixlab[ii],"Sans 8");
   }

   LMclick = RMclick = 0;
   Fpaint2();
   return;
}


//  Keyboard function
//  KB arrow keys tweak position of last touched tetragon corner.

void perspective::KBfunc(int key)
{
   using namespace perspective;

   int      ii, xstep, ystep;

   xstep = ystep = 0;
   if (key == GDK_KEY_Left) xstep = -1;
   if (key == GDK_KEY_Right) xstep = +1;
   if (key == GDK_KEY_Up) ystep = -1;
   if (key == GDK_KEY_Down) ystep = +1;

   ii = PSP_corner;
   if (ii < 0 || ii > 3) return;                                                 //  last corner touched, 0-3

   PSP_pixel[ii][0] += xstep;
   PSP_pixel[ii][1] += ystep;

   mousefunc();
   return;
}


//  perspective warp function - make input tetragon into a rectangle

void perspective::warpfunc(void)
{
   using namespace perspective;

   int         ii, jj, tempx, tempy, vstat;
   float       px3, py3, trpx[4], trpy[4];
   float       sqpx0, sqpy0, sqpx1, sqpy1, sqpx2, sqpy2, sqpx3, sqpy3;
   float       cdx0, cdy0, cdx1, cdy1, cdx2, cdy2, cdx3, cdy3;
   float       px1, py1, dispx, dispy, sqww, sqhh;
   float       f0, f1, f2, f3;
   float       vpix1[4], *pix3;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);

   if (PSP_npix != 4) {
      zmessageACK(Mwin,ZTX("must have 4 corners"));
      return;
   }

   for (ii = 0; ii < 4; ii++) {                                                  //  get 4 selected tetragon points
      trpx[ii] = PSP_pixel[ii][0];
      trpy[ii] = PSP_pixel[ii][1];
   }

   //  sort 4 points in clockwise order NW, NE, SE, SW

   for (ii = 0; ii < 4; ii++) {                                                  //  sort top to bottom (y order)
      for (jj = ii; jj < 4; jj++) {
         if (trpy[jj] < trpy[ii]) {
            tempx = trpx[ii];
            tempy = trpy[ii];
            trpx[ii] = trpx[jj];
            trpy[ii] = trpy[jj];
            trpx[jj] = tempx;
            trpy[jj] = tempy;
         }
      }
   }

   if (trpx[1] < trpx[0]) {                                                      //  sort upper two left, right
      tempx = trpx[0];
      tempy = trpy[0];
      trpx[0] = trpx[1];
      trpy[0] = trpy[1];
      trpx[1] = tempx;
      trpy[1] = tempy;
   }

   if (trpx[2] < trpx[3]) {                                                      //  sort lower two right, left
      tempx = trpx[2];
      tempy = trpy[2];
      trpx[2] = trpx[3];
      trpy[2] = trpy[3];
      trpx[3] = tempx;
      trpy[3] = tempy;
   }

   if (trpx[0] < trpx[3]) sqpx0 = sqpx3 = trpx[0];                               //  rectangle enclosing tetragon
   else  sqpx0 = sqpx3 = trpx[3];
   if (trpx[1] > trpx[2]) sqpx1 = sqpx2 = trpx[1];
   else  sqpx1 = sqpx2 = trpx[2];
   if (trpy[0] < trpy[1]) sqpy0 = sqpy1 = trpy[0];
   else  sqpy0 = sqpy1 = trpy[1];
   if (trpy[2] > trpy[3]) sqpy2 = sqpy3 = trpy[2];
   else  sqpy2 = sqpy3 = trpy[3];

/***
   sqpx0 = sqpx3 = 0.5 * (trpx[0] + trpx[3]);                                    //  rectangle bisecting tetragon sides
   sqpx1 = sqpx2 = 0.5 * (trpx[1] + trpx[2]);
   sqpy0 = sqpy1 = 0.5 * (trpy[0] + trpy[1]);
   sqpy2 = sqpy3 = 0.5 * (trpy[2] + trpy[3]);
***/

   cdx0 = sqpx0 - trpx[0];                                                       //  displavement of tetragon corner
   cdy0 = sqpy0 - trpy[0];                                                       //    to corresponding rectangle corner
   cdx1 = sqpx1 - trpx[1];
   cdy1 = sqpy1 - trpy[1];
   cdx2 = sqpx2 - trpx[2];
   cdy2 = sqpy2 - trpy[2];
   cdx3 = sqpx3 - trpx[3];
   cdy3 = sqpy3 - trpy[3];

   sqww = 1.0 / (sqpx1 - sqpx0);                                                 //  rectangle width and height
   sqhh = 1.0 / (sqpy3 - sqpy0);

   for (py3 = 0; py3 < E3pxm->hh; py3++)                                         //  loop all output pixels
   for (px3 = 0; px3 < E3pxm->ww; px3++)
   {
      f0 = (1.0 - (px3 - sqpx0) * sqww) * (1.0 - (py3 - sqpy0) * sqhh);
      f1 = (px3 - sqpx0) * sqww * (1.0 - (py3 - sqpy0) * sqhh);
      f2 = (px3 - sqpx0) * sqww * (py3 - sqpy0) * sqhh;
      f3 = (1.0 - (px3 - sqpx0) * sqww) * (py3 - sqpy0) * sqhh;

      dispx = cdx0 * f0 + cdx1 * f1 + cdx2 * f2 + cdx3 * f3;
      dispy = cdy0 * f0 + cdy1 * f1 + cdy2 * f2 + cdy3 * f3;

      px1 = px3 - dispx;                                                         //  input virtual pixel for px3/py3
      py1 = py3 - dispy;

      pix3 = PXMpix(E3pxm,int(px3),int(py3));                                    //  output pixel

      vstat = vpixel(E1pxm,px1,py1,vpix1);                                       //  output pixel = input virtual pixel
      if (vstat) memcpy(pix3,vpix1,pcc);
      else memset(pix3,0,pcc);                                                   //  voided pixel
   }

   CEF->Fmods++;                                                                 //  image is modified
   CEF->Fsaved = 0;
   Fpaint2();                                                                    //  update window
   return;
}


/********************************************************************************/

//  warp/distort area - select image area and pull with mouse

float       *WarpAx, *WarpAy;                                                    //  memory for pixel warp vectors
int         WarpAcc, WarpAnew;
int         WarpA_started;
int         WarpA_areanumber;

editfunc    EFwarpA;

void  WarpA_init();
void  WarpA_warpfunc(float mdx, float mdy, float mdw, float mdh, int acc);
void  WarpA_warpfunc2(float mdx, float mdy, float mdw, float mdh, int acc);
void  WarpA_mousefunc();
void  WarpA_expand();
void  WarpA_edgeblend();


//  menu function

void m_warp_area(GtkWidget *, cchar *)
{
   int      WarpA_dialog_event(zdialog *zd, cchar *event);

   cchar  *WarpA_message = ZTX(
             " Select an area to warp using select area function. \n"
             " Press [start warp] and pull area with mouse. \n"
             " Make multiple mouse pulls until satisfied. \n"
             " When finished, select another area or press [done].");

   F1_help_topic = "warp_area";

   EFwarpA.menufunc = m_warp_area;
   EFwarpA.funcname = "warp-area";
   EFwarpA.Farea = 2;                                                            //  select area usable
   EFwarpA.mousefunc = WarpA_mousefunc;                                          //  mouse function
   if (! edit_setup(EFwarpA)) return;                                            //  setup edit

   zdialog *zd = zdialog_new(ZTX("Warp area"),Mwin,Bdone,Bcancel,null);
   EFwarpA.zd = zd;

   zdialog_add_widget(zd,"label","lab1","dialog",WarpA_message,"space=3");
   zdialog_add_widget(zd,"hbox","hb1","dialog",0,"space=5");
   zdialog_add_widget(zd,"button","start","hb1",ZTX("start warp"),"space=5");
   zdialog_add_widget(zd,"button","reset","hb1",Breset,"space=5");

   WarpAcc = E3pxm->ww * E3pxm->hh * sizeof(float);
   WarpAx = (float *) zmalloc(WarpAcc);                                          //  get memory for pixel warp vectors
   WarpAy = (float *) zmalloc(WarpAcc);
   memset(WarpAx,0,WarpAcc);
   memset(WarpAy,0,WarpAcc);

   WarpA_started = 0;
   WarpA_areanumber = 0;
   WarpAnew = 0;

   zdialog_run(zd,WarpA_dialog_event,"save");                                    //  run dialog, parallel
   return;
}


//  dialog event and completion callback function

int WarpA_dialog_event(zdialog * zd, cchar *event)
{
   int      init = 0;

   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (zd->zstat)
   {
      if (zd->zstat == 1) edit_done(0);                                          //  commit edit
      else edit_cancel(0);                                                       //  discard edit

      zfree(WarpAx);                                                             //  release memory
      zfree(WarpAy);
      return 1;
   }

   if (! sa_validate()) init = 1;                                                //  area invalid for curr. image file

   if (sa_stat != 3 || sa_mode == mode_image)                                    //  no select area active
      init = 1;

   if (WarpA_started && WarpA_areanumber != areanumber)                          //  select area changed
      init = 1;

   if (init) {
      memset(WarpAx,0,WarpAcc);
      memset(WarpAy,0,WarpAcc);
      WarpA_started = 0;
      WarpA_areanumber = 0;
      return 1;
   }

   if (strmatch(event,"start"))                                                  //  start warp
   {
      if (sa_stat != 3 || sa_mode == mode_image) {                               //  no select area active
         zmessageACK(Mwin,ZTX("no active Select Area"));
         return 1;
      }
      sa_edgecalc();                                                             //  calculate area edge distances
      takeMouse(WarpA_mousefunc,dragcursor);                                     //  connect mouse function
      WarpA_started = 1;
      WarpA_areanumber = areanumber;
   }

   if (strmatch(event,"focus"))                                                  //  toggle mouse capture
      takeMouse(WarpA_mousefunc,dragcursor);

   if (strmatch(event,"reset")) {                                                //  undo all warps
      edit_reset();
      memset(WarpAx,0,WarpAcc);                                                  //  clear warp vectors
      memset(WarpAy,0,WarpAcc);
   }

   if (strmatch(event,"blendwidth"))
      WarpA_edgeblend();

   return 1;
}


//  warp mouse function

void  WarpA_mousefunc()
{
   static float   mdx, mdy, mdw, mdh;
   static int     warped = 0;

   if (! WarpA_started) return;

   if (Mxdrag || Mydrag)                                                         //  mouse drag underway
   {
      mdx = Mxdown;                                                              //  drag origin, image coordinates
      mdy = Mydown;
      mdw = Mxdrag - Mxdown;                                                     //  drag increment
      mdh = Mydrag - Mydown;
      WarpA_warpfunc(mdx,mdy,mdw,mdh,0);                                         //  warp image
      warped = 1;
      Mxdrag = Mydrag = 0;
      return;
   }

   else if (warped)
   {
      warped = 0;
      WarpA_warpfunc(mdx,mdy,mdw,mdh,1);                                         //  drag done, update warp vectors
      if (WarpAnew) WarpA_expand();
   }

   return;
}


//  warp image and accumulate warp memory
//  mdx/y = mouse initial position
//  mdw/h = mouse drag vector

void WarpA_warpfunc(float mdx, float mdy, float mdw, float mdh, int acc)
{
   int         ii, px, py, ww, hh, vstat;
   float       ddx, ddy, dpe, dpm1, dpm2, dpm;
   float       mag, dispx, dispy;
   float       vpix[4], *pix3;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);

   if (sa_stat != 3) return;                                                     //  area erased

   ww = E1pxm->ww;
   hh = E1pxm->hh;

   ii = mdy * ww + mdx;
   if (! sa_pixmap[ii]) {                                                        //  if mouse drag outside select area
      WarpA_warpfunc2(mdx,mdy,mdw,mdh,acc);                                      //    then use alternate function
      return;
   }

   for (py = sa_miny; py < sa_maxy; py++)                                        //  loop all pixels in area
   for (px = sa_minx; px < sa_maxx; px++)
   {
      ii = py * ww + px;
      dpe = sa_pixmap[ii];                                                       //  distance from area edge
      if (dpe < 1) continue;                                                     //  outside area
      dpe -= 1;                                                                  //  rebase edge = 0

      ddx = (px - mdx);                                                          //  pixel distance to mouse position
      ddy = (py - mdy);                                                          //    before drag
      dpm1 = sqrt(ddx*ddx + ddy*ddy);

      ddx -= mdw;                                                                //  pixel distance to mouse position
      ddy -= mdh;                                                                //    after drag
      dpm2 = sqrt(ddx*ddx + ddy*ddy);

      dpm = 0.5 * (dpm1 + dpm2);                                                 //  mean

      mag = (dpe / (dpm + dpe));                                                 //  1...0 for pixel at mouse...edge

      dispx = -mdw * mag;                                                        //  pixel movement from drag movement
      dispy = -mdh * mag;

      dispx += WarpAx[ii];                                                       //  add this warp to prior
      dispy += WarpAy[ii];

      if (acc) {                                                                 //  mouse drag done,
         WarpAx[ii] = dispx;                                                     //    accumulate warp memory
         WarpAy[ii] = dispy;
         continue;
      }

      vstat = vpixel(E1pxm,px+dispx,py+dispy,vpix);                              //  input virtual pixel
      if (vstat) {
         pix3 = PXMpix(E3pxm,px,py);                                             //  output pixel
         memcpy(pix3,vpix,pcc);
      }
   }

   ww = sa_maxx - sa_minx;                                                       //  update window
   hh = sa_maxy - sa_miny;
   Fpaint3(sa_minx,sa_miny,ww,hh);

   CEF->Fmods++;
   CEF->Fsaved = 0;
   return;
}


//  warp function when mouse is dragged outside select area

void  WarpA_warpfunc2(float mdx, float mdy, float mdw, float mdh, int acc)
{
   int         ii, jj, px, py, dpx, dpy, dpe;
   int         xlo, xhi, ylo, yhi;
   int         ww, hh, vstat;
   float       ddx, ddy, dpm1, dpm2, dpm;
   float       mag, dispx, dispy;
   float       mdmax, mdmax2;
   float       vpix[4], *pix1, *pix3;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);

   #define distance(px1,py1,px2,py2) \
            sqrtf((px2-px1)*(px2-px1)+(py2-py1)*(py2-py1))

   WarpAnew = 1;                                                                 //  note select area will expand

   ww = E1pxm->ww;
   hh = E1pxm->hh;

   ylo = sa_miny;                                                                //  select area enclosing rectangle
   yhi = sa_maxy;
   xlo = sa_minx;
   xhi = sa_maxx;

   if (mdy < ylo) ylo = mdy - 200;                                               //  extend area to mouse posn + more
   if (mdy > yhi) yhi = mdy + 200;
   if (mdx < xlo) xlo = mdx - 200;
   if (mdx > xhi) xhi = mdx + 200;

   if (ylo < 0) ylo = 0;
   if (yhi > hh) yhi = hh;
   if (xlo < 0) xlo = 0;
   if (xhi > ww) xhi = ww;

   mdmax = distance(mdx,mdy,xlo,ylo);                                            //  find max. mouse distance from
   mdmax2 = distance(mdx,mdy,xhi,ylo);                                           //    enclosing rectangle corners
   if (mdmax2 > mdmax) mdmax = mdmax2;
   mdmax2 = distance(mdx,mdy,xhi,yhi);
   if (mdmax2 > mdmax) mdmax = mdmax2;
   mdmax2 = distance(mdx,mdy,xlo,yhi);
   if (mdmax2 > mdmax) mdmax = mdmax2;

   for (py = ylo; py < yhi; py++)                                                //  loop all pixels in enclosing rectangle
   for (px = xlo; px < xhi; px++)
   {
      ddx = (px - mdx);                                                          //  pixel distance to mouse initial position
      ddy = (py - mdy);
      dpm1 = sqrt(ddx*ddx + ddy*ddy);

      ddx -= mdw;                                                                //  pixel distance to mouse curr. position
      ddy -= mdh;
      dpm2 = sqrt(ddx*ddx + ddy*ddy);

      dpm = 0.5 * (dpm1 + dpm2);                                                 //  mean

      mag = 1.0 - dpm / mdmax;                                                   //  1...0 for pixel-mouse 0...mdmax

      dispx = -mdw * mag;
      dispy = -mdh * mag;                                                        //  pixel movement from drag movement

      ii = py * ww + px;
      dispx += WarpAx[ii];                                                       //  add this warp to prior
      dispy += WarpAy[ii];

      dpx = px + dispx;                                                          //  source virtual pixel
      dpy = py + dispy;                                                          //  (nearest real pixel)
      if (dpx < 0) dpx = 0;
      if (dpx > ww-1) dpx = ww-1;
      if (dpy < 0) dpy = 0;
      if (dpy > hh-1) dpy = hh-1;
      jj = dpy * ww + dpx;
      dpe = sa_pixmap[jj];                                                       //  distance from area edge

      if (dpe < 1)                                                               //  outside area
         dispx = dispy = 0;

      if (acc) {                                                                 //  mouse drag done,
         WarpAx[ii] = dispx;                                                     //    accumulate warp memory
         WarpAy[ii] = dispy;
         continue;
      }

      if (dpe) {                                                                 //  source pixel inside area
         vstat = vpixel(E1pxm,px+dispx,py+dispy,vpix);                           //   = input virtual pixel
         if (vstat) {
            pix3 = PXMpix(E3pxm,px,py);                                          //  output pixel
            memcpy(pix3,vpix,pcc);
            continue;
         }
      }

      pix1 = PXMpix(E1pxm,px,py);                                                //  pixel is unchanged
      pix3 = PXMpix(E3pxm,px,py);
      memcpy(pix3,pix1,pcc); 
   }

   ww = xhi - xlo;                                                               //  update window
   hh = yhi - ylo;
   Fpaint3(xlo,ylo,ww,hh);

   CEF->Fmods++;
   CEF->Fsaved = 0;
   return;
}


//  expand select area if pulled outside the original bounds

void  WarpA_expand()
{
   int      ii, px, py, ww, hh;
   int      blend = sa_blend;

   ww = E1pxm->ww;
   hh = E1pxm->hh;

   for (py = 0; py < hh; py++)                                                   //  loop all pixels
   for (px = 0; px < ww; px++)
   {
      ii = py * ww + px;                                                         //  find pixels that have moved
      if (! WarpAx[ii] && ! WarpAy[ii]) continue;
      sa_pixmap[ii] = 1;                                                         //  mark within area
   }

   sa_map_pixels();
   sa_finish_auto();
   sa_edgecalc();                                                                //  recalculate edge distances

   sa_blend = blend;
   WarpA_edgeblend();

   WarpA_areanumber = areanumber;
   WarpAnew = 0;

   return;
}


//  blend area edges according to sa_blend

void  WarpA_edgeblend()
{
   int      ii, px, py, ww, hh, dist, vstat;
   float    *pix1, *pix3, vpix[4];
   float    dold, dnew, dispx, dispy;

   if (! sa_blend) return;

   ww = E1pxm->ww;
   hh = E1pxm->hh;

   for (py = 0; py < hh; py++)                                                   //  loop output pixels
   for (px = 0; px < ww; px++)
   {
      ii = py * ww + px;                                                         //  blend changes at edge
      dist = sa_pixmap[ii];
      if (! dist) continue;                                                      //  outside area
      if (dist > sa_blend) continue;                                             //  beyond blend area

      dispx = WarpAx[ii];                                                        //  warp
      dispy = WarpAy[ii];
      vstat = vpixel(E1pxm,px+dispx,py+dispy,vpix);                              //  input virtual pixel (new)
      if (! vstat) continue;
      pix1 = PXMpix(E1pxm,px,py);                                                //  input pixel (old)
      pix3 = PXMpix(E3pxm,px,py);                                                //  output pixel
      dnew = sa_blendfunc(dist);                                                 //  16.08
      dold = 1.0 - dnew;
      pix3[0] = dnew * vpix[0] + dold * pix1[0];
      pix3[1] = dnew * vpix[1] + dold * pix1[1];
      pix3[2] = dnew * vpix[2] + dold * pix1[2];
   }

   ww = sa_maxx - sa_minx;                                                       //  update window
   hh = sa_maxy - sa_miny;
   Fpaint3(sa_minx,sa_miny,ww,hh);

   return;
}


/********************************************************************************/

//  Unwarp closeup face photo - shrink magnified areas closest to camera

int         unwarpCU_started;
int         unwarpCU_areanumber;
float       unwarpCU_warpval;
int         unwarpCU_cx, unwarpCU_cy;
editfunc    EFunwarpCU;

int   unwarpCU_dialog_event(zdialog *zd, cchar *event);
void  unwarpCU_warpfunc();
void  unwarpCU_mousefunc();

//  menu function

void m_unwarp_closeup(GtkWidget *, cchar *)                                      //  16.11
{
   cchar  *unwarpCU_message = ZTX(
             " Use Select Area to select a face. \n"
             " Click on the tip of the nose. \n"
             " Move the slider. \n");

   F1_help_topic = "unwarp_closeup";

   EFunwarpCU.menufunc = m_unwarp_closeup;
   EFunwarpCU.funcname = "unwarp-closeup";
   EFunwarpCU.Farea = 2;                                                         //  select area usable
   EFunwarpCU.mousefunc = unwarpCU_mousefunc;                                    //  mouse function
   if (! edit_setup(EFunwarpCU)) return;                                         //  setup edit

/***
       ______________________________________
      |           Unwarp Closeup             |
      |                                      |
      |  Use Select Area to select a face.   |
      |  Click on the tip of the nose.       |
      |  Move the slider.                    |
      |                                      |
      |  [===============[]================] |
      |                                      |
      |                      [done] [cancel] |
      |______________________________________|

***/

   zdialog *zd = zdialog_new(ZTX("Unwarp Closeup"),Mwin,Bdone,Bcancel,null);
   EFunwarpCU.zd = zd;

   zdialog_add_widget(zd,"label","lab1","dialog",unwarpCU_message,"space=3");
   zdialog_add_widget(zd,"hbox","hbw","dialog",0,"space=5");
   zdialog_add_widget(zd,"hscale","warpval","hbw","0.0|1.0|0.01|0.0","space=5|expand");

   takeMouse(unwarpCU_mousefunc,dragcursor);                                     //  connect mouse function

   unwarpCU_started = 0;
   unwarpCU_areanumber = 0;

   zdialog_run(zd,unwarpCU_dialog_event,"save");                                 //  run dialog, parallel
   return;
}


//  dialog event and completion callback function

int unwarpCU_dialog_event(zdialog * zd, cchar *event)
{
   int      init = 0;

   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (zd->zstat)
   {
      if (zd->zstat == 1) edit_done(0);                                          //  commit edit
      else edit_cancel(0);                                                       //  discard edit
      return 1;
   }

   if (! sa_validate()) init = 1;                                                //  area invalid for curr. image file

   if (sa_stat != 3 || sa_mode == mode_image)                                    //  no select area active
      init = 1;

   if (unwarpCU_started && unwarpCU_areanumber != areanumber)                    //  select area changed
      init = 1;

   if (init) {
      unwarpCU_started = 0;
      unwarpCU_areanumber = 0;
      return 1;
   }
   
   if (strmatch(event,"focus"))                                                  //  reconnect mouse
      takeMouse(unwarpCU_mousefunc,dragcursor);

   if (strmatch(event,"warpval")) {                                              //  slider movement
      zdialog_fetch(zd,"warpval",unwarpCU_warpval);
      unwarpCU_warpfunc();
   }

   return 1;
}


//  mouse function

void unwarpCU_mousefunc()
{
   if (LMclick) {
      unwarpCU_cx = Mxclick;                                                     //  capture central point
      unwarpCU_cy = Myclick;
   }
   else if (Mxdrag || Mydrag) {
      unwarpCU_cx = Mxdrag;
      unwarpCU_cy = Mydrag;
   }
   else return;
   
   LMclick = Mxdrag = Mydrag = 0;

   if (sa_stat != 3 || sa_mode == mode_image) {                                  //  no select area active
      zmessageACK(Mwin,ZTX("no active Select Area"));
      unwarpCU_started = 0;
      Mdrag = 0;
      return;
   }
   
   unwarpCU_started = 1;                                                         //  unwarp can proceed
   unwarpCU_areanumber = areanumber;
   sa_edgecalc();                                                                //  calculate area edge distances
   
   unwarpCU_warpfunc();
   return;
}


//  warp image according to slider position

void unwarpCU_warpfunc()
{
   int         ii, px, py, ww, vstat;
   float       hsize, vsize, hpos, vpos;
   float       warpval = unwarpCU_warpval;
   float       ed, cx, cy, dx, dy;
   float       vpix[4], *pix3;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);

   if (sa_stat != 3) return;                                                     //  area erased
   
   hsize = sa_maxx - sa_minx;                                                    //  area size
   vsize = sa_maxy - sa_miny;
   
   cx = unwarpCU_cx;                                                             //  unwarp center (nose tip)
   cy = unwarpCU_cy;

   ww = E1pxm->ww;

   for (py = sa_miny; py < sa_maxy; py++)                                        //  loop all pixels in area
   for (px = sa_minx; px < sa_maxx; px++)
   {
      ii = py * ww + px;
      if (sa_pixmap[ii] < 1) continue;                                           //  pixel outside area
      
      hpos = 2.0 * (px - cx) / hsize;                                            //  horizontal pixel position, -1 ... +1
      vpos = 2.0 * (py - cy) / vsize;                                            //  vertical pixel position, -1 ... +1

      ed = sa_pixmap[ii];                                                        //  pixel edge distance

      hpos = hpos * ed / hsize;                                                  //  scale pixel position 
      vpos = vpos * ed / vsize;
      
      dx = sinf(PI * hpos);                                                      //  pixel displacement, -1 ... +1
      dy = sinf(PI * vpos);
      
      dx = dx * 0.1 * warpval * hsize;                                           //  pixel displacement, -10% ... +10%
      dy = dy * 0.1 * warpval * vsize;
      
      vstat = vpixel(E1pxm,px+dx,py+dy,vpix);                                    //  input virtual pixel
      if (vstat) {
         pix3 = PXMpix(E3pxm,px,py);                                             //  output pixel
         memcpy(pix3,vpix,pcc);
      }
   }

   Fpaint3(sa_minx,sa_miny,hsize,vsize);

   CEF->Fmods++;
   CEF->Fsaved = 0;

   return;
}


/********************************************************************************/

//  warp/distort whole image with a curved transform
//  fix perspective problems (e.g. curved walls, leaning buildings)

namespace warpC_names
{
   float       *WarpCx, *WarpCy;                                                 //  memory of all dragged pixels
   float       WarpCmem[5][100];                                                 //  undo memory, last 100 drags
   int         NWarpC;                                                           //  WarpCmem count
   int         WarpCdrag;
   int         E3ww, E3hh;

   float       $mdx, $mdy, $mdw, $mdh;                                           //  warpC_warpfunc() and warpC_wthread()
   float       $D, $span;                                                        //    use these $ args
   int         $acc;

   editfunc    EFwarpC;

   int    WarpC_dialog_event(zdialog *zd, cchar *event);
   void   WarpC_warpfunc();
   void   WarpC_mousefunc(void);
   void * WarpC_wthread(void *arg);
}

//  menu function

void m_warp_curved(GtkWidget *, cchar *)
{
   using namespace warpC_names;

   cchar  *WarpC_message = ZTX(
             " Pull an image position using the mouse. \n"
             " Make multiple mouse pulls until satisfied. \n"
             " When finished, press [done].");

   int         px, py, ii;

   F1_help_topic = "warp_curved";

   EFwarpC.menufunc = m_warp_curved;
   EFwarpC.funcname = "warp-curved";
   EFwarpC.FprevReq = 1;                                                         //  use preview
   EFwarpC.mousefunc = WarpC_mousefunc;                                          //  mouse function
   if (! edit_setup(EFwarpC)) return;                                            //  setup edit

   PXM_addalpha(E0pxm);                                                          //  15.09
   PXM_addalpha(E1pxm);
   PXM_addalpha(E3pxm);

   zdialog *zd = zdialog_new(ZTX("Warp curved"),Mwin,Bdone,Bcancel,null);
   EFwarpC.zd = zd;

   zdialog_add_widget(zd,"label","lab1","dialog",WarpC_message,"space=3");
   zdialog_add_widget(zd,"hbox","hb1","dialog",0,"space=8");
   zdialog_add_widget(zd,"button","undolast","hb1",Bundolast,"space=8");
   zdialog_add_widget(zd,"button","undoall","hb1",Bundoall,"space=2");
   zdialog_add_widget(zd,"button","grid","hb1",Bgrid,"space=15");
   zdialog_add_widget(zd,"hbox","hb2","dialog",0,"space=4");
   zdialog_add_widget(zd,"label","lab2","hb2",ZTX("warp span"),"space=8");
   zdialog_add_widget(zd,"spin","span","hb2","0.00|1.0|0.01|0.1","space=1");

   currgrid = 3;                                                                 //  use warp-curved grid

   NWarpC = WarpCdrag = 0;                                                       //  no drag data

   int cc = E3pxm->ww * E3pxm->hh * sizeof(float);
   WarpCx = (float *) zmalloc(cc);                                               //  get memory for pixel displacements
   WarpCy = (float *) zmalloc(cc);

   for (py = 0; py < E3pxm->hh; py++)                                            //  no pixel displacements
   for (px = 0; px < E3pxm->ww; px++)
   {
      ii = py * E3pxm->ww + px;
      WarpCx[ii] = WarpCy[ii] = 0.0;
   }

   E3ww = E3pxm->ww;                                                             //  preview dimensions
   E3hh = E3pxm->hh;

   zdialog_restore_inputs(zd);                                                   //  restore previous inputs
   zdialog_fetch(zd,"span",$span);                                               //  save span value

   zdialog_run(zd,WarpC_dialog_event,"save");                                    //  run dialog, parallel
   takeMouse(WarpC_mousefunc,dragcursor);                                        //  connect mouse function
   return;
}


//  dialog event and completion callback function

int warpC_names::WarpC_dialog_event(zdialog * zd, cchar *event)
{
   using namespace warpC_names;

   int         px, py, ii;
   int         fpx, fpy, epx, epy, vstat;
   float       scale, dispx, dispy;
   float       vpix[4], *pix3;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);

   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (zd->zstat) goto complete;

   if (strmatch(event,"undolast"))
   {
      if (NWarpC == 1) event = "undoall";
      else if (NWarpC) {                                                         //  undo most recent drag
         ii = --NWarpC;
         $mdx = WarpCmem[0][ii];
         $mdy = WarpCmem[1][ii];
         $mdw = -WarpCmem[2][ii];
         $mdh = -WarpCmem[3][ii];
         $span = WarpCmem[4][ii];
         zdialog_stuff(zd,"span",$span);
         $acc = 0;
         WarpC_warpfunc();                                                       //  undrag image
         $acc = 1;
         WarpC_warpfunc();                                                       //  undrag memory
      }
   }

   if (strmatch(event,"undoall"))                                                //  undo all drags
   {
      NWarpC = 0;                                                                //  erase undo memory

      for (py = 0; py < E3pxm->hh; py++)                                         //  reset pixel displacements
      for (px = 0; px < E3pxm->ww; px++)
      {
         ii = py * E3pxm->ww + px;
         WarpCx[ii] = WarpCy[ii] = 0.0;
      }
      edit_reset();                                                              //  restore image 1
   }

   if (strmatch(event,"grid")) m_gridlines(0,"grid 3");                          //  grid settings dialog

   if (strmatch(event,"span"))
      zdialog_fetch(zd,"span",$span);

   return 1;

complete:

   currgrid = 0;                                                                 //  restore normal grid settings

   if (zd->zstat)
   {
      if (NWarpC == 0) zd->zstat = 2;                                            //  no warps, cancel

      if (zd->zstat == 1)                                                        //  done
      {
         edit_fullsize();                                                        //  get full size image

         scale = 1.0 * (E3pxm->ww + E3pxm->hh) / (E3ww + E3hh);

         for (fpy = 0; fpy < E3pxm->hh; fpy++)                                   //  scale net pixel displacements
         for (fpx = 0; fpx < E3pxm->ww; fpx++)                                   //    to full image size
         {
            epx = E3ww * fpx / E3pxm->ww;
            epy = E3hh * fpy / E3pxm->hh;
            ii = epy * E3ww + epx;
            dispx = WarpCx[ii] * scale;
            dispy = WarpCy[ii] * scale;
            vstat = vpixel(E1pxm,fpx+dispx,fpy+dispy,vpix);                      //  input virtual pixel
            pix3 = PXMpix(E3pxm,fpx,fpy);                                        //  output pixel
            if (vstat) memcpy(pix3,vpix,pcc);
            else memset(pix3,0,pcc);                                             //  voided pixel
         }
         signal_thread();
         edit_done(0);                                                           //  commit edit
      }

      else edit_cancel(0);                                                       //  discard edit
      zfree(WarpCx);                                                             //  release memory
      zfree(WarpCy);
   }

   return 1;
}


//  WarpC mouse function

void warpC_names::WarpC_mousefunc(void)
{
   using namespace warpC_names;

   int      ii;

   if (Mxdrag || Mydrag)                                                         //  mouse drag underway
   {
      $mdx = Mxdown;                                                             //  drag origin
      $mdy = Mydown;
      $mdw = Mxdrag - Mxdown;                                                    //  drag increment
      $mdh = Mydrag - Mydown;
      $acc = 0;
      WarpC_warpfunc();                                                          //  drag image
      WarpCdrag = 1;
      Mxdrag = Mydrag = 0;
      return;
   }

   else if (WarpCdrag)
   {
      WarpCdrag = 0;
      $acc = 1;
      WarpC_warpfunc();                                                          //  drag done, add to memory

      if (NWarpC == 100)                                                         //  if full, throw away oldest
      {
         NWarpC = 99;
         for (ii = 0; ii < NWarpC; ii++)
         {
            WarpCmem[0][ii] = WarpCmem[0][ii+1];
            WarpCmem[1][ii] = WarpCmem[1][ii+1];
            WarpCmem[2][ii] = WarpCmem[2][ii+1];
            WarpCmem[3][ii] = WarpCmem[3][ii+1];
            WarpCmem[4][ii] = WarpCmem[4][ii+1];
         }
      }

      ii = NWarpC;
      WarpCmem[0][ii] = $mdx;                                                    //  save drag for undo
      WarpCmem[1][ii] = $mdy;
      WarpCmem[2][ii] = $mdw;
      WarpCmem[3][ii] = $mdh;
      WarpCmem[4][ii] = $span;
      NWarpC++;
   }

   return;
}


//  warp image and accumulate warp memory
//  mouse at (mx,my) is moved (mw,mh) pixels

void warpC_names::WarpC_warpfunc()
{
   using namespace warpC_names;

   float       D, d1, d2, d3, d4;

   d1 = ($mdx-0) * ($mdx-0) + ($mdy-0) * ($mdy-0);                               //  distance, mouse to 4 corners
   d2 = (E3pxm->ww-$mdx) * (E3pxm->ww-$mdx) + ($mdy-0) * ($mdy-0);
   d3 = (E3pxm->ww-$mdx) * (E3pxm->ww-$mdx) + (E3pxm->hh-$mdy) * (E3pxm->hh-$mdy);
   d4 = ($mdx-0) * ($mdx-0) + (E3pxm->hh-$mdy) * (E3pxm->hh-$mdy);

   D = d1;
   if (d2 > D) D = d2;                                                           //  find greatest corner distance
   if (d3 > D) D = d3;
   if (d4 > D) D = d4;

   $D = D * $span;

   for (int ii = 0; ii < NWT; ii++)                                              //  start worker threads
      start_wthread(WarpC_wthread,&Nval[ii]);
   wait_wthreads();                                                              //  wait for completion

   CEF->Fmods++;
   CEF->Fsaved = 0;
   Fpaint2();                                                                    //  update window
   return;
}


//  working thread to process the pixels

void * warpC_names::WarpC_wthread(void *arg)
{
   using namespace warpC_names;

   int      index = *((int *) arg);
   int      ii, px, py, vstat;
   float    d, mag, dispx, dispy;
   float    vpix[4], *pix3;
   int      nc = E1pxm->nc, pcc = nc * sizeof(float);

   for (py = index; py < E3pxm->hh; py += NWT)                                   //  process all pixels
   for (px = 0; px < E3pxm->ww; px++)
   {
      d = (px-$mdx)*(px-$mdx) + (py-$mdy)*(py-$mdy);
      mag = (1.0 - d / $D);
      if (mag < 0) continue;

      mag = mag * mag;                                                           //  faster than pow(mag,4);
      mag = mag * mag;

      dispx = -$mdw * mag;                                                       //  displacement = drag * mag
      dispy = -$mdh * mag;

      ii = py * E3pxm->ww + px;

      if ($acc) {                                                                //  drag done, accumulate drag sum
         WarpCx[ii] += dispx;
         WarpCy[ii] += dispy;
         continue;
      }

      dispx += WarpCx[ii];                                                       //  add this drag to prior sum
      dispy += WarpCy[ii];

      vstat = vpixel(E1pxm,px+dispx,py+dispy,vpix);                              //  input virtual pixel
      pix3 = PXMpix(E3pxm,px,py);                                                //  output pixel
      if (vstat) memcpy(pix3,vpix,pcc);
      else memset(pix3,0,pcc);                                                   //  voided pixel
   }

   exit_wthread();                                                               //  exit thread
   return 0;                                                                     //  not executed, avoid gcc warning
}


/********************************************************************************/

//  warp/distort whole image with a linear transform (almost)
//  fix perspective problems (e.g. leaning buildings)

namespace warpL_names
{
   float       *WarpLx, *WarpLy;                                                 //  memory of all dragged pixels
   float       WarpLmem[4][100];                                                 //  undo memory, last 100 drags
   int         NWarpL;                                                           //  WarpLmem count
   int         WarpLdrag;
   int         E3ww, E3hh;
   float       $mdx, $mdy, $mdw, $mdh;                                           //  warpL_warpfunc() and warpL_wthread()
   int         $D, $Dx, $Dy, $acc;                                               //    use these $ args

   editfunc    EFwarpL;

   int    WarpL_dialog_event(zdialog *zd, cchar *event);
   void   WarpL_mousefunc(void);
   void   WarpL_warpfunc();
   void * WarpL_wthread(void *arg);
}


//  menu function

void m_warp_linear(GtkWidget *, cchar *)
{
   using namespace warpL_names;

   cchar  *WarpL_message = ZTX(
             " Pull an image position using the mouse. \n"
             " Make multiple mouse pulls until satisfied. \n"
             " When finished, press [done].");

   int         px, py, ii;

   F1_help_topic = "warp_linear";

   EFwarpL.menufunc = m_warp_linear;
   EFwarpL.funcname = "warp-linear";
   EFwarpL.FprevReq = 1;                                                         //  use preview
   EFwarpL.mousefunc = WarpL_mousefunc;                                          //  mouse function
   if (! edit_setup(EFwarpL)) return;                                            //  setup edit

   PXM_addalpha(E0pxm);                                                          //  15.09
   PXM_addalpha(E1pxm);
   PXM_addalpha(E3pxm);

   zdialog *zd = zdialog_new(ZTX("Warp linear"),Mwin,Bdone,Bcancel,null);
   EFwarpL.zd = zd;

   zdialog_add_widget(zd,"label","lab1","dialog",WarpL_message,"space=3");
   zdialog_add_widget(zd,"hbox","hb1","dialog",0,"space=10");
   zdialog_add_widget(zd,"button","undolast","hb1",Bundolast,"space=5");
   zdialog_add_widget(zd,"button","undoall","hb1",Bundoall,"space=5");
   zdialog_add_widget(zd,"hbox","hb2","dialog",0,"space=6");
   zdialog_add_widget(zd,"button","grid","hb2",Bgrid,"space=10");

   currgrid = 4;                                                                 //  use warp-linear grid

   NWarpL = WarpLdrag = 0;                                                       //  no drag data

   int cc = E3pxm->ww * E3pxm->hh * sizeof(float);
   WarpLx = (float *) zmalloc(cc);                                               //  get memory for pixel displacements
   WarpLy = (float *) zmalloc(cc);

   for (py = 0; py < E3pxm->hh; py++)                                            //  no pixel displacements
   for (px = 0; px < E3pxm->ww; px++)
   {
      ii = py * E3pxm->ww + px;
      WarpLx[ii] = WarpLy[ii] = 0.0;
   }

   E3ww = E3pxm->ww;                                                             //  preview dimensions
   E3hh = E3pxm->hh;

   zdialog_run(zd,WarpL_dialog_event,"save");                                    //  run dialog, parallel

   takeMouse(WarpL_mousefunc,dragcursor);                                        //  connect mouse function
   return;
}


//  dialog event and completion callback function

int warpL_names::WarpL_dialog_event(zdialog * zd, cchar *event)
{
   using namespace warpL_names;

   int         px, py, ii;
   int         fpx, fpy, epx, epy, vstat;
   float       scale, dispx, dispy;
   float       vpix[4], *pix3;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);

   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (zd->zstat) goto complete;

   if (strmatch(event,"undolast"))
   {
      if (NWarpL == 1) event = "undoall";
      else if (NWarpL) {                                                         //  undo most recent drag
         ii = --NWarpL;
         $mdx = WarpLmem[0][ii];
         $mdy = WarpLmem[1][ii];
         $mdw = -WarpLmem[2][ii];
         $mdh = -WarpLmem[3][ii];
         $acc = 0;
         WarpL_warpfunc();                                                       //  undrag image
         $acc = 1;
         WarpL_warpfunc();                                                       //  undrag memory
      }
   }

   if (strmatch(event,"undoall"))                                                //  undo all drags
   {
      NWarpL = 0;                                                                //  erase undo memory

      for (py = 0; py < E3pxm->hh; py++)                                         //  reset pixel displacements
      for (px = 0; px < E3pxm->ww; px++)
      {
         ii = py * E3pxm->ww + px;
         WarpLx[ii] = WarpLy[ii] = 0.0;
      }
      edit_reset();                                                              //  restore image 1
   }

   if (strmatch(event,"grid")) m_gridlines(0,"grid 4");                          //  grid settings dialog

   return 1;

complete:

   currgrid = 0;                                                                 //  restore normal grid settings

   if (zd->zstat != 1 || NWarpL == 0)                                            //  cancel or no warps made
      edit_cancel(0);

   else
   {
      edit_fullsize();                                                           //  get full-size E1/E3

      scale = 1.0 * (E3pxm->ww + E3pxm->hh) / (E3ww + E3hh);

      for (fpy = 0; fpy < E3pxm->hh; fpy++)                                      //  scale net pixel displacements
      for (fpx = 0; fpx < E3pxm->ww; fpx++)                                      //    to full image size
      {
         epx = E3ww * fpx / E3pxm->ww;
         epy = E3hh * fpy / E3pxm->hh;
         ii = epy * E3ww + epx;
         dispx = WarpLx[ii] * scale;
         dispy = WarpLy[ii] * scale;

         vstat = vpixel(E1pxm,fpx+dispx,fpy+dispy,vpix);                         //  input virtual pixel
         pix3 = PXMpix(E3pxm,fpx,fpy);                                           //  output pixel
         if (vstat) memcpy(pix3,vpix,pcc);
         else memset(pix3,0,pcc);                                                //  voided pixel
      }

      signal_thread();
      edit_done(0);
   }

   zfree(WarpLx);                                                                //  release memory
   zfree(WarpLy);
   return 1;
}


//  WarpL mouse function

void warpL_names::WarpL_mousefunc(void)
{
   using namespace warpL_names;

   int      ii;

   if (Mxdrag || Mydrag)                                                         //  mouse drag underway
   {
      $mdx = Mxdown;                                                             //  drag origin, window coordinates
      $mdy = Mydown;
      $mdw = Mxdrag - Mxdown;                                                    //  drag increment
      $mdh = Mydrag - Mydown;
      $acc = 0;
      WarpL_warpfunc();                                                          //  drag image
      WarpLdrag = 1;
      Mxdrag = Mydrag = 0;
      return;
   }

   else if (WarpLdrag)
   {
      WarpLdrag = 0;
      $acc = 1;
      WarpL_warpfunc();                                                          //  drag done, add to memory

      if (NWarpL == 100)                                                         //  if full, throw away oldest
      {
         NWarpL = 99;
         for (ii = 0; ii < NWarpL; ii++)
         {
            WarpLmem[0][ii] = WarpLmem[0][ii+1];
            WarpLmem[1][ii] = WarpLmem[1][ii+1];
            WarpLmem[2][ii] = WarpLmem[2][ii+1];
            WarpLmem[3][ii] = WarpLmem[3][ii+1];
         }
      }

      ii = NWarpL;
      WarpLmem[0][ii] = $mdx;                                                    //  save drag for undo
      WarpLmem[1][ii] = $mdy;
      WarpLmem[2][ii] = $mdw;
      WarpLmem[3][ii] = $mdh;
      NWarpL++;
   }

   return;
}


//  warp image and accumulate warp memory
//  mouse at ($mdx,$mdy) is moved ($mdw,$mdh) pixels

void warpL_names::WarpL_warpfunc()
{
   using namespace warpL_names;

   float       d1, d2, d3, d4;

   d1 = ($mdx-0) * ($mdx-0) + ($mdy-0) * ($mdy-0);                               //  distance, mouse to 4 corners
   d2 = (E3pxm->ww-$mdx) * (E3pxm->ww-$mdx) + ($mdy-0) * ($mdy-0);
   d3 = (E3pxm->ww-$mdx) * (E3pxm->ww-$mdx) + (E3pxm->hh-$mdy) * (E3pxm->hh-$mdy);
   d4 = ($mdx-0) * ($mdx-0) + (E3pxm->hh-$mdy) * (E3pxm->hh-$mdy);

   $D = d1;
   if (d2 > $D) $D = d2;                                                         //  find greatest corner distance
   if (d3 > $D) $D = d3;
   if (d4 > $D) $D = d4;

   if ($D == d1) {                                                               //  NW corner
      $D = 1;
      $Dx = $mdx;                                                                //  x/y distance, mouse to edges
      $Dy = $mdy;
   }

   if ($D == d2) {                                                               //  NE
      $D = 2;
      $Dx = E3ww - $mdx;
      $Dy = $mdy;
   }

   if ($D == d3) {                                                               //  SE
      $D = 3;
      $Dx = E3ww - $mdx;
      $Dy = E3hh - $mdy;
   }

   if ($D == d4) {                                                               //  SW
      $D = 4;
      $Dx = $mdx;
      $Dy = E3hh - $mdy;
   }

   for (int ii = 0; ii < NWT; ii++)                                              //  start worker threads
      start_wthread(WarpL_wthread,&Nval[ii]);
   wait_wthreads();                                                              //  wait for completion

   CEF->Fmods++;
   CEF->Fsaved = 0;
   Fpaint2();                                                                    //  update window
   return;
}


//  working thread to process the pixels

void * warpL_names::WarpL_wthread(void *arg)                                     //  algorithm change
{
   using namespace warpL_names;

   int      index = *((int *) arg);
   int      ii, px, py, vstat;
   float    dx, dy, mag, dispx, dispy;
   float    vpix[4], *pix3;
   int      nc = E1pxm->nc, pcc = nc * sizeof(float);

   for (py = index; py < E3pxm->hh; py += NWT)                                   //  process all pixels
   for (px = 0; px < E3pxm->ww; px++)
   {
      if ($D == 1) {
         dx = $mdx - px;                                                         //  x/y distance, pixel to mouse
         dy = $mdy - py;
      }

      else if ($D == 2) {
         dx = px - $mdx;
         dy = $mdy - py;
      }

      else if ($D == 3) {
         dx = px - $mdx;
         dy = py - $mdy;
      }

      else /* $D == 4 */ {
         dx = $mdx - px;
         dy = py - $mdy;
      }

      mag = (1.0 - dx / $Dx) * (1.0 - dy / $Dy);                                 //  pixel movement / mouse drag

      dispx = -$mdw * mag;                                                       //  displacement = drag * mag
      dispy = -$mdh * mag;

      ii = py * E3pxm->ww + px;

      if ($acc) {                                                                //  drag done, accumulate drag sum
         WarpLx[ii] += dispx;
         WarpLy[ii] += dispy;
         continue;
      }

      dispx += WarpLx[ii];                                                       //  add this drag to prior sum
      dispy += WarpLy[ii];

      vstat = vpixel(E1pxm,px+dispx,py+dispy,vpix);                              //  input virtual pixel
      pix3 = PXMpix(E3pxm,px,py);                                                //  output pixel
      if (vstat) memcpy(pix3,vpix,pcc);
      else memset(pix3,0,pcc);                                                   //  voided pixel
   }

   exit_wthread();                                                               //  exit thread
   return 0;                                                                     //  not executed, avoid gcc warning
}


/********************************************************************************/

//  warp/distort whole image using affine transform
//  (straight lines remain straight)

float       WarpF_old[3][2];                                                     //  3 original image points
float       WarpF_new[3][2];                                                     //  corresponding warped points
float       WarpF_coeff[6];                                                      //  transform coefficients
float       WarpF_Icoeff[6];                                                     //  inverse transform coefficients
int         WarpF_ftf;                                                           //  first time flag

editfunc    EFwarpF;

void  WarpF_warpfunc();                                                          //  image warp function
void  WarpF_mousefunc(void);
void  WarpF_affine(float po[3][2], float pn[3][2], float coeff[6]);              //  compute affine transform coefficients
void  WarpF_invert(float coeff[6], float Icoeff[6]);                             //  compute reverse transform coefficients


//  menu function

void m_warp_affine(GtkWidget *, cchar *)
{
   int      WarpF_dialog_event(zdialog *zd, cchar *event);

   cchar  *WarpF_message = ZTX(
             " Pull on an image corner using the mouse. \n"
             " Make multiple mouse pulls until satisfied. \n"
             " When finished, press [done].");

   F1_help_topic = "warp_affine";

   EFwarpF.menufunc = m_warp_affine;
   EFwarpF.funcname = "warp-affine";
   EFwarpF.FprevReq = 1;                                                         //  use preview
   EFwarpF.mousefunc = WarpF_mousefunc;                                          //  mouse function
   if (! edit_setup(EFwarpF)) return;                                            //  setup edit

   PXM_addalpha(E0pxm);                                                          //  15.09
   PXM_addalpha(E1pxm);
   PXM_addalpha(E3pxm);

   zdialog *zd = zdialog_new(ZTX("Warp affine"),Mwin,Bdone,Bcancel,null);
   EFwarpF.zd = zd;

   zdialog_add_widget(zd,"label","lab1","dialog",WarpF_message,"space=3");
   zdialog_add_widget(zd,"hbox","hb2","dialog",0,"space=6");
   zdialog_add_widget(zd,"button","grid","hb2",Bgrid,"space=10");

   currgrid = 5;                                                                 //  use warp-affine grid
   WarpF_ftf = 1;                                                                //  1st warp flag

   zdialog_run(zd,WarpF_dialog_event,"save");                                    //  run dialog, parallel

   takeMouse(WarpF_mousefunc,dragcursor);                                        //  connect mouse function
   return;
}


//  dialog event and completion callback function

int WarpF_dialog_event(zdialog *zd, cchar *event)
{
   float       scale;
   int         ww, hh;

   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (strmatch(event,"grid")) m_gridlines(0,"grid 5");                          //  grid settings dialog

   if (! zd->zstat) return 1;                                                    //  wait for completion

   currgrid = 0;                                                                 //  restore normal grid settings

   if (zd->zstat != 1 || ! CEF->Fmods) {
      edit_cancel(0); 
      return 1;
   }

   ww = E3pxm->ww;                                                               //  preview image dimensions
   hh = E3pxm->hh;

   edit_fullsize();                                                              //  get full-size images

   scale = 1.0 * (E3pxm->ww + E3pxm->hh) / (ww + hh);                            //  preview to full-size scale factor

   WarpF_old[0][0] = WarpF_old[0][0] * scale;                                    //  re-scale new and old points
   WarpF_old[0][1] = WarpF_old[0][1] * scale;
   WarpF_old[1][0] = WarpF_old[1][0] * scale;
   WarpF_old[1][1] = WarpF_old[1][1] * scale;
   WarpF_old[2][0] = WarpF_old[2][0] * scale;
   WarpF_old[2][1] = WarpF_old[2][1] * scale;

   WarpF_new[0][0] = WarpF_new[0][0] * scale;
   WarpF_new[0][1] = WarpF_new[0][1] * scale;
   WarpF_new[1][0] = WarpF_new[1][0] * scale;
   WarpF_new[1][1] = WarpF_new[1][1] * scale;
   WarpF_new[2][0] = WarpF_new[2][0] * scale;
   WarpF_new[2][1] = WarpF_new[2][1] * scale;

   WarpF_warpfunc();                                                             //  warp full-size image

   edit_done(0);
   return 1;
}


//  WarpF mouse function

void  WarpF_mousefunc(void)
{
   int      mdx1, mdy1, mdx2, mdy2;
   float    x1o, y1o, x2o, y2o, x3o, y3o;
   float    x1n, y1n, x2n, y2n, x3n, y3n;
   float    a, b, c, d, e, f;

   if (! Mxdrag && ! Mydrag) return;

   mdx1 = Mxdown;                                                                //  mouse drag origin
   mdy1 = Mydown;
   mdx2 = Mxdrag;                                                                //  mouse drag position
   mdy2 = Mydrag;

   Mxdown = Mxdrag;                                                              //  reset origin for next time
   Mydown = Mydrag;

   x1n = mdx1;                                                                   //  point 1 = drag origin
   y1n = mdy1;
   x2n = E3pxm->ww - x1n;                                                        //  point 2 = mirror of point1
   y2n = E3pxm->hh - y1n;
   x3n = E3pxm->ww * (y2n / E3pxm->hh);
   y3n = E3pxm->hh * (1.0 - (x2n / E3pxm->ww));

   if (WarpF_ftf)                                                                //  first warp
   {
      WarpF_ftf = 0;
      x1o = x1n;                                                                 //  old = current positions
      y1o = y1n;
      x2o = x2n;
      y2o = y2n;
      x3o = x3n;
      y3o = y3n;
   }
   else
   {
      WarpF_invert(WarpF_coeff,WarpF_Icoeff);                                    //  get inverse coefficients
      a = WarpF_Icoeff[0];
      b = WarpF_Icoeff[1];
      c = WarpF_Icoeff[2];
      d = WarpF_Icoeff[3];
      e = WarpF_Icoeff[4];
      f = WarpF_Icoeff[5];

      x1o = a * x1n + b * y1n + c;                                               //  compute old from current positions
      y1o = d * x1n + e * y1n + f;
      x2o = a * x2n + b * y2n + c;
      y2o = d * x2n + e * y2n + f;
      x3o = a * x3n + b * y3n + c;
      y3o = d * x3n + e * y3n + f;
   }

   WarpF_old[0][0] = x1o;                                                        //  set up 3 old points and corresponding
   WarpF_old[0][1] = y1o;                                                        //    new points for affine translation
   WarpF_old[1][0] = x2o;
   WarpF_old[1][1] = y2o;
   WarpF_old[2][0] = x3o;
   WarpF_old[2][1] = y3o;

   x1n = mdx2;                                                                   //  point 1 new position = drag position
   y1n = mdy2;
   x2n = E3pxm->ww - x1n;                                                        //  point 2 new = mirror of point1 new
   y2n = E3pxm->hh - y1n;

   WarpF_new[0][0] = x1n;                                                        //  3 new points
   WarpF_new[0][1] = y1n;
   WarpF_new[1][0] = x2n;
   WarpF_new[1][1] = y2n;
   WarpF_new[2][0] = x3n;
   WarpF_new[2][1] = y3n;

   WarpF_warpfunc();                                                             //  do the warp

   Mxdrag = Mydrag = 0;
   return;
}


//  warp image and accumulate warp memory

void  WarpF_warpfunc()
{
   void * WarpF_wthread(void *);

   WarpF_affine(WarpF_old, WarpF_new, WarpF_coeff);                              //  get coefficients for forward transform
   WarpF_invert(WarpF_coeff, WarpF_Icoeff);                                      //  get coefficients for reverse transform

   for (int ii = 0; ii < NWT; ii++)                                              //  start worker threads
      start_wthread(WarpF_wthread,&Nval[ii]);
   wait_wthreads();                                                              //  wait for completion

   CEF->Fmods++;
   CEF->Fsaved = 0;
   Fpaint2();                                                                    //  update window
   return;
}


//  working thread to process the pixels

void * WarpF_wthread(void *arg)
{
   int      index = *((int *) arg);
   float    a, b, c, d, e, f;
   int      px3, py3, vstat;
   float    px1, py1;
   float    vpix1[4], *pix3;
   int      nc = E1pxm->nc, pcc = nc * sizeof(float);

   a = WarpF_Icoeff[0];                                                          //  coefficients to map output pixels
   b = WarpF_Icoeff[1];                                                          //    to corresponding input pixels
   c = WarpF_Icoeff[2];
   d = WarpF_Icoeff[3];
   e = WarpF_Icoeff[4];
   f = WarpF_Icoeff[5];

   for (py3 = index; py3 < E3pxm->hh; py3 += NWT)                                //  process all pixels
   for (px3 = 0; px3 < E3pxm->ww; px3++)
   {
      px1 = a * px3 + b * py3 + c;                                               //  corresponding input pixel
      py1 = d * px3 + e * py3 + f;

      vstat = vpixel(E1pxm,px1,py1,vpix1);                                       //  input virtual pixel
      pix3 = PXMpix(E3pxm,px3,py3);                                              //  output pixel
      if (vstat) memcpy(pix3,vpix1,pcc);
      else memset(pix3,0,pcc);                                                   //  voided pixel
   }

   exit_wthread();                                                               //  exit thread
   return 0;                                                                     //  not executed, avoid gcc warning
}


/********************************************************************************

   Compute affine transformation of an image (warp image).

   Given 3 new (warped) positions for 3 image points, derive the
   coefficients of the translation function to warp the entire image.

   Inputs:
      pold[3][2]  (x,y) coordinates for 3 points in original image
      pnew[3][2]  (x,y) coordinates for same points in warped image

   Output:
      coeff[6]  coefficients of translation function which can be used
                to convert all image points to their warped positions

   If coeff[6] = (a, b, c, d, e, f) then the following formula
   can be used to convert an image point to its warped position:

      Xnew = a * Xold + b * Yold + c
      Ynew = d * Xold + e * Yold + f

*********************************************************************************/

void WarpF_affine(float pold[3][2], float pnew[3][2], float coeff[6])
{
   float    x11, y11, x12, y12, x13, y13;                                        //  original points
   float    x21, y21, x22, y22, x23, y23;                                        //  moved points
   float    a, b, c, d, e, f;                                                    //  coefficients
   float    A1, A2, B1, B2, C1, C2;

   x11 = pold[0][0];
   y11 = pold[0][1];
   x12 = pold[1][0];
   y12 = pold[1][1];
   x13 = pold[2][0];
   y13 = pold[2][1];

   x21 = pnew[0][0];
   y21 = pnew[0][1];
   x22 = pnew[1][0];
   y22 = pnew[1][1];
   x23 = pnew[2][0];
   y23 = pnew[2][1];

   A1 = x11 - x12;
   A2 = x12 - x13;
   B1 = y11 - y12;
   B2 = y12 - y13;
   C1 = x21 - x22;
   C2 = x22 - x23;

   a = (B1 * C2 - B2 * C1) / (A2 * B1 - A1 * B2);
   b = (A1 * C2 - A2 * C1) / (A1 * B2 - A2 * B1);
   c = x23 - a * x13 - b * y13;

   C1 = y21 - y22;
   C2 = y22 - y23;

   d = (B1 * C2 - B2 * C1) / (A2 * B1 - A1 * B2);
   e = (A1 * C2 - A2 * C1) / (A1 * B2 - A2 * B1);
   f = y23 - d * x13 - e * y13;

   coeff[0] = a;
   coeff[1] = b;
   coeff[2] = c;
   coeff[3] = d;
   coeff[4] = e;
   coeff[5] = f;

   return;
}


/********************************************************************************

   Invert affine transform

   Input:
      coeff[6]  coefficients of translation function to convert
                image points to their warped positions
   Output:
      Icoeff[6]  coefficients of translation function to convert
                 warped image points to their original positions

   If Icoeff[6] = (a, b, c, d, e, f) then the following formula can be
      used to translate a warped image point to its original position:

      Xold = a * Xnew + b * Ynew + c
      Yold = d * Xnew + e * Ynew + f

*********************************************************************************/

void WarpF_invert(float coeff[6], float Icoeff[6])
{
   float    a, b, c, d, e, f, Z;

   a = coeff[0];
   b = coeff[1];
   c = coeff[2];
   d = coeff[3];
   e = coeff[4];
   f = coeff[5];

   Z = 1.0 / (a * e - b * d);

   Icoeff[0] = e * Z;
   Icoeff[1] = - b * Z;
   Icoeff[2] = Z * (b * f - c * e);
   Icoeff[3] = - d * Z;
   Icoeff[4] = a * Z;
   Icoeff[5] = Z * (c * d - a * f);

   return;
}


/********************************************************************************/

//  Flatten a photographed book page.
//  Compensate for page curvature at the center binding.

editfunc    EFflatbook;                                                          //  edit function data

namespace flatbook {
   int      Tmx[20], Tmy[20], Bmx[20], Bmy[20];                                  //  top/bottom mouse click points
   int      Tnm = 0, Bnm = 0;                                                    //  top/bottom click points counts
   int      Tbase, Bbase;                                                        //  top/bottom points, low values
   int      ww, hh;                                                              //  image dimensions
   int      E3warped;                                                            //  flag, E3 image is warped
   double   *Tfy, *Bfy;                                                          //  derived top/bottom y-shifts [ww]
   double   *Tfx, *Bfx;                                                          //  derived top/bottom x-shifts [ww]
   double   Tstretch, Bstretch;                                                  //  top/bottom x-shift stretch factors
}


//  menu function

void m_flatbook(GtkWidget *, const char *)
{
   using namespace flatbook;

   int    flatbook_dialog_event(zdialog* zd, const char *event);
   void * flatbook_thread(void *);
   void   flatbook_mousefunc();
   void   flatbook_draw();

   cchar    *title = ZTX("Flatten Book Page Photo");
   cchar    *guide = ZTX("Trim image to isolate one page. \n"
                         "Map top and bottom edges with \n"
                         "4+ mouse clicks, then flatten: ");
   cchar    *stretch = ZTX("Stretch curved-down surfaces:");
   zdialog  *zd;

   F1_help_topic = "flatten_book";

   EFflatbook.menufunc = m_flatbook;
   EFflatbook.funcname = "flatbook";                                             //  func name, no preview, no area
   EFflatbook.threadfunc = flatbook_thread;
   EFflatbook.mousefunc = flatbook_mousefunc;
   if (! edit_setup(EFflatbook)) return;                                         //  setup edit

   int cc = E1pxm->ww * sizeof(double);                                          //  allocate memory
   Tfy = (double *) zmalloc(cc);
   Bfy = (double *) zmalloc(cc);
   Tfx = (double *) zmalloc(cc);
   Bfx = (double *) zmalloc(cc);

   if (ww != E1pxm->ww) Tnm = Bnm = 0;                                           //  clear prior points if image
   if (hh != E1pxm->hh) Tnm = Bnm = 0;                                           //    size is different

   ww = E1pxm->ww;                                                               //  set image size
   hh = E1pxm->hh;
   E3warped = 0;                                                                 //  no E3 warp yet

/***
           _________________________________
          |    Flatten Book Page Photo      |
          |                                 |
          | Trim image to isolate one page. |
          | Map top and bottom edges with   |
          | 4+ mouse clicks, then flatten:  |
          |   [clear]  [flatten]  [undo]    |
          |                                 |
          | Stretch curved-down surfaces:   |
          |  Top: ========[]==============  |
          |  Bottom: =========[]==========  |
          |                                 |
          |              [done] [cancel]    |
          |_________________________________|

***/

   zd = zdialog_new(title,Mwin,Bdone,Bcancel,null);                              //  flatbook dialog
   EFflatbook.zd = zd;

   zdialog_add_widget(zd,"hbox","hbg","dialog",0,"space=3");
   zdialog_add_widget(zd,"label","labg","hbg",guide,"space=8");
   zdialog_add_widget(zd,"hbox","hbf","dialog");
   zdialog_add_widget(zd,"button","clear","hbf",Bclear,"space=10");
   zdialog_add_widget(zd,"button","flatten","hbf",Bflatten,"space=10");
   zdialog_add_widget(zd,"button","undo","hbf",Bundo,"space=10");
   zdialog_add_widget(zd,"hbox","space","dialog",0,"space=5");
   zdialog_add_widget(zd,"hbox","hbs1","dialog");
   zdialog_add_widget(zd,"label","labs1","hbs1",stretch,"space=8");
   zdialog_add_widget(zd,"hbox","hbs2","dialog");
   zdialog_add_widget(zd,"label","labs2","hbs2",ZTX("Top:"),"space=8");
   zdialog_add_widget(zd,"hscale","top","hbs2","1|30|0.01|1","expand|space=5");
   zdialog_add_widget(zd,"hbox","hbs3","dialog");
   zdialog_add_widget(zd,"label","labs3","hbs3",ZTX("Bottom:"),"space=8");
   zdialog_add_widget(zd,"hscale","bottom","hbs3","1|30|0.01|1","expand|space=5");

   zdialog_restore_inputs(zd);                                                   //  restore prior inputs

   zdialog_resize(zd,300,0);
   zdialog_run(zd,flatbook_dialog_event,"save");                                 //  run dialog - parallel

   takeMouse(flatbook_mousefunc,dragcursor);                                     //  connect mouse function

   if (Tnm + Bnm > 0) flatbook_draw();
   return;
}


//  flatbook dialog event and completion function

int flatbook_dialog_event(zdialog *zd, const char *event)                        //  flatbook dialog event function
{
   using namespace flatbook;

   void flatbook_draw();

   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (zd->zstat) {
      zfree(Tfy);                                                                //  free memory
      zfree(Bfy);
      zfree(Tfx);
      zfree(Bfx);
      if (zd->zstat == 1) edit_done(0);                                          //  commit edit
      else edit_cancel(0);                                                       //  cancel or destroy
      return 1;
   }

   if (strmatch(event,"clear")) {                                                //  clear all mouse points
      Tnm = Bnm = 0;
      edit_undo();
      flatbook_draw();
   }

   if (strstr("flatten top bottom",event)) {
      if (Tnm < 4 || Bnm < 4) return 1;
      zdialog_fetch(zd,"top",Tstretch);
      zdialog_fetch(zd,"bottom",Bstretch);
      signal_thread();                                                           //  trigger update thread
   }

   if (strmatch(event,"undo")) {
      edit_undo();
      E3warped = 0;
      flatbook_draw();
   }

   if (strmatch(event,"line_color")) flatbook_draw();                            //  refresh lines after window redraw

   return 1;
}


//  flatbook mouse function

void flatbook_mousefunc()
{
   using namespace flatbook;

   void flatbook_draw();

   int      ii, jj, mx, my, dist;
   int      Fadd = 0, close;

   if (! (LMclick || RMclick || Mxdrag || Mydrag)) return;                       //  ignore mouse movement

   if (E3warped) return;                                                         //  E3 image warped, ignore mouse

   if (LMclick || RMclick) {                                                     //  left or right mouse click
      mx = Mxclick;
      my = Myclick;
      if (LMclick) Fadd = 1;
      LMclick = RMclick = 0;
   }
   else {                                                                        //  mouse drag
      mx = Mxdrag;
      my = Mydrag;
      Mxdrag = Mydrag = 0;
      Fadd = 1;
   }

   close = 0.01 * ww;                                                            //  1% of image size

   if (my < hh/2) {
      for (ii = 0; ii < Tnm; ii++) {                                             //  compare mouse position with
         dist = abs(mx-Tmx[ii]);                                                 //    existing top points
         if (dist < close) {
            Tnm--;                                                               //  delete any too close
            for (jj = ii--; jj < Tnm; jj++) {
               Tmx[jj] = Tmx[jj+1];
               Tmy[jj] = Tmy[jj+1];
            }
         }
      }
   }
   else {
      for (ii = 0; ii < Bnm; ii++) {                                             //  same for bottom points
         dist = abs(mx-Bmx[ii]);
         if (dist < close) {
            Bnm--;
            for (jj = ii--; jj < Bnm; jj++) {
               Bmx[jj] = Bmx[jj+1];
               Bmy[jj] = Bmy[jj+1];
            }
         }
      }
   }

   if (Fadd)                                                                     //  add new mouse position
   {
      if (my < hh/2) {                                                           //  add to top points
         if (Tnm == 20) return;
         for (ii = 0; ii < Tnm; ii++)
            if (mx < Tmx[ii]) break;
         for (jj = Tnm; jj > ii; jj--) {
            Tmx[jj] = Tmx[jj-1];
            Tmy[jj] = Tmy[jj-1];
         }
         Tmx[ii] = mx;
         Tmy[ii] = my;
         Tnm++;
      }
      else {                                                                     //  add to bottom points
         if (Bnm == 20) return;
         for (ii = 0; ii < Bnm; ii++)
            if (mx < Bmx[ii]) break;
         for (jj = Bnm; jj > ii; jj--) {
            Bmx[jj] = Bmx[jj-1];
            Bmy[jj] = Bmy[jj-1];
         }
         Bmx[ii] = mx;
         Bmy[ii] = my;
         Bnm++;
      }
   }

   edit_undo();                                                                  //  erase and redraw
   flatbook_draw();

   return;
}


//  generate spline curves from mouse points and draw on the image

void flatbook_draw()
{
   using namespace flatbook;

   float    Tnx[20], Tny[20], Bnx[20], Bny[20];                                  //  top/bottom spline nodes
   float    *ppix3;
   int      topmax, topmin, bottmax, bottmin;                                    //  limits of image update areas
   int      ii, px, py, qx, qy, npq;

   topmax = topmin = bottmax = bottmin = 0;
   npq = 3;                                                                      //  pixel block for drawing nodes

   if (Tnm > 0) topmax = topmin = Tmy[0];

   for (ii = 0; ii < Tnm; ii++)                                                  //  draw top mouse points
   {
      px = Tmx[ii];
      py = Tmy[ii];

      for (qy = py-npq; qy <= py+npq; qy++)
      for (qx = px-npq; qx <= px+npq; qx++) {
         if (qx < 0 || qx > ww-1) continue;
         if (qy < 0 || qy > hh-1) continue;
         ppix3 = PXMpix(E3pxm,qx,qy);
         ppix3[0] = LINE_COLOR[0];
         ppix3[1] = LINE_COLOR[1];
         ppix3[2] = LINE_COLOR[2];
         if (qy > topmax) topmax = qy;
         if (qy < topmin) topmin = qy;
      }
   }

   if (Bnm > 0) bottmax = bottmin = Bmy[0];

   for (ii = 0; ii < Bnm; ii++)                                                  //  draw bottom mouse points
   {
      px = Bmx[ii];
      py = Bmy[ii];

      for (qy = py-npq; qy <= py+npq; qy++)
      for (qx = px-npq; qx <= px+npq; qx++) {
         if (qx < 0 || qx > ww-1) continue;
         if (qy < 0 || qy > hh-1) continue;
         ppix3 = PXMpix(E3pxm,qx,qy);
         ppix3[0] = LINE_COLOR[0];
         ppix3[1] = LINE_COLOR[1];
         ppix3[2] = LINE_COLOR[2];
         if (qy > bottmax) bottmax = qy;
         if (qy < bottmin) bottmin = qy;
      }
   }

   if (Tnm > 3)                                                                  //  top mouse points
   {
      Tbase = Tmy[0];                                                            //  find lowest top point
      for (ii = 0; ii < Tnm; ii++)
         if (Tmy[ii] < Tbase) Tbase = Tmy[ii];

      for (ii = 0; ii < Tnm; ii++) {                                             //  copy top points to spline nodes
         Tnx[ii] = Tmx[ii];                                                      //    and convert to 0 base
         Tny[ii] = Tmy[ii] - Tbase;
      }

      spline1(Tnm,Tnx,Tny);                                                      //  generate spline curve

      for (px = 0; px < ww; px++)                                                //  generate top curve y shifts
         Tfy[px] = spline2(px);

      for (px = 0; px < ww; px++) {                                              //  draw top curve over image
         py = Tfy[px] + Tbase;
         if (py < 0 || py > hh-1) continue;
         ppix3 = PXMpix(E3pxm,px,py);
         ppix3[0] = LINE_COLOR[0];
         ppix3[1] = LINE_COLOR[1];
         ppix3[2] = LINE_COLOR[2];
         if (py < topmin) topmin = py;
         if (py > topmax) topmax = py;
      }
   }

   if (Bnm > 3)                                                                  //  bottom mouse points
   {
      Bbase = Bmy[0];                                                            //  find lowest bottom point
      for (ii = 0; ii < Bnm; ii++)
         if (Bmy[ii] < Bbase) Bbase = Bmy[ii];

      for (ii = 0; ii < Bnm; ii++) {                                             //  copy bottom points to spline nodes
         Bnx[ii] = Bmx[ii];                                                      //    and convert to 0 base
         Bny[ii] = Bmy[ii] - Bbase;
      }

      spline1(Bnm,Bnx,Bny);                                                      //  generate spline curve

      for (px = 0; px < ww; px++)                                                //  generate bottom curve y shifts
         Bfy[px] = spline2(px);

      for (px = 0; px < ww; px++) {                                              //  draw bottom curve over image
         py = Bfy[px] + Bbase;
         if (py < 0 || py > hh-1) continue;                                      //  15.09
         ppix3 = PXMpix(E3pxm,px,py);
         ppix3[0] = LINE_COLOR[0];
         ppix3[1] = LINE_COLOR[1];
         ppix3[2] = LINE_COLOR[2];
         if (py < bottmin) bottmin = py;
         if (py > bottmax) bottmax = py;
      }
   }

   if (Tnm) Fpaint3(0,topmin,ww,topmax-topmin);                                  //  paint modified areas
   if (Bnm) Fpaint3(0,bottmin,ww,bottmax-bottmin);

   if (Tnm + Bnm) CEF->Fmods++;                                                  //  necessary
   CEF->Fsaved = 0;

   return;
}


//  flatbook thread function - warp image based on input parameters

void * flatbook_thread(void *)
{
   using namespace flatbook;

   void  * flatbook_wthread(void *arg);                                          //  worker thread

   int      px;
   double   Tslope, Bslope, slope, Rt, Rb;

   while (true)
   {
      thread_idle_loop();

      E3warped = 1;                                                              //  flag, E3 image is warped

      Tfx[0] = Bfx[0] = 0;

      for (px = 1; px < ww; px++)                                                //  loop page width
      {
         Tslope = fabs(Tfy[px] - Tfy[px-1]);                                     //  slope of top function
         Bslope = fabs(Bfy[px] - Bfy[px-1]);                                     //  slope of bottom function

         slope = Tslope * Tstretch;                                              //  mean slope, stretched              15.09
         slope = cos(atan(slope));
         Tfx[px] = Tfx[px-1] + slope;                                            //  resulting pixel x-shift

         slope = Bslope * Bstretch;
         slope = cos(atan(slope));
         Bfx[px] = Bfx[px-1] + slope;
      }

      Rt = (ww-1) / Tfx[ww-1];                                                   //  make overall width the same
      Rb = (ww-1) / Bfx[ww-1];
      for (px = 0; px < ww; px++) {
         Tfx[px] = Rt * Tfx[px];                                                 //  (using separate loops causes
         Bfx[px] = Rb * Bfx[px];                                                 //      GCC optimization bug)
      }

      for (int ii = 0; ii < NWT; ii++)                                           //  start worker threads
         start_wthread(flatbook_wthread,&Nval[ii]);
      wait_wthreads();                                                           //  wait for completion

      CEF->Fmods++;                                                              //  image is modified
      CEF->Fsaved = 0;
      Fpaint2();                                                                 //  update window
   }

   return 0;
}


void * flatbook_wthread(void *arg)
{
   using namespace flatbook;

   int      index = *((int *) (arg));
   int      px, py, vstat;
   double   Rt, Rb;
   float    *pix3, vpix[4], vpx, vpy;
   int      nc = E1pxm->nc, pcc = nc * sizeof(float);

   for (py = index; py < hh; py += NWT)                                          //  loop from top to bottom
   {
      Rt = 1.0 - 1.0 * py / hh;                                                  //  Rt from 1 to 0  - top weight
      Rb = 1.0 - Rt;                                                             //  Rb from 0 to 1  - bottom weight

      for (px = 0; px < ww; px++)                                                //  loop page width
      {
         pix3 = PXMpix(E3pxm,px,py);                                             //  output pixel

         vpx = Rt * Tfx[px] + Rb * Bfx[px];                                      //  source pixel for (px,py)
         vpy = py + Rt * Tfy[int(vpx)] + Rb * Bfy[int(vpx)];

         vstat = vpixel(E1pxm,vpx,vpy,vpix);                                     //  get source pixel
         if (vstat) memcpy(pix3,vpix,pcc);
         else memset(pix3,0,pcc);                                                //  off page, void pixel
      }
   }

   exit_wthread();
   return 0;
}


/********************************************************************************/

//  project image on to a sphere with adjustable radius (flatness)

namespace sphere_names
{
   int         E3ww, E3hh;                                                       //  image dimensions
   int         Xcen, Ycen, Dia;                                                  //  center and diameter of sphere
   float       flatten = 0;                                                      //  flatten parameter
   float       magnify = 1;                                                      //  magnify parameter
   editfunc    EFsphere;                                                         //  edit function data
}


//  menu function

void m_sphere(GtkWidget *, const char *)                                         //  15.11
{
   using namespace sphere_names;

   int    sphere_dialog_event(zdialog* zd, const char *event);
   void   sphere_mousefunc(void);
   void * sphere_thread(void *);
   
   cchar  *title = ZTX("Spherical Projection");

   F1_help_topic = "sphere";

   EFsphere.menufunc = m_sphere;
   EFsphere.funcname = "sphere";                                                 //  function name
   EFsphere.FprevReq = 1;                                                        //  use preview edit mode
   EFsphere.mousefunc = sphere_mousefunc;                                        //  mouse function
   EFsphere.threadfunc = sphere_thread;                                          //  thread function
   if (! edit_setup(EFsphere)) return;                                           //  setup edit

   PXM_addalpha(E0pxm);
   PXM_addalpha(E1pxm);
   PXM_addalpha(E3pxm);

   E3ww = E3pxm->ww;
   E3hh = E3pxm->hh;
   
   Xcen = E3ww / 2;                                                              //  set dialog defaults
   Ycen = E3hh / 2;
   Dia = E3ww;
   if (E3hh < E3ww) Dia = E3hh;
   flatten = 0;
   magnify = 1;

/***
       _____________________________
      |    Spherical Projection     |
      |                             |
      | Flatten ======[]==========  |
      | Magnify ==========[]======  |
      |                             |
      |            [done] [cancel]  |
      |_____________________________|

***/

   zdialog *zd = zdialog_new(title,Mwin,Bdone,Bcancel,null);                     //  sphere dialog
   CEF->zd = zd;

   zdialog_add_widget(zd,"hbox","hbflat","dialog");
   zdialog_add_widget(zd,"label","labflat","hbflat",Bflatten,"space=5");
   zdialog_add_widget(zd,"hscale","flatten","hbflat","0.0|0.999|0.001|0.0","expand");
   zdialog_add_widget(zd,"hbox","hbmag","dialog");
   zdialog_add_widget(zd,"label","labmag","hbmag",ZTX("Magnify"),"space=5");
   zdialog_add_widget(zd,"hscale","magnify","hbmag","1.0|2.0|0.001|1.0","expand");

   zdialog_resize(zd,250,0);
   zdialog_run(zd,sphere_dialog_event,"save");                                   //  run dialog - parallel
   takeMouse(sphere_mousefunc,0);                                                //  connect mouse function
   signal_thread();                                                              //  trigger update thread
   return;
}


//  sphere dialog event and completion function

int sphere_dialog_event(zdialog *zd, const char *event)
{
   using namespace sphere_names;
   
   float    scale;

   if (strmatch(event,"done")) zd->zstat = 1;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 1;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 2;                                  //  from f_open()

   if (zd->zstat)
   {
      if (zd->zstat == 1)                                                        //  done
      {
         edit_fullsize();                                                        //  get full size image
         scale = 1.0 * (E3pxm->ww + E3pxm->hh) / (E3ww + E3hh);                  //  scale up the parameters
         E3ww = E3pxm->ww;
         E3hh = E3pxm->hh;
         Xcen = scale * Xcen;
         Ycen = scale * Ycen;
         Dia = scale * Dia;
         signal_thread();                                                        //  recalculate image
         edit_done(0);                                                           //  commit edit
      }
      else edit_cancel(0);                                                       //  discard edit
      return 1;
   }
   
   if (strmatch(event,"flatten")) {
      zdialog_fetch(zd,"flatten",flatten);
      signal_thread();                                                           //  trigger update thread
   }

   if (strmatch(event,"magnify")) {
      zdialog_fetch(zd,"magnify",magnify);
      signal_thread();                                                           //  trigger update thread
   }

   return 1;
}


//  mouse function - get new center from mouse click or drag

void sphere_mousefunc(void)
{
   using namespace sphere_names;
   
   int      edist;
   
   if (! LMclick && ! Mdrag) return;
   LMclick = 0;
   
   if (Mxposn < 0.1 * E3ww || Mxposn > 0.9 * E3ww) return;                       //  ignore if near image edge
   if (Myposn < 0.1 * E3hh || Myposn > 0.9 * E3hh) return;

   Xcen = Mxposn;                                                                //  new center
   Ycen = Myposn;
   
   edist = Xcen;                                                                 //  find nearest edge distance
   if (Ycen < edist) edist = Ycen;
   if (E3ww - Xcen < edist) edist = E3ww - Xcen;
   if (E3hh - Ycen < edist) edist = E3hh - Ycen;
   
   Dia = 2 * edist;                                                              //  new sphere diameter
   signal_thread();

   return;
}


//  thread function - multiple working threads to update image

void * sphere_thread(void *)
{
   void  * sphere_wthread(void *arg);                                            //  worker thread

   while (true)
   {
      thread_idle_loop();                                                        //  wait for work or exit request

      for (int ii = 0; ii < NWT; ii++)                                           //  start worker threads
         start_wthread(sphere_wthread,&Nval[ii]);
      wait_wthreads();                                                           //  wait for completion

      CEF->Fmods++;                                                              //  image modified
      CEF->Fsaved = 0;                                                           //  not saved

      Fpaint2();                                                                 //  update window
   }

   return 0;                                                                     //  not executed, stop warning
}


void * sphere_wthread(void *arg)                                                 //  worker thread function
{
   using namespace sphere_names;

   int         index = *((int *) (arg));
   int         px3, py3, dx, dy, vstat;
   float       px1, py1, *pix3, vpix[4];
   float       s1, s2, D, T;
   int         nc = E1pxm->nc, pcc = nc * sizeof(float);
   
   D = Dia / (1.0 - flatten);
   
   for (py3 = index; py3 < E3hh; py3 += NWT)                                     //  loop output pixels
   for (px3 = 0; px3 < E3ww; px3++)
   {
      pix3 = PXMpix(E3pxm,px3,py3);                                              //  output pixel
      dx = px3 - Xcen;
      dy = py3 - Ycen;
      s1 = sqrtf(dx*dx + dy*dy);                                                 //  dist. from center to output pixel
      T = s1 * PI / D / magnify;                                                 //  sine of subtended angle
      if (T > 1.0) {
         memset(pix3,0,pcc);                                                     //  pixel = black, alpha = 0
         continue;
      }
      s2 = D / PI * asinf(T);                                                    //  corresp. dist. on sphere
      px1 = Xcen + dx * s2 / s1;                                                 //  input v.pixel
      py1 = Ycen + dy * s2 / s1;
      vstat = vpixel(E1pxm,px1,py1,vpix);
      if (vstat) memcpy(pix3,vpix,pcc);
      else memset(pix3,0,pcc);                                                   //  pixel = black, alpha = 0
   }

   exit_wthread();                                                               //  exit thread
   return 0;                                                                     //  not executed, avoid gcc warning
}


/********************************************************************************/

//  Rescale an image while leaving selected areas unchanged.

namespace SLrescale_names
{
   editfunc    EFslrescale;

   int      Fsetups = 0;
   int      dragx, dragy;
   int      E3ww, E3hh;
   char     *sqrow, *sqcol;
   int      Nsqrow, Nsqcol;
   int      *npx, *npy;

   int    dialog_event(zdialog *zd, cchar *event);
   void   setups();
   void   cleanups();
   void   mousefunc();
   void   warpfunc();
   void   *warpthread(void *);
}


//  menu function

void m_selective_rescale(GtkWidget *, cchar *)                                   //  16.05
{
   using namespace SLrescale_names;

   cchar  *message = ZTX(" Select areas to remain unchanged. \n"
                         " Pull image from upper left corner. \n"
                         " When finished, press [done].");

   F1_help_topic = "selective_rescale";

   EFslrescale.menufunc = m_selective_rescale;
   EFslrescale.funcname = "selective rescale";
   EFslrescale.Farea = 2;                                                        //  select area usable
   EFslrescale.mousefunc = mousefunc;                                            //  mouse function
   if (! edit_setup(EFslrescale)) return;                                        //  setup edit

   PXM_addalpha(E0pxm);
   PXM_addalpha(E1pxm);
   PXM_addalpha(E3pxm);

   zdialog *zd = zdialog_new(ZTX("Selective Rescale"),Mwin,Bproceed,Bdone,Bcancel,null);
   EFslrescale.zd = zd;
   zdialog_add_widget(zd,"label","lab1","dialog",message,"space=3");

   zdialog_run(zd,dialog_event,"save");                                          //  run dialog, parallel
   return;
}


//  dialog event and completion callback function

int SLrescale_names::dialog_event(zdialog * zd, cchar *event)
{
   using namespace SLrescale_names;

   if (strmatch(event,"done")) zd->zstat = 2;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 2;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 3;                                  //  from f_open()
   
   if (! zd->zstat) return 1;                                                    //  wait for completion
   
   if (zd->zstat == 1) {                                                         //  [proceed]
      zd->zstat = 0;                                                             //  keep dialog active
      if (sa_stat != 3)
         zmessageACK(Mwin,ZTX("select areas first"));
      else setups();                                                             //  start drag and warp
      return 1;
   }
   
   if (zd->zstat != 2 || dragx + dragy == 0) {                                   //  [cancel] or no change
      edit_cancel(0);
      cleanups();
      return 1;
   }

   edit_done(0);                                                                 //  [done]
   cleanups();
   m_trimrotate(0,"auto");                                                       //  set up automatic trim
   return 1;
}


//  do setups based on select area data

void SLrescale_names::setups()
{
   int      ii, spx, spy, sum;

   cleanups();                                                                   //  free prior if any

   dragx = dragy = 0;                                                            //  no drag data

   E3ww = E3pxm->ww;                                                             //  image dimensions
   E3hh = E3pxm->hh;

   sqrow = (char *) zmalloc(E3hh);                                               //  maps squishable rows/cols
   sqcol = (char *) zmalloc(E3ww);
   memset(sqrow,1,E3hh);                                                         //  mark all rows/cols squishable
   memset(sqcol,1,E3ww);
   
   for (spy = 0; spy < E3hh; spy++)                                              //  loop all source pixels
   for (spx = 0; spx < E3ww; spx++)
   {
      ii = spy * E3ww + spx;                                                     //  pixel within area?
      if (sa_pixmap[ii]) sqrow[spy] = sqcol[spx] = 0;                            //  mark row/col non-squishable
   }
   
   Nsqrow = Nsqcol = 0;
   for (spy = 0; spy < E3hh; spy++)                                              //  count total squishable rows/cols
      Nsqrow += sqrow[spy];
   for (spx = 0; spx < E3ww; spx++)
      Nsqcol += sqcol[spx];

   npx = (int *) zmalloc(E3ww * sizeof(int));                                    //  count of squishable rows/cols
   npy = (int *) zmalloc(E3hh * sizeof(int));                                    //    predeeding a given row/col
   
   for (sum = spx = 0; spx < E3ww; spx++)
   {
      if (sqcol[spx]) sum++;
      npx[spx] = sum;
   }

   for (sum = spy = 0; spy < E3hh; spy++)
   {
      if (sqrow[spy]) sum++;
      npy[spy] = sum;                                                            //  squishable rows < spy
   }

   Fsetups = 1;
   sa_unselect();                                                                //  delete area
   takeMouse(mousefunc,dragcursor);                                              //  connect mouse function
   return;
}


//  free allocated memory

void SLrescale_names::cleanups()
{
   if (! Fsetups) return;
   Fsetups = 0;
   zfree(sqrow);
   zfree(sqcol);   
   zfree(npx);
   zfree(npy);
   return;
}


//  mouse function

void SLrescale_names::mousefunc()
{
   using namespace SLrescale_names;
   
   float    R;
   
   if (Mxdrag || Mydrag)                                                         //  mouse drag underway
   {
      R = 1.0 * Mxdown / E3ww;                                                   //  ignore drag not from NW corner
      if (R > 0.2) return;
      R = 1.0 * Mydown / E3hh;
      if (R > 0.2) return;
      dragx = Mxdrag - Mxdown;                                                   //  drag amount
      dragy = Mydrag - Mydown;
      warpfunc();                                                                //  drag image
      Mxdrag = Mydrag = 0;
   }

   return;
}


//  warp function

void SLrescale_names::warpfunc()
{
   using namespace SLrescale_names;

   for (int ii = 0; ii < NWT; ii++)                                              //  start worker threads
      start_wthread(warpthread,&Nval[ii]);
   wait_wthreads();                                                              //  wait for completion

   CEF->Fmods++;
   CEF->Fsaved = 0;
   Fpaint2();                                                                    //  update window

   return;
}


//  warp thread

void * SLrescale_names::warpthread(void *arg)
{
   using namespace SLrescale_names;

   int      index = *((int *) (arg));
   int      spx, spy, dpx, dpy;
   float    Rx, Ry;
   float    *spix, *dpix;
   int      nc = E1pxm->nc, pcc = nc * sizeof(float);
   
   for (spy = index; spy < E3hh; spy += NWT)                                     //  loop all source pixels
   for (spx = 0; spx < E3ww; spx++)
   {
      if (spx < dragx || spy < dragy) {                                          //  pixels < dragx/dragy:
         spix = PXMpix(E3pxm,spx,spy);                                           //    black, transparent
         memset(spix,0,pcc);
      }

      Rx = 1.0 * npx[spx] / Nsqcol;                                              //  squishable pixel ratios, 0 - 1.0
      Ry = 1.0 * npy[spy] / Nsqrow;
      
      dpx = spx + dragx * (1.0 - Rx);                                            //  destination pixel
      dpy = spy + dragy * (1.0 - Ry);

      if (dpx < 0 || dpx > E3ww-1) continue;                                     //  necessary, why?  ///////////
      if (dpy < 0 || dpy > E3hh-1) continue;
      
      dpix = PXMpix(E3pxm,dpx,dpy);                                              //  source pixel >> destination pixel
      spix = PXMpix(E1pxm,spx,spy);
      memcpy(dpix,spix,pcc);
   }

   exit_wthread();
   return 0;
}


/********************************************************************************/

//  make waves menu function
//  distort the image with a wave pattern
//  independent horizontal and vertical wavelengths
//  variance option: make waves more or less irregular
//  perspective option: wavelengths lengthen going down

namespace waves_names
{
   editfunc EFwaves;
   int      ww, hh;                          //  image dimensions
   int      WLV, WLH;                        //  vertical and horizontal wavelengths
   int      AMPV, AMPH;                      //  vertical and horizontal amplitudes
   int      VARV, VARH;                      //  vertical and horizontal variance
   int      PERSP;                           //  perspective 0-100
}


//  menu function

void m_waves(GtkWidget *, const char *)
{
   using namespace waves_names;

   int    waves_dialog_event(zdialog* zd, const char *event);
   void * waves_thread(void *);

   F1_help_topic = "make_waves";

   EFwaves.menufunc = m_waves;
   EFwaves.funcname = "make_waves";
   EFwaves.Farea = 2;                                                            //  select area usable
   EFwaves.threadfunc = waves_thread;

   if (! edit_setup(EFwaves)) return;

/***
          ___________________________________
         |           Make waves             |
         |                                   |
         |             Horizontal  Vertical  |
         |  wavelength [___|-+]    [___|-+]  |
         |  amplitude  [___|-+]    [___|-+]  |
         |  variance   [___|-+]    [___|-+]  |
         |                                   |
         |  perspective [___|-+]             |
         |                                   |
         |          [apply] [done] [cancel]  |
         |___________________________________|

***/

   ww = E3pxm->ww;
   hh = E3pxm->hh;

   zdialog *zd = zdialog_new(ZTX("Make Waves"),Mwin,Bapply,Bdone,Bcancel,null);
   EFwaves.zd = zd;
   zdialog_add_widget(zd,"hbox","hbw","dialog",0,"space=3");
   zdialog_add_widget(zd,"vbox","vbw1","hbw",0,"space=3");
   zdialog_add_widget(zd,"vbox","vbw2","hbw",0,"space=3");
   zdialog_add_widget(zd,"vbox","vbw3","hbw",0,"space=3");
   zdialog_add_widget(zd,"label","space","vbw1"," ","space=1");
   zdialog_add_widget(zd,"label","labwl","vbw1",ZTX("wavelength"),"expand");
   zdialog_add_widget(zd,"label","labamp","vbw1",ZTX("amplitude"),"expand");
   zdialog_add_widget(zd,"label","labamp","vbw1",ZTX("variance"),"expand");
   zdialog_add_widget(zd,"label","labh","vbw2",ZTX("horizontal"),"space=1");
   zdialog_add_widget(zd,"spin","wlh","vbw2","3|500|1|50","expand");
   zdialog_add_widget(zd,"spin","amph","vbw2","0|100|1|20","expand");
   zdialog_add_widget(zd,"spin","varh","vbw2","0|100|1|20","expand");
   zdialog_add_widget(zd,"label","labh","vbw3",ZTX("vertical"),"space=1");
   zdialog_add_widget(zd,"spin","wlv","vbw3","3|500|1|50","expand");
   zdialog_add_widget(zd,"spin","ampv","vbw3","0|100|1|20","expand");
   zdialog_add_widget(zd,"spin","varv","vbw3","0|100|1|20","expand");
   zdialog_add_widget(zd,"hsep","sepp","dialog",0,"space=3");
   zdialog_add_widget(zd,"hbox","hbp","dialog",0,"space=2");
   zdialog_add_widget(zd,"label","labp","hbp",ZTX("perspective"),"space=3");
   zdialog_add_widget(zd,"spin","persp","hbp","0|100|1|0","space=5");

   zdialog_restore_inputs(zd);                                                   //  restore previous inputs
   zdialog_run(zd,waves_dialog_event,"save");                                    //  run dialog - parallel
   return;
}


//  waves dialog event and completion function

int waves_dialog_event(zdialog *zd, const char *event)                           //  waves dialog event function
{
   using namespace waves_names;

   if (strmatch(event,"done")) zd->zstat = 2;                                    //  from edit_setup() or f_save()
   if (strmatch(event,"enter")) zd->zstat = 2;                                   //  KB input
   if (strmatch(event,"cancel")) zd->zstat = 3;                                  //  from f_open()

   if (zd->zstat)
   {
      if (zd->zstat == 1) {                                                      //  apply
         zd->zstat = 0;                                                          //  keep dialog active
         edit_reset();
         zdialog_fetch(zd,"wlv",WLV);                                            //  get user inputs
         zdialog_fetch(zd,"wlh",WLH);
         zdialog_fetch(zd,"ampv",AMPV);
         zdialog_fetch(zd,"amph",AMPH);
         zdialog_fetch(zd,"varv",VARV);
         zdialog_fetch(zd,"varh",VARH);
         zdialog_fetch(zd,"persp",PERSP);
         signal_thread();                                                        //  calculate
         return 1;
      }
      if (zd->zstat == 2) edit_done(0);                                          //  commit edit
      else edit_cancel(0);                                                       //  discard edit
   }

   return 1;
}


//  thread function - multiple working threads to update image

void * waves_thread(void *)
{
   using namespace waves_names;

   void * waves_wthread(void *);

   while (true)
   {
      thread_idle_loop();                                                        //  wait for work or exit request

      for (int ii = 0; ii < NWT; ii++)                                           //  start worker threads
         start_wthread(waves_wthread,&Nval[ii]);
      wait_wthreads();                                                           //  wait for completion

      CEF->Fmods++;                                                              //  image modified
      CEF->Fsaved = 0;                                                           //  not saved
      Fpaint2();                                                                 //  update window
   }

   return 0;                                                                     //  not executed, stop warning
}


void * waves_wthread(void *arg)                                                  //  working threads
{
   using namespace waves_names;

   int      index = *((int *) (arg));
   float    wlv, wlh, ampv, amph, varv, varh, persp;
   int      py, px, pylo, pyhi;
   int      ii, edist = 0;
   float    dx, dy, *pix1, *pix3, vpix[4];
   float    f1, f2, ffv, ffh, red, green, blue;

   wlv = WLV;                                                                    //  vertical and horizontal wavelengths
   wlh = WLH;
   ampv = 0.01 * AMPV * wlv;                                                     //  vertical and horizoneal amplitudes
   amph = 0.01 * AMPH * wlh;
   varv = 0.01 * VARV;                                                           //  vertical and horizontal variance
   varh = 0.01 * VARH;
   persp = 0.01 * PERSP;                                                         //  perspective 0 to 1.0

   pylo = 0;                                                                     //  py range
   pyhi = hh;

   if (sa_stat == 3) {
      pylo = sa_miny - wlv;                                                      //  rescale to select area height
      if (pylo < 0) pylo = 0;
      pyhi = sa_maxy + wlv;
      if (pyhi > hh) pyhi = hh;
   }

   for (py = pylo + index; py < pyhi; py += NWT)
   {
      if (persp == 0) ffv = ffh = 1.0;
      else {
         ffv = 2.0 * persp * (py - pylo) / (pyhi - pylo);                        //  0 ... 2.0
         ffv = 1.0 + ffv * ffv;                                                  //  1 ... 5.0
         ffh = sqrtf(ffv);                                                       //  1 ... 2.23
      }

      dx = ampv * sin(2.0 * PI * py / (ffv * wlv));

      if (varv > 0)
         dx += ampv * varv * sin(2.0 * PI * py / (1.37 * ffv * wlv));            //  make irregular wave

      for (px = 0; px < ww; px++)
      {
         if (sa_stat == 3) {                                                     //  select area active
            ii = py * ww + px;
            edist = sa_pixmap[ii];                                               //  distance from edge
            if (! edist) continue;                                               //  pixel outside area
         }

         dy = amph * sin(2.0 * PI * px / (ffh * wlh));

         if (varh > 0)
            dy += amph * varh * sin(2.0 * PI * px / (1.37 * ffh * wlh));         //  make irregular wave

         vpixel(E1pxm,px+dx,py+dy,vpix);                                         //  source pixel, displaced
         red = vpix[0];
         green = vpix[1];
         blue = vpix[2];

         if (sa_stat == 3 && edist < sa_blend) {                                 //  blend edges
            f1 = sa_blendfunc(edist);                                            //  16.08
            f2 = 1.0 - f1;
            pix1 = PXMpix(E1pxm,px,py);
            red = f1 * red + f2 * pix1[0];
            green = f1 * green + f2 * pix1[1];
            blue = f1 * blue + f2 * pix1[2];
         }

         pix3 = PXMpix(E3pxm,px,py);                                             //  destination pixel
         pix3[0] = red;
         pix3[1] = green;
         pix3[2] = blue;
      }
   }

   exit_wthread();                                                               //  exit thread
   return 0;                                                                     //  not executed, stop gcc warning
}