/* See figs1d.h. */
/* Last edited on 2011-09-15 18:00:42 by stolfilocal */

#define _GNU_SOURCE
#include <stdio.h>
#include <math.h>
#include <stdlib.h>
#include <string.h>
#include <values.h>

#include <set32.h>
#include <interval.h>
#include <pswr.h>
#include <jsmath.h>

// #include <dg_grid.h>
// #include <dg_locus.h>
// #include <dg_pulse.h>
// #include <dg_tent.h>
// #include <dg_tent_basis.h>

#include <plot1d.h>
#include <figtools.h>
#include <somegrids.h>
  
#include <figs1d.h>

#define DIM 1
#define MAXTAGSZ 200 
  
void figs_1d_mother_pulse_catalog 
  ( PlotOptions *o,
    dg_pulse_kind_t pkind,
    dg_cont_t cMin, 
    dg_cont_t cMax,
    dg_degree_t gOver
  )
  {  
    char pkch = dg_pulse_kind_to_char(pkind);
    dg_cont_t c;
    for (c = cMin; c <= cMax; c++)
      { /* The variable {gMin} holds the min degree valid for {pkind,c} found so far: */
        dg_degree_t gMin = dg_PULSE_MAX_DEGREE + 1;  /* For now. */
        /* Loop over all degrees {g}: */
        dg_degree_t g;
        for (g = 0; g <= gMin + gOver; g++)
          { dg_pulse_family_t fam = dg_pulse_family(pkind, c, g);
            if (fam.nmp > 0)
              { /* Found a valid degree {g} for the pair {pkind,c}. */
                if (g < gMin) { /* First valid {g}: */ gMin = g; }
                /* Get the maximum mother pulse support size {mszMax} for this family: */
                dg_grid_size_t mszMax = dg_pulse_max_supp_count(&fam);
                /* Loop over grid sizes to show effect of self-overlap: */
                dg_rank_t r = 0; /* Grid level. */
                while(TRUE) 
                  { /* Get the grid size {gsz} along each axis: */
                    dg_grid_size_t gsz;
                    dg_grid_size(DIM, r, &gsz);
                    /* Get the min grid size {gszMin} along any axis: */
                    dg_grid_size_t gszMin = ipow(2, r/DIM);
                    /* Plot one figure with all mother pulses of {fam} on {gsz}. */
                    /* Use the same arg and value ranges for all plots with same {c,gsz}: */
                    interval_t fr[2];
                    compute_1d_common_pulse_ranges(c, gsz, 1.0, /*out*/ fr);
                    char *tag = NULL;
                    asprintf(&tag, "d1-k%c-g%d-c%d-r%02d", pkch, g, c, r);
                    dg_pulse_family_t fam = dg_pulse_family(pkind, c, g); 
                    fig_1d_mother_pulse_set(o, tag, &fam, gsz, fr);
                    free(tag);
                    if (gszMin > mszMax) { break; }
                    r++;
                  }
              }
          }
      }
  }

void fig_1d_mother_pulse
  ( PlotOptions *o,
    char *tag,                   /* Tag for figure/file name. */
    dg_pulse_family_t *fam,      /* Pulse family. */
    dg_pulse_mother_index_t pix, /* Index of mother element. */
    dg_grid_size_t gsz,          /* Domain period (number of grid intervals). */
    interval_t fr[]              /* Arg and function ranges of interest. */
  )
  { if (! check_pulse_index(fam, pix)) { return; }
    
    dg_pulse_t p = dg_pulse(fam, pix, gsz, 0);   /* The specific pulse. */
    interval_t R[1];
    R[0] = (interval_t){{ 0, (double)gsz }};      /* The root interval. */
    
    /* Nominal plot scale (client coords to mm): */
    double scale = o->scale;
    
    double mag[2];    /* Graphing scales (client uns per function un). */
    double tstep[2];  /* Major tick spacing. */
    double sstep[2];  /* Minor tick spacing. */
    mag[0] = 1.0; tstep [0]= 1.0; sstep [0]= 0.5; /* For function argument. */
    mag[1] = 0.5; tstep [1]= 1.0; sstep [1]= 0.5; /* For function values. */
    
    interval_t gr[2];   /* Slightly widened ranges for plot/clip. */
    interval_t win[2];  /* Plot window (in client coords). */
    
    /* Start figure, plot axes etc: */
    compute_1d_plot_scales(ps, o, fr, scale, mag, /* OUT */ gr, win);
    start_figure(o, tag, tag, win, scale);
    plot_1d_graph_decoration
      ( ps, o, scale, 
        mag, fr, tstep, sstep, gr, win
      );
    
    /* Plot the function: */
    if (o->color)
      { pswr_set_pen(ps, 0.0,0.0,0.5, 0.15, 0.0, 0.0); }
    else
      { pswr_set_pen(ps, 0.0,0.0,0.0, 0.10, 0.0, 0.0); }
    
    auto double func(double x);

    double func(double x)
      { double fx;
        dg_pulse_eval(fam, &p, &(R[0]), x, 0, &fx);
        return fx;
      }
    
    plot_1d_function(ps, o, func,
      0.002, gr[0], mag[0], gr[1], mag[1]
    );

    finish_figure(o);
  }
  
void fig_1d_mother_pulse_set
  ( PlotOptions *o, 
    char *tag,              /* Tag for figure/file name. */
    dg_pulse_family_t *fam, /* The pulse's family. */
    dg_grid_size_t gsz,     /* Domain period (number of grid intervals). */
    interval_t fr[]         /* Arg and function ranges of interest. */
  )
  { 
    if (! check_pulse_index(fam, 0)) { return; }
    /* Get max valid pulse index {mix}: */
    int mix = fam->nmp - 1;
    
    /* Nominal plot scale (client coords to mm): */
    double scale = o->scale;
    
    interval_t R[1];
    R[0] = (interval_t){{ 0, (double)gsz }};      /* The root interval. */
    
    double mag[2];    /* Graphing scales (client uns per function un). */
    double tstep[2];  /* Major tick spacing. */
    double sstep[2];  /* Minor tick spacing. */
    mag[0] = 1.0; tstep [0]= 1.0; sstep [0]= 0.5; /* For function argument. */
    mag[1] = 0.5; tstep [1]= 1.0; sstep [1]= 0.5; /* For function values. */
    
    interval_t gr[2];      /* Slightly widened ranges for plot/clip. */
    interval_t win[2];     /* Plot window (in client coords). */
    
    /* Start figure, plot axes etc: */
    compute_1d_plot_scales(ps, o, fr, scale, mag, /* OUT */ gr, win);
    start_figure(o, tag, tag, win, scale);
    plot_1d_graph_decoration
      ( ps, o, scale, 
        mag, fr, tstep, sstep, gr, win
      );
    
    /* Plot the function: */
    if (o->color)
      { pswr_set_pen(ps, 0.0,0.0,0.5, 0.15, 0.0, 0.0); }
    else
      { pswr_set_pen(ps, 0.0,0.0,0.0, 0.10, 0.0, 0.0); }
    dg_pulse_mother_index_t pix;
    for (pix = 0; pix <= mix; pix++)
      {
        dg_pulse_t p = dg_pulse(fam, pix, gsz, 0);   /* The specific pulse. */

        auto double func(double x);

        double func(double x)
          { double fx;
            dg_pulse_eval(fam, &p, &(R[0]), x, 0, &fx);
            return fx;
          }

        plot_1d_function(ps, o, func, 
          0.002, gr[0], mag[0], gr[1], mag[1]
        );
      }

    finish_figure(o);
  }

void figs_1d_some_grid_both_bases 
  ( PlotOptions *o,
    dg_rank_t rMin,
    dg_rank_t rMax,
    dg_pulse_kind_t pkind,
    dg_cont_t c,
    dg_degree_t g
  )
  {
    interval_t rootCell[1];
    dg_cell_box_canonical(1, dg_ROOT_CELL, rootCell);

    /* Build a finite dyadic grid: */
    char *tree_tag = NULL;
    asprintf(&tree_tag, "%02d:%02d", rMin, rMax);
    dg_tree_node_t *G = somegrid_varied(1, rMin, rMax, NULL);
    
    /* Choose a pulse family for the basis: */
    dg_pulse_family_t fam = dg_pulse_family(pkind, c, g);

    /* Plot basis: */
    figs_1d_basis(o, tree_tag, G, rootCell, &fam, FALSE);
    figs_1d_basis(o, tree_tag, G, rootCell, &fam, TRUE);
    
    free(tree_tag);
  }

void figs_1d_basis
  ( PlotOptions *o, 
    char *tree_tag,          /* A character strings that identifies the tree. */
    dg_tree_node_t *G,       /* Root node of grid. */
    interval_t rootCell[],   /* Root interval. */
    dg_pulse_family_t *fam,  /* Pulse family. */
    bool_t superior          /* TRUE for maximal-support elements. */
  )
  {
    char pkch = dg_pulse_kind_to_char(fam->pkind);
    char *si_tag = (superior ? "SUP" : "INF");
    /* Generate the filename tag: */
    char *basis_tag = NULL;
    asprintf(&basis_tag, "d1-t%s-k%c-g%d-c%d-%s", tree_tag, pkch, fam->g, fam->c, si_tag );
    
    /* Get max valid pulse index {mix}: */
    if (fam->nmp <= 0)
      { fprintf(stderr, "** family has no mothers");
        dg_pulse_family_print(stderr, fam);
        fprintf(stderr, "\n");
        return;
      }
    
    int i;
    
    /* Assuming that there is no leaf at level 0: */   
    dg_grid_size_t szMin = (superior ? 1 : 2);
    double wd = 1.0;
    
    /* Obtain the function arg/value ranges of interest: */
    interval_t fr[2];
    compute_1d_common_pulse_ranges(fam->c, szMin, wd, /*out*/ fr);
    fr[0] = rootCell[0];
      
    dg_tent_vec_t tv = dg_tent_basis_get(DIM, fam, G, superior);
    fprintf(stderr, "  basis has %d tents\n", tv.nel);
    
    /* Ensure that the basis starts on a new page: */
    pswr_sync_canvas(ps, NULL);
    
    /* Plot the tents: */
    for (i = 0; i < tv.nel; i++)
      { dg_tent_t t = tv.el[i];
      
        dg_pulse_t p[DIM];
        dg_tent_unpack(DIM, fam, &t, p);
        dg_rank_t r = dg_cell_rank(t.ck);
        dg_grid_pos_t pos;
        dg_cell_position(DIM, t.ck, &pos);
      
        fprintf(stderr, "[%03d] ", i); 
        dg_tent_print(stderr, DIM, fam, p);
        fprintf(stderr, "\n"); 
        
        char *elTag = NULL;
        asprintf(&elTag, "%s-e%05d", basis_tag, i);
        fig_1d_basis_element(o, elTag, G, rootCell, fam, p[0].pix, r, pos, fr);
        free(elTag);
      }

    /* Ensure that subsequent pictures will go to the next canvas: */
    pswr_fill_canvas(ps);
    free(basis_tag);
  }

void fig_1d_basis_element
  ( PlotOptions *o,
    char *tag,                   /* Tag for figure/file name. */
    dg_tree_node_t *G,           /* Root node of grid. */
    interval_t rootCell[],       /* The root cell. */
    dg_pulse_family_t *fam,      /* Pulse family. */
    dg_pulse_mother_index_t pix, /* The mother pulse index */
    dg_rank_t r,                 /* The pulse's rank. */
    dg_grid_pos_t pos,           /* The pulse's shift. */
    interval_t fr[]              /* Arg and function ranges of interest. */
  )
  { if (! check_pulse_index(fam, pix)) { return; }
    
    /* Nominal plot scale (client coords to mm): */
    double scale = o->scale;
    
    double mag[2];    /* Graphing scales (client uns per function un). */
    double tstep[2];  /* Major tick spacing. */
    double sstep[2];  /* Minor tick spacing. */
    mag[0] = 1.00; tstep [0]= 1.00; sstep [0]= 0.125; /* For function argument. */
    mag[1] = 0.25; tstep [1]= 1.00; sstep [1]= 0.500; /* For function values. */
    
    interval_t gr[2];   /* Slightly widened ranges for plot/clip. */
    interval_t win[2];  /* Plot window (in client coords). */
    
    /* Start figure, plot axes etc: */
    compute_1d_plot_scales(ps, o, fr, scale, mag, /* OUT */ gr, win);
    start_figure(o, tag, tag, win, scale);
    plot_1d_graph_decoration
      ( ps, o, scale, 
        mag, fr, tstep, sstep, gr, win
      );
    
    /* Plot the leaf cell boundaries: */
    if (o->color)
      { pswr_set_pen(ps, 0.75,0.75,0.75, 0.15, 0.0,0.0); }
    else
      { pswr_set_pen(ps, 0.0,0.0,0.0, 0.10, 0.75,0.75); }
    plot_1d_tree_leaves(ps, o, G, rootCell, gr[0], mag[0], gr[1], mag[1]);
    
    /* Get the grid's size: */
    dg_grid_size_t gsz;
    dg_grid_size(DIM, r, &gsz);
    
    /* Create the 1-dimensional tent: */
    dg_pulse_t p[DIM];
    p[0] = dg_pulse(fam, pix, gsz, 0);

    if ((p[0].msz % 2) == 0) 
      { /* Even support size, plot the tent's center: */
        dg_grid_pos_t ctrpos = (pos + p[0].msz/2) % gsz;
        dg_cell_index_t ucell = dg_cell_from_position(DIM, r, &ctrpos);
        dg_locus_t E = dg_locus(set32_range(0,DIM-1), ucell); 
        plot_1d_locus(ps, o, E, G, rootCell, gr[0], mag[0]); 
      }

    /* Derivative orders desired in evaluation: */
    dg_degree_t ord[DIM];
    ord[0] = 0;
    
    /* Plot the function: */

    auto double func(double x);

    double func(double x)
      { double fx;
        dg_tent_eval(DIM, fam, p, rootCell, &x, ord, &fx);
        return fx;
      }

    if (o->color)
      { switch (pix)
          {
            case 0: pswr_set_pen(ps, 0.0,0.0,0.6, 0.15, 0.0, 0.0); break;
            case 1: pswr_set_pen(ps, 0.0,0.3,0.0, 0.15, 0.0, 0.0); break;
            case 2: pswr_set_pen(ps, 0.5,0.0,0.0, 0.15, 0.0, 0.0); break;
            default: pswr_set_pen(ps, 0.0,0.0,0.0, 0.15, 0.0, 0.0);
          }
      }
    else
      {  switch(pix)
          {
            case 0: pswr_set_pen(ps, 0.0,0.0,0.0, 0.10, 0.0, 0.0); break;
            case 1: pswr_set_pen(ps, 0.0,0.0,0.0, 0.10, 0.5, 1.0); break;
            case 2: pswr_set_pen(ps, 0.0,0.0,0.0, 0.10, 0.5, 0.5); break;
            default: pswr_set_pen(ps, 0.0,0.0,0.0, 0.10, 0.0, 0.0);
          }
      }
    plot_1d_function(ps, o, func, 
      0.002, gr[0], mag[0], gr[1], mag[1]
    );

    finish_figure(o);
  }