/* Last edited on 2012-12-12 01:09:53 by stolfilocal */

    if (it > 0)
      { pswr_set_fill_color(*psP, 1.000, 0.000, 0.000); 
        pswr_set_label_font(*psP, "Helvetica", 14);
        char *title = NULL;
        asprintf(&title, "After %d iterations", it);
        pswr_label(*psP, title, 0.5*(xmin+xmax), ymax + 0.5*ymrg_top, 0.5, 0.5);
        free(title);
      }
    

void figs_poisson_convergence_plot
  ( PSStream **psP, 
    plot_options_t *o,
    double Y[],
    int nit, 
    char *tag,
    int level,
    int freq
  );
  /* Plots {(it,Y[it])} for {it} in {0..nit-1}, in log-Y scale. */


void figs_poisson_convergence_plot
  ( PSStream **psP, 
    plot_options_t *o,
    double Y[],
    int nit, 
    char *tag,
    int level,
    int freq
  )
  {
 
    /* Define figure dimensions: */
    double xsz = 160;              /* X width of actual plot area. */
    double ysz = 120;              /* Y height actual plot area. */
    
    double xmrg_left = 20;         /* Margin at left side of plot. */
    double xmrg_rite =  8;         /* Margin at right side of plot. */
    double ymrg_bot =  12;         /* Margin at bottom of plot. */
    double ymrg_top =   6;         /* Margin at top of plot. */

    /* Define graph area proper {win[0..1]}: */
    double xmin = xmrg_left;       /* X start of plot area. */     
    double xmax = xmin + xsz;      /* X end of plot area. */  
    double ymin = ymrg_bot;        /* Y start of plot area. */
    double ymax = ymin + ysz;      /* Y end of plot area. */  
    interval_t win[2] = { (interval_t){{ xmin, xmax }}, (interval_t){{ ymin, ymax }} };
    
    double xtot = xsz + xmrg_left + xmrg_rite;  /* Width of figure. */
    double ytot = ysz + ymrg_bot + ymrg_top;    /* Height of figure. */
    interval_t bbox[2] = { (interval_t){{ 0, xtot }}, (interval_t){{ 0, ytot }} };

    /* Open Postscript file: */
    char *fulltag = NULL;
    asprintf(&fulltag, "%s-v%02d-f%03d", tag, level, freq);
    start_figure(psP, o, fulltag, NULL, bbox, 1.0);

    /* Convert error magnitude {Y} to {log_2} scale: */
    double eps = 1.0e-100; /* A small value outside of the plot range. */
    double tol = 1.0e-8;   /* A very small error. */
    double Xplot[nit];
    double Yplot[nit];
    int it;
    double YplotMin = log(tol)/M_LN10;
    double YplotMax = -INF;
    for (it = 0; it < nit; it++)
      { Xplot[it] = it+1;
        Yplot[it] = log(fmax(Y[it],eps))/M_LN10;
        if (Yplot[it] > YplotMax) { YplotMax = Yplot[it]; }
        if ((it < 64) || ((it % 500) == 0))
          { fprintf(stderr, "%6d %10.5e  %8.0f %8.5f\n", it, Y[it], Xplot[it], Yplot[it]); }
      }
      
    /* Choose plot scales: */
    double XplotMin = floor(Xplot[0]);
    double XplotMax = ceil(Xplot[nit-1]);
    YplotMin = floor(YplotMin);
    YplotMax = ceil(YplotMax);
    
    fprintf(stderr, "Xplot = [ %8.0f _ %8.0f ]\n", XplotMin, XplotMax);
    fprintf(stderr, "Yplot = [ %8.5f _ %8.5f ]\n", YplotMin, YplotMax);
    
    /* Intervals for ticks: */
    interval_t fr[2] = { (interval_t){{ XplotMin, XplotMax }}, (interval_t){{ YplotMin, YplotMax }} };
    double tstep[2] = { 100, 1 };
    double sstep[2] = { 10, 1 };
    
    /* Adjust intervals for graph clipping: */
    double Xpmrg = 0.01*(XplotMax-XplotMin);
    double Ypmrg = 0.01*(YplotMax-YplotMin);
    
    XplotMin = XplotMin - 0.01*Xpmrg;
    XplotMax = XplotMax + 0.01*Xpmrg;
    YplotMin = YplotMin - 0.01*Ypmrg;
    YplotMax = YplotMax + 0.01*Ypmrg;
    
    /* Plot scales {mag[0..1]} and function clip interval {gr[0..1]}: */
    double mag[2] = { xsz/(XplotMax - XplotMin), ysz/(YplotMax - YplotMin) };
    interval_t gr[2] = { (interval_t){{ XplotMin, XplotMax }}, (interval_t){{ YplotMin, YplotMax }} };
    
    /* Plot graph: */
    pswr_set_pen(*psP, 0,0,0, 0.15, 0,0);
    double xa = 0, ya = 0;
    for (it = 0; it < nit; it++)
      { double xb = xmin + (Xplot[it] - XplotMin)*mag[0];
        double yb = ymin + (Yplot[it] - YplotMin)*mag[1];
        if (it > 0) { pswr_segment(*psP, xa, ya, xb, yb); }
        xa = xb; ya = yb;
      }

    finish_figure(psP, o);
  }