/* Last edited on 2023-02-22 12:28:26 by stolfi */

#include <pz_types.h>
#include <pz_proc.h>
#include <pz_plot.h>
#include <pz_plot_gauge.h>

#include <r3.h>
#include <r2.h>
#include <pswr.h>
#include <argparser.h>
#include <interval.h>

#include <math.h>
#include <stdio.h>
#include <string.h>
#include <values.h>

#define pt_per_mm (72.0/25.4)

r2_t pz_plot_gauge_graph_point( int i, double step, double lambda, double height );

r2_t pz_plot_gauge_graph_point( int i, double step, double lambda, double height )
  { double zp = ((double)i)*step;
    double xp = zp;
    double tp = zp/lambda;
    double yp = height * exp(-tp*tp/2.0);
    return (r2_t){{xp, yp}};
  }

void pz_plot_gaussian_kernel
  ( PSStream *f,
    double lambda,
    double x,
    double y,
    double width,
    double height,
    double lineWidth,
    double dotStep,   /* (:= 0.0) */
    double dotSize    /* (:= 0.0) */
  )
  {
    r2_t a;  /* Previous point on graph */
    int i;
    
    if (lineWidth > 0.0)
      { pswr_set_pen(f, 0.0,0.0,0.5, lineWidth, 0.0,0.0);
        double plotStep = lambda/12.0;
        int nSteps = (int)floor((width/2.0) / plotStep);
        /* Draw the graph: */
        for (i = -nSteps; i <= +nSteps ; i++)
          { r2_t p = pz_plot_gauge_graph_point(i, plotStep, lambda, height);
            if (i > -nSteps){ pswr_segment(f, x+a.c[0], y+a.c[1], x+p.c[0], y+p.c[1]); }
            a = p;
          }
      }

    if ((dotStep > 0.0) && (dotSize > 0.0))
      { /* Draw the dots: */
        int nSteps = (int)floor((width/2.0) / dotStep);
        for (i = -nSteps; i <= +nSteps ; i++)
          { r2_t p = pz_plot_gauge_graph_point(i, dotStep, lambda, height);
            pswr_dot(f, x+p.c[0], y+p.c[1], dotSize, TRUE, FALSE);
          }
      }
  }

r2_t pz_plot_gauge_circle_point ( int i, double step, double radius );

r2_t pz_plot_gauge_circle_point ( int i, double step, double radius )
  { double theta = ((double)i)*step/radius;
    double xp = radius*sin(theta);
    double yp = radius*cos(theta);
    return (r2_t){{xp, yp}};
  }

void pz_plot_critical_circle
  ( PSStream *f,
    double radius,
    double x,
    double y,
    double width,
    double height,
    double lineWidth,
    double dotStep,   /* (:= 0.0) */
    double dotSize    /* (:= 0.0) */
  )
  {
    double xc, yc; /* Circle center */
    r2_t a;        /* Previous point on graph */
    int i;
    
    auto bool_t pz_plot_gauge_circle_inside( r2_t *p );
    bool_t pz_plot_gauge_circle_inside( r2_t *p )
      { return 
          (fabs(p->c[0] - x) <= width/2.0) && 
          (p->c[1] >= y) && 
          (p->c[1] <= y + height);
      }

    double halfHeight = height/2.0;
    
    /* Assumes width >= height; will not look good if width < height. */
    xc = x;
    if (radius <= halfHeight)
      { yc = y + halfHeight; }
    else
      { yc = y + height - radius; }
    
    if (lineWidth > 0.0)
      {
        pswr_set_pen(f, 0.0,0.0,0.5, lineWidth, 0.0,0.0);
        int nSteps = 32;
        double plotStep = Pi*radius/((double)nSteps);
        /* Draw the graph: */
        for (i = -nSteps; i <= +nSteps ; i++)
          { r2_t p = pz_plot_gauge_circle_point(i, plotStep, radius);
            bool_t a_inside = pz_plot_gauge_circle_inside(&a);
            bool_t p_inside = pz_plot_gauge_circle_inside(&p);
            if ((i > -nSteps) && a_inside && p_inside)
              { pswr_segment(f, xc+a.c[0], yc+a.c[1], xc+p.c[0], yc+p.c[1]); }
            a = p;
          }
      }

    if ((dotStep > 0.0) && (dotSize > 0.0))
      { /* Draw the dots: */
        int nSteps = (int)floor(M_PI*radius / dotStep);
        for (i = -nSteps; i <= +nSteps ; i++)
          { r2_t p = pz_plot_gauge_circle_point(i, dotStep, radius);
            pswr_dot(f, xc+p.c[0], yc+p.c[1], dotSize, TRUE,FALSE);
          }
      }
  }

void pz_plot_scale_bar(
    PSStream *f,
    double barStep,
    unsigned barCount,
    double x,
    double y,
    double lineWidth,
    double tickLength 
  )
  { if (lineWidth > 0.0)
      { pswr_set_pen(f, 0.0,0.0,0.5, lineWidth, 0.0,0.0);
        double tickExtent = ((double)tickLength/2.0 * pt_per_mm) / f->yScale;
        double halfWidth = ((double)barCount)*barStep/2.0;
        int i;
        pswr_segment(f, x - halfWidth, y, x + halfWidth, y);
        for (i = 0; i <= barCount ; i++)
          { double xTick = x - halfWidth + ((double)i)*barStep;
            pswr_segment(f, xTick, y - tickExtent, xTick, y + tickExtent);
          }
      }
  }

#define CriticalAlpha (0.8577638846)  /* sqrt(2/e) */ 

void pz_plot_gauge ( PSStream *f, double lambda, pz_plot_gauge_options_t *o )
  {
    double iconWidth, iconHeight;  /* Dimensions of icon (minus scale bar). */ 
    double iconLift;               /* Y distance from base of bar to base of icon. */
    double barWidth, barHeight;    /* Dimensions of scale bar. */
    double totWidth, totHeight;    /* Dimensions of icon plus scale bar. */ 
    if ((o->style == NoGauge) && (o->barCount == 0)) { return; }

    /* Scale factors on each axis (pt/client_unit): */
    double xScale = f->xScale;
    double yScale = f->yScale;
    /* Plotting area (client coords): */
    interval_t xRange = (interval_t){{f->xMin, f->xMax}};
    interval_t yRange = (interval_t){{f->yMin, f->yMax}};
    /* Leave 2mm of blank space around drawing: */
    double xMargin = 2.0 * pt_per_mm / xScale;
    double yMargin = 2.0 * pt_per_mm / yScale;
    /* Leave another 2mm of distance from the drawing's border: */
    double xSafety = 2.0 * pt_per_mm / xScale;
    double ySafety = 2.0 * pt_per_mm / yScale;
    
    /* Compute height, depth, and width of scale bar: */
    if (o->barCount > 0) 
      { /* The tick marks of the scale bar will extend below "y": */
        barHeight = ((double)o->tickLength) / yScale;
        barWidth = (0.5 + ((double)o->barCount)) * lambda;
      }
    else
      { barHeight = 0.0; barWidth = 0.0; }

    /* Compute ideal height and width of icon (exclusive of scale bar): */
    switch(o->style)
      {
        case NoGauge:
          iconWidth = 0.0;
          iconHeight = 0.0;
          iconLift = 0.0;
          break;
        case KernelGauge: 
          iconWidth = fmax(5.0, ((double)o->barCount)) * lambda;
          iconHeight = 8.0 / yScale;
          iconLift = 0.5 * barHeight;
          break;
        case CircleGauge: 
          iconWidth = 2.0 * (lambda/CriticalAlpha);
          iconHeight = 8.0 / yScale;
          iconLift = barHeight + 0.5 * yMargin;
          break;
        default:
          demand(FALSE, "invalid gauge style");
          iconWidth =  0.0;
          iconHeight = 0.0;
          iconLift = 0.0;
      }

    /* Compute plotting icon and scale bar width: */
    /* Don't allow the them to be too big: */

    { double totWidthMax = (xRange.end[1] - xRange.end[0] - 2.0*(xSafety+xMargin))/3.0;
      double totHeightMax = (yRange.end[1] - yRange.end[0] - 2.0*(ySafety+yMargin))/3.0;
      totWidth = fmin(totWidthMax, fmax(iconWidth, barWidth));
      totHeight = fmin(totHeightMax, fmax(iconLift + iconHeight, barHeight));
    }

    /* Plot something only if one of them fits: */
    if
      ( ( (iconWidth <= totWidth) && (iconHeight <= totHeight - iconLift) ) || 
        ( (barWidth <= totWidth) && (barHeight <= totHeight) )
      )
      {
        /* Compute midpoint of bottom side of icon/bar area: */
        double x = xRange.end[1] - xSafety - xMargin - totWidth/2.0;
        double y = yRange.end[0] + ySafety + yMargin;
        /* Blank out the background: */
        pswr_set_fill_color(f, 1.0,1.0,1.0);
        double dx = totWidth/2.0 + xMargin;
        pswr_rectangle(f, x - dx, x + dx, y - yMargin, y + totHeight + yMargin, TRUE,FALSE);
        if ((iconWidth <= totWidth) && (iconHeight <= totHeight - iconLift))
          { /* Plot the kernel or circle, with lines: */
            pswr_set_pen(f, 0.0,0.0,0.0, 0.15, 0.0,0.0);
            switch(o->style)
              {
                case NoGauge: 
                  /* OK */
                  break;
                case KernelGauge: 
                  pz_plot_gaussian_kernel
                    ( f, 
                      /* lambda */    lambda,
                      /* x */         x, 
                      /* y */         y + iconLift,
                      /* width */     iconWidth, 
                      /* height */    totHeight - iconLift,
                      /* lineWidth */ 0.15,
                      /* dotStep */   o->dotStep, 
                      /* dotSize */   2.0
                    );
                  break;
                case CircleGauge: 
                  pz_plot_critical_circle
                    ( f, 
                      /* radius */    lambda / CriticalAlpha,
                      /* x */         x, 
                      /* y */         y + iconLift,
                      /* width */     totWidth, 
                      /* height */    totHeight - iconLift,
                      /* lineWidth */ 0.15,
                      /* dotStep */   o->dotStep,
                      /* dotSize */   2.0
                    );
                break;
              }
          }
        if ((barWidth <= totWidth) && (barHeight <= totHeight))
          { pswr_set_pen(f, 0.0,0.0,0.0, 0.15, 0.0,0.0);
            pz_plot_scale_bar
              ( f, 
                /* barStep */    lambda, 
                /* barCount */   o->barCount,
                /* x */          x,
                /* y */          y + 0.5 * barHeight,
                /* lineWidth */  0.15,
                /* tickLength */ o->tickLength
              );
          }
      }
  }

pz_plot_gauge_options_t *pz_plot_parse_gauge_options ( argparser_t *pp )
  {
    pz_plot_gauge_options_t *o = malloc(sizeof(pz_plot_gauge_options_t));
    if (argparser_keyword_present(pp, "-scaleBar"))
      { o->style = NoGauge;
        o->barCount = argparser_get_next_int(pp, 0, 100);
        o->dotStep = 0.0;
      }
    else if (argparser_keyword_present(pp, "-scaleIcon"))
      { char *s = argparser_get_next(pp);
        if (strcmp(s, "None") == 0)
          { o->style = NoGauge; o->barCount = 0; }
        else if (strcmp(s, "Bar") == 0)
          { o->style = NoGauge; o->barCount = 1; }
        else if (strcmp(s, "Circle") == 0)
          { o->style = CircleGauge; o->barCount = 1; }
        else if (strcmp(s, "Kernel") == 0)
          { o->style = KernelGauge; o->barCount = 2; }
        else
          { argparser_error(pp, "bad scale scale icon style"); }
        if (argparser_keyword_present(pp, "-scaleDotStep"))
          { o->dotStep = argparser_get_next_double(pp, 0.0, DBL_MAX); }
        else
          { o->dotStep = 0.0; }
      }
    else
      { o->style = NoGauge;
        o->barCount = 0;    /* No scale bar */
        o->dotStep = 0.0;
      }
    o->tickLength = 1.5;
    return o;
  }