;+
;
; Project   : s3drs reconstruction software
;                   
; Name      : p_standardplots
;               
; Purpose   : part of Pixon, part of the s3drs GUI
;               
; Explanation: part of the Pixon wrappers used to access and
; manipulate the various arcane Pixon structures and objects.
;               
; our standard output format
; Usually just needs p_standardplots,pxwrapper,oPixon
; but you can alternate say p_standardplots,filein=filein
; to have it ingest an IDL save file that contains pxwrapper & oPixon
;
; Use       : typically called from 's3drs.pro' 
;    
; Inputs    : (various)
;               
; Outputs   : (various)
;
; Keywords  : (various)
;               
; Common    : none
;               
; Restrictions: none
;               
; Side effects: none
;               
; Category    : reconstruction, pixon
;               
; Prev. Hist. : None.
;
; Written     : Sandy Antunes, NRC, March-April 2005-2009
;               
;-            
; our standard output format
; Usually just needs p_standardplots,pxwrapper,oPixon
; but you can alternate say p_standardplots,filein=filein
; to have it ingest an IDL save file that contains pxwrapper & oPixon
;
; You can also pass it a source image 'synthetic=image' for problems
; where the original source is known, and it will include it in the plots.
;

; now turned off
;function data2photonslabel,data,oom,istd
;  return,''
;  dmax=max(data)
;  dmin=dmax/(10.0^oom)
;  dmax=dmax*istd
;  dmin=dmin*istd
;  sdmax=strtrim(string(dmax),2)
;  sdmin=strtrim(string(dmin),2)
;  sdmax=strmid(sdmax,0,strpos(sdmax,'.')+3) + strmid(sdmax,strpos(sdmax,'e'))
;;;     sdmin=strmid(sdmin,0,strpos(sdmin,'.')+3) + $
;;;       strmid(sdmin,strpos(sdmin,'e'))
;  title='['+sdmin+','+sdmax+'] photons sec!U-1!N cm!U-2!N sr!U-1!N!C'
;  return,title
;end


; note extra wrapper for 'generate postscript for all Pixon
; save files in the given directory'
PRO p_standardplots_all,best=best,shellbest=shellbest,fluxbest=fluxbest

  pixonlist=file_search('.',"*solved.sav",count=icount)

  for ifile=0,icount-1 do begin
      p_standardplots,filein=pixonlist[ifile],/noscreen,/quiet,$
        best=best,shellbest=shellbest,fluxbest=fluxbest
  end

END


PRO p_sptv,plotset,id,tag,header,data,saved,image=image,syn=syn,$
           hardlimits=hardlimits

;     ybuffer=0;256*8

  image=KEYWORD_SET(image)

  Nd=n_elements(data[0,0,*])

  blankrow=0
  limits=[!VALUES.F_NAN,!VALUES.F_NAN]; placeholder
  if (image) then begin
      Nl=3
      lset=['output X','output Y','output Z','input X','input Y','input Z']
      threeview=1
      if (n_elements(hardlimits) ne 0) then limits=hardlimits[2:3]
  end else begin
      Nl=plotset.Nd
      lset=plotset.lset
      ; also handle 'singlets' which are not 'input vs output'
      if (Nd eq plotset.Nd) then begin
          lset=lset[(plotset.Nd*2):(plotset.Nd*3)-1]
          blankrow=1
      end
      threeview=0
      if (n_elements(hardlimits) ne 0) then limits=hardlimits[0:1]
  end

  fname1=plotset.tag+'_'+strtrim(string(id),2)+'a_'+tag+'_lin'
  fname2=plotset.tag+'_'+strtrim(string(id),2)+'b_'+tag+'_log'

  hA=strtrim(string(id),2)+'a) '+header
  hB=strtrim(string(id),2)+'b) log('+header+')'

  if (plotset.noscreen eq 0) then begin
;      print,tag
;      tv_multi,quiet=plotset.quiet,data[*,*,0:plotset.Nd-1],$

      if (plotset.gif) then begin
          ; only 'slideshow' if gifs
          if (threeview) then winid=plotset.winid3 else winid=plotset.winid
      end

      tv_multi,quiet=plotset.quiet,data,$
        labels=lset,line=Nl,$
        fixed=256,title=plotset.title, header=hA,caption=plotset.desc,$
        oom=plotset.oom,/show,autocap=plotset.cap,threeview=threeview,$
        blankrow=blankrow,minval=limits[0],maxval=limits[1],$
        inverse=plotset.inverse,gif=plotset.gif,filename=fname1+'.gif',$
        winid=winid,slideshow=plotset.gif

      if (plotset.gif) then begin
          ; only 'slideshow' if gifs
          if (threeview) then plotset.winid3=winid else plotset.winid=winid
      end else begin
          winid=-1
      end

;      tv_multi,quiet=plotset.quiet,data[*,*,0:plotset.Nd-1],$
      tv_multi,quiet=plotset.quiet,data,$
        labels=lset,line=Nl,$
        fixed=256,title=plotset.title, header=hB,caption=plotset.desc,$
        oom=plotset.oom,/log,/show,autocap=plotset.cap,threeview=threeview,$
        blankrow=blankrow,minval=limits[0],maxval=limits[1],$
        inverse=plotset.inverse,gif=plotset.gif,filename=fname2+'.gif',$
        winid=winid,slideshow=plotset.gif

  end

;  tv_multi,quiet=plotset.quiet,data[*,*,0:plotset.Nd-1],labels=lset,line=Nl,$
  tv_multi,quiet=plotset.quiet,data,labels=lset,line=Nl,$
    fixed=256,title=plotset.title, header=hA,caption=plotset.desc,$
    oom=plotset.oom,/border,/ps,/landscape,filename=fname1+'.ps',$
    autocap=plotset.cap,threeview=threeview,blankrow=blankrow,$
    minval=limits[0],maxval=limits[1],inverse=plotset.inverse,$
    nocloseps=plotset.mega,appendps=plotset.mega

  
;  tv_multi,quiet=plotset.quiet,data[*,*,0:plotset.Nd-1],labels=lset,line=Nl,$
  tv_multi,quiet=plotset.quiet,data,labels=lset,line=Nl,$
    fixed=256,title=plotset.title, header=hB,caption=plotset.desc,$
    oom=plotset.oom,/log,/border,/ps,/landscape,filename=fname2+'.ps',$
    autocap=plotset.cap,threeview=threeview,blankrow=blankrow,$
    minval=limits[0],maxval=limits[1],inverse=plotset.inverse,$
    nocloseps=plotset.mega,appendps=plotset.mega

  saved=[saved,fname1,fname2]

END

PRO p_sp_lset,pxwrapper,lset,captions,runname,saved,Nd=Nd,Nstart=Nstart,$
              noises=noises

   if (n_elements(Nstart) eq 0) then Nstart=0
   if (n_elements(Nd) eq 0) then Nd=pxwrapper.sr.Nd

     ; now do automatic label making
   lset=strarr(Nd*3)
   for i=0,Nd-1 do lset[i]='Output '
   for i=Nd,(Nd*2)-1 do lset[i]='Input '
   for i=(Nd*2),(Nd*3)-1 do lset[i]=''

   front = strarr(Nd*3)
   rhostr = strarr(Nd)
   phistr = strarr(Nd)
   for i=0,Nd-1 do rhostr[i] = $
     strtrim(string(float(pxwrapper.sr.rho[i+Nstart])),2)
   for i=0,Nd-1 do rhostr[i]=strmid(rhostr[i],0,strpos(rhostr[i],'.')+2)
   for i=0,Nd-1 do phistr[i] = $
     strtrim(string(float(pxwrapper.sr.phi[i+Nstart])),2)
   for i=0,Nd-1 do phistr[i]=strmid(phistr[i],0,strpos(phistr[i],'.')+2)

   satF=strmid(pxwrapper.sr.sat_name,0,1)
   satB=strmid(pxwrapper.sr.sat_name,0,1,/reverse)
   sat_name=satB; assumes STEREO_A or STEREO_B then checks for others
   for i=0,Nd-1 do $
     if (sat_name[i] eq 'O') then sat_name[i]='L' else $
     if (sat_name[i] ne 'A' and sat_name[i] ne 'B') then sat_name[i]=satB[i]
   
   for i=0,Nd-1 do front[i] = pxwrapper.sr.detector[i+Nstart] + $
     sat_name[i+Nstart] + ' (' + rhostr[i] + '!U0!N,' + phistr[i] + '!U0!N)'
   
   ; if polar angles are not all 0 and thus relevant
   if (total(pxwrapper.sr.fan) ne 0) then $
     for i=0,Nd-1 do front[i]=front[i] + ' polar ' + $
     strtrim(string(floor(pxwrapper.sr.fan[i+Nstart])),2) + '!Uo!N'
                                ; now duplicate for set 2
   for i=Nd,(Nd*2)-1 do front[i]=front[i-Nd]
   for i=(Nd*2),(Nd*3)-1 do front[i]=front[i-Nd]
   
   for i=0,(Nd*3)-1 do lset[i]=lset[i] + front[i]

   mydate=pxwrapper.sr.date_obs[0]
   adate=strsplit(mydate,' ',/extract)
   if (n_elements(adate) gt 3) then $
     newdate=adate[4] + '-' + adate[1] + '-' + adate[2] + 'T' + adate[3] $
   else $
     newdate=adate[0] + ' ' + (strsplit(adate[1],':',/extract))[0]
   
;     los_utc=anytim2utc(newdate); needs year month day hour
   
   runname2=strjoin(strsplit(runname,' ',/extract),'_')
   
   fname=runname2+'.txt'
   p_report,pxwrapper,/nosave,returntext=captionarr,/terse,runname=runname,$
     noises=noises
;     p_report,pxwrapper,/nosave,filename=fname,runname=runname,$
;       returntext=captionlongarr

   p_report,pxwrapper,/nosave,runname=runname,returntext=captionlongarr,$
     noises=noises
   saved=[fname]

   captionarr[2]='   '+captionarr[2]

   caption=STRJOIN(captionarr,'!C')
   caption='!C!C'+caption ; leading line allows for 2-line 'title' below
   
   captionlong=STRJOIN(captionlongarr,'!C')
   
   rsun=(1.0/pxwrapper.sr.L0tau)
   rsun=[rsun,rsun]
   
;     l1=data2photonslabel(odataB,1.0,pxwrapper.sr.istd[0+Nstart])
;     l4=data2photonslabel(odataB,4.0,pxwrapper.sr.istd[0+Nstart])
   
;     l1title=  date + ' cor2 data'
;     l4title=  date + ' cor2 data'
;     lrtitle=  date + ' fractional error'
;     lstitle=  date + ' noise'
;     lsntitle=  date + ' S/N'
;     lwttitle=  date + ' weighting'
;     lmasktitle=  date + ' masks'

   captions=strarr(2)
   captions[0]=caption
   captions[1]=captionlong
   
   ;print,'caption: ',caption
   ;print,'captionlong: ',captionlong
   ;print,'sc.note: ',pxwrapper.sc.note

END

PRO p_sp_cover,plotset,pxwrapper,saved,id=id

  if (n_elements(id) eq 0) then id=0
;  runname2=strjoin(strsplit(runname,' ',/extract),'_')

  lshort=pxwrapper.sr.sat_name

;  fname='0_'+runname2+'_orbit.ps'
  if (plotset.mega) then begin
      fname=strjoin(strsplit(plotset.title,' ',/extract),'_') ;+'_all.ps'
  end else begin
      fname=plotset.tag+'_'+strtrim(string(id),2)+'_cover_orbit'
  end
  if (plotset.noscreen eq 0) then $
    sat_los_view,pxwrapper.sr.rho,pxwrapper.sr.LL,caption=plotset.banner,$
    labelset=lshort,phiset=pxwrapper.sr.phi,/fpa,/newwindow,$
    filename=fname+'.gif',writegif=plotset.gif,inverse=plotset.gif

  sat_los_view,pxwrapper.sr.rho,pxwrapper.sr.LL,$
    /ps,filename=fname+'.ps',labelset=lshort,/fpa,$
    caption=plotset.banner,phiset=pxwrapper.sr.phi,nocloseps=plotset.mega
  saved=[saved,fname]
        
END


FUNCTION p_genflyby,img,tag,noisy=noisy

   ; also create a 360-view flyby plot
   flybysave=tag+'_flyby'
   flyby,img,/grandtour,/cap,/addaxis,/mask,savefile=flybysave,noisy=noisy
   print_flyby,flybysave,/noscreen,skip=2,nper=4; make bigger
   return,flybysave+'.ps'

END
    
    
; optionally can set 'imgmax' and 'imgmin' arrays for threeview plots,
; also use keyword /showmaxmin to show the current values

; suggested hardlimits=[5e3,5e5,1e-13,1e-11]
;
; /megafile prints all postscript to one file, rather than individual ones.
;
; /gif dumps the windows out as gifs (conflicts with /mega and /ps)


PRO p_standardplots,pxwrapper,oPixon,runname,noisy=noisy,quiet=quiet,$
                    masks=masks,date=date,noscreen=noscreen,oom=oom,$
                    inputonly=inputonly,filein=filein,hardlimits=hardlimits,$
                    quickcheck=quickcheck,cap=cap,midpoint=midpoint,$
                    synthetic=synthetic,fullset=fullset,extras=extras,$
                    slice=slice,best=best,Nend=Nend,Nstart=Nstart,$
                    flyby=flyby,shellbest=shellbest,fluxbest=fluxbest,$
                    imgmax=imgmax,imgmin=imgmin,proof=proof,mxit=mxit,$
                    datamax=datamax,datamin=datamin,megafile=megafile,$
                    nocut=nocut,cor1=cor1,inverse=inverse,addtag=addtag,$
                    multips=multips,gif=gif,c3=c3,noises=noises

    ; multips 0 = usual, 1 = do cover but no image and leave open,
    ;         2 = skip cover and close when done
    if (n_elements(multips) eq 0) then multips=0
    megafile=KEYWORD_SET(megafile)
    if (megafile) then begin
        noscreen=1;
        ps=1;
        quiet=1;
        proof=0;
        if (multips eq 2) then current_device='X' else current_device=!d.name
;        set_plot,'PS'
;        device,filename='mega.ps',/color,bits_per_pixel=8,/landscape,$
;          font_size=9,/inches
    end

    heap_gc; we allocate a lot of temporary data items, so free things first
    
    if (n_elements(filein) ne 0) then begin
        restore,filein
        ; if file has a standardly-named input model 'timage', use it
        if (n_elements(synthetic) eq 0 and n_elements(timage) ne 0) then $
          synthetic=timage
    end

    pixon_env,pxwrapper.object ; ensures proper pixon env is loaded, if needed
    
; *** PARSE INPUT OPTIONS
    
    best=KEYWORD_SET(best)
    cor1=KEYWORD_SET(cor1)
    c3=KEYWORD_SET(c3)
    if (best) then begin
        if (n_elements(cap) eq 0) then cap=1
        if (n_elements(oom) eq 0) then oom=2
        if (n_elements(datamax) eq 0) then datamax=5E-11
        if (n_elements(imgmax) eq 0) then imgmax=1E6
; removed instrument-specific limits, they got confusing to plot-readers
;        if (n_elements(datamax) eq 0) then $
;          if (cor1) then datamax=1E-10 else $
;          if (c3) then datamax=5E-11 else datamax=1E-8
;        if (n_elements(imgmax) eq 0) then $
;          if (cor1) then imgmax=1E6 else $
;          if (c3) then imgmax=1E7 else imgmax=5E7
    end
;    print,'imagemax ',imgmax,' and datamax ',datamax

    shellbest=KEYWORD_SET(shellbest)
    if (shellbest) then begin
        if (n_elements(cap) eq 0) then cap=1
        if (n_elements(oom) eq 0) then oom=1
        if (n_elements(datamax) eq 0) then datamax=1E-10
        if (n_elements(imgmax) eq 0) then imgmax=5E6
    end
    fluxbest=KEYWORD_SET(fluxbest)
    if (fluxbest) then begin
        if (n_elements(cap) eq 0) then cap=1
        if (n_elements(oom) eq 0) then oom=2
        if (n_elements(datamax) eq 0) then datamax=5E-11
        if (n_elements(imgmax) eq 0) then imgmax=1e6
    end

    proof=KEYWORD_SET(proof)
    zbuffer=KEYWORD_SET(zbuffer)
    ps=KEYWORD_SET(ps)
    gif=KEYWORD_SET(gif)
    inverse=KEYWORD_SET(inverse)
    nocut=KEYWORD_SET(nocut)
    noisy=KEYWORD_SET(noisy); default is now to be quiet unless told otherwise
    if (n_elements(quiet) ne 0) then quiet=KEYWORD_SET(quiet) else quiet=~ noisy
    flyby=KEYWORD_SET(flyby)
    cap=KEYWORD_SET(cap) ; whether to remove bright points before plotting or not
    extras=KEYWORD_SET(extras)
    if (n_elements(quickcheck) eq 0) then quickcheck=0
    showmaxmin=KEYWORD_SET(showmaxmin)
    noscreen=KEYWORD_SET(noscreen) ; only makes PS, no on-screen plotting
    inputonly=KEYWORD_SET(inputonly) ; ignore output
    fullset=KEYWORD_SET(fullset) ; give everything, not just input and result
    slice=KEYWORD_SET(slice);
    if (n_elements(oom) eq 0) then oom=-1 ; no oom if not given
    if (n_elements(date) eq 0) then date=''
    if (n_elements(runname) eq 0) then runname=pxwrapper.sc.gname
    if (n_elements(addtag) ne 0) then runname += addtag
    if (n_elements(oom) eq 0) then oom=4
    if (n_elements(hardlimits) eq 0) then $
      hardlimits=[!VALUES.F_NAN,!VALUES.F_NAN,!VALUES.F_NAN,!VALUES.F_NAN]

    p_update_runstats,pxwrapper,oPixon,/alt
    
    if (n_elements(Nend) eq 0) then Nend=pxwrapper.sr.Nd
    if (n_elements(Nstart) eq 0) then Nstart=0
    Nd=Nend-Nstart

    if (n_elements(datamax) ne 0) then $
      pxwrapper.sc.note += '!Cdatamax='+int2str(datamax*1E11)+'E-11'
    if (n_elements(imgmax) ne 0) then $
      pxwrapper.sc.note += ' , imagemax='+int2str(imgmax)
    if (tag_exist(pxwrapper.sc,'datacut')) then begin
        if (pxwrapper.sc.datacut ne 0) then pxwrapper.sc.note += $
          ' , datacut='+strtrim(string(pxwrapper.sc.datacut),2)
    end

; *** EXTRACT DATA
    
    indata=stereo_fetch(pxwrapper,oPixon,/data)
    indata=indata[*,*,Nstart:Nend-1]; optional cutting back
    render=stereo_fetch(pxwrapper,oPixon,/render)
    render=render[*,*,Nstart:Nend-1]; optional cutting back
    s2n=stereo_fetch(pxwrapper,oPixon,/sn)
    s2n=s2n[*,*,Nstart:Nend-1]; optional cutting back

    if (n_elements(masks) eq 0) then masks=wt2mask(oPixon) else $
      masks=bytarr(pxwrapper.sr.Ncube,pxwrapper.sr.Ncube,pxwrapper.sr.Nd)+1

    img=stereo_fetch(pxwrapper,oPixon,/image)
    m=p_getmask(pxwrapper.sr.detector[0],pxwrapper.sr.ncube,rpixels=rpixels)
    mask3D = cube_sphere(pxwrapper.sr.ncube,rpixels[0],r_outer=rpixels[1])
    if (nocut eq 0) then img=img*mask3D

    threeview,img,threexyz,/noplot,/pixon,xyangles=pxwrapper.sr.rho


; *** GENERATE LABELS
    
    tag=strmid(runname,0,21)

; cursory raw plot for quick sanity checks


    if (proof) then begin


        if (pxwrapper.sr.nd gt 2) then iele=3 else iele=2
        ; proof sheet is four lines: input_data,output_data,sigma,imgcube

        l1=indata[*,*,0:iele-1]*masks[*,*,0:iele-1]
        if (n_elements(datamax) eq 0) then datamax=max(l1)

        l2=render[*,*,0:iele-1]*masks[*,*,0:iele-1]

        qmap=stereo_fetch(pxwrapper,oPixon,/qmap)
        l3=qmap[*,*,0:iele-1]*masks[*,*,0:iele-1]
        l3= (l3 / max(l3)) * datamax

        if (iele eq 2) then begin
            blank=fltarr(pxwrapper.sr.ncube,pxwrapper.sr.ncube)
            l1=[[[l1]],[[blank]]]
            l2=[[[l2]],[[blank]]]
            l3=[[[l3]],[[blank]]]
        end

        if (n_elements(imgmax) eq 0) then imgmax=max(threexyz)
        threexyz=threexyz * datamax/imgmax

        rset=['Image XYZ','QMap','Render','Input data']
        if (n_elements(oom) eq 0) then oom=2
        if (noscreen eq 0) then $
          tv_multi,[[[threexyz]],[[l3]],[[l2]],[[l1]]],oom=oom,line=3,$
          /log,title='Proof sheet '+pxwrapper.sc.gname,rowtitle=rset,$
          /autocap,max=datamax

        tv_multi,[[[threexyz]],[[l3]],[[l2]],[[l1]]],oom=oom,line=3,$
          /log,title='Proof sheet '+pxwrapper.sc.gname,rowtitle=rset,/ps,$
          file='proof_'+runname+'.ps',/autocap,max=datamax

        return

    end


; more processing

    p_sp_lset,pxwrapper,lset,captions,runname,saved,Nd=Nd,Nstart=Nstart,$
      noises=noises
    
                                ; save standard plotables here for clarity
    plotset={desc:captions[0],quiet:quiet,Nd:Nd,oom:oom,lset:lset,ybuffer:0,$
             noscreen:noscreen,title:runname,banner:captions[1],cap:cap,$
             tag:tag,inverse:inverse,mega:megafile,gif:gif,winid:-1,winid3:-1}


    mm=masks[*,*,0]
;    if (n_elements(synthetic) ne 0) then m0=mm else 
    m0=mm*0.0+1.0

    masksix=[[[mm]],[[mm]],[[mm]],[[m0]],[[m0]],[[m0]]]
    
;    if (cap) then begin
;        sm1=smooth(render*masks,[8,8,0])
;        sm2=smooth(indata*masks,[8,8,0])
;        datamax=max([max(sm1),max(sm2)])
;    end
    
    
    mset=[[[render]],[[indata]]]
    
    if (n_elements(datamin) ne 0) then begin
;        mset=mset>(datamin)
        hardlimits[0]=datamin
    end
    if (n_elements(datamax) ne 0) then begin
;        mset=mset<(datamax)
        hardlimits[1]=datamax
    end

    if (inputonly) then begin
        mset=indata
        lset=lset[Nd:(Nd*2)-1]
    end

    if (n_elements(imgmin) ne 0) then begin
        ;imgxyz=imgxyz>imgmin
        hardlimits[2]=imgmin
    end else begin
        imgmin=min(threexyz)
    end

    if (n_elements(imgmax) ne 0) then begin
        ;imgxyz=imgxyz<imgmax
        hardlimits[3]=imgmax
    end else begin
        imgmax=max(threexyz)
    end
        
    if (n_elements(synthetic) ne 0) then begin
        threeview,synthetic,synthree,/noplot,/pixon,xyangles=pxwrapper.sr.rho
    end else begin
        if (inverse) then synthree=threexyz*0.0+imgmin else $
          synthree=threexyz*0.0+imgmax
    end

    imgxyz=[[[threexyz]],[[synthree]]]


; additional items


    if (quickcheck eq 0) then begin
        ; for speed, only generate these for full plotting calls

        sigma=stereo_fetch(pxwrapper,oPixon,/sigma)
        sigma=sigma[*,*,Nstart:Nend-1] ; optional cutting back

        qmap=stereo_fetch(pxwrapper,oPixon,/qmap)
        qmap=qmap[*,*,Nstart:Nend-1] ; optional cutting back

        fov=pr_fov(pxwrapper)
        threeview,fov,fovthree,/noplot,/pixon,xyangles=pxwrapper.sr.rho
    end
    
    if (fullset) then begin
        ; optional data items

        pxnmap=oPixon->get(/pxnidx,/nop)

        wt=stereo_fetch(pxwrapper,oPixon,/wt)
        if (size(wt,/n_dim) eq 1 and datatype(wt[0]) ne 'PTR') then begin
                                ; special handling for scalar masks
            twt=indata*0.0
            for i=0,n_elements(wt)-1 do twt[*,*,i]=twt[*,*,i] + wt[i]
            wt=twt
        end
    end

; *** IF NEEDED, DO 'QUICKCHECK' ITEMS ONLY

    if (quickcheck) then begin

;        p_sptv,plotset,0,'inout','Input (masked) & Output (unmasked) data',$
;          mset*[[[masks]],[[masks*0.0]]],saved,hardlimits=hardlimits
        tv_multi,mset
        p_sptv,plotset,0,'inout','Input & Output data, unmasked',$
          mset,saved,hardlimits=hardlimits

        p_sptv,plotset,0,'signoi','Signal/Noise ratio > 1 < 100',$
          (s2n > 1.0 < 100.0),saved

;        p_sptv,plotset,0,'image','Image Cubes (masked)',$
;          imgxyz*masksix,saved,/image,hardlimits=hardlimits
        p_sptv,plotset,0,'image','Image Cubes (unmasked)',$
          imgxyz,saved,/image,hardlimits=hardlimits

;        saved=[saved,p_genflyby(img,plotset.tag)]

        return
    end
    
; *** MAKE COVER SHEET

    id=0
    if (multips ne 2) then p_sp_cover,id=id++,plotset,pxwrapper,saved

; *** AND NOW PLOT DATA ITEMS
    
;    p_sptv,plotset,id++,'inout','Input (masked) & Output (unmasked) data',$
;      mset*[[[masks*0.0+1.0]],[[masks]]],saved,hardlimits=hardlimits
    p_sptv,plotset,id++,'inout','Input (masked) & Output (masked) data',$
      mset*[[[masks]],[[masks]]],saved,hardlimits=hardlimits

    p_sptv,plotset,id++,'signoi','Signal/Noise ratio > 1 < 100',$
      (s2n > 1.0 < 100.0),saved,hardlimits=[1,100]
    


; *** DONE DATA ITEMS, NOW DO GLOBALS AND IMAGE
    p_sptv,plotset,id++,'qmap','abs[Input-Output]/Sigma',qmap,saved,$
      hardlimits=[0.01,10]
    
    if (multips ne 1) then $
      p_qplots,id=id++,oPixon,noscreen=plotset.noscreen,tag=plotset.tag,saved,$
      appendps=plotset.mega,mxit=mxit

;    p_sptv,plotset,id++,'image','Image Cubes (masked)',$
;      imgxyz*masksix,saved,/image,hardlimits=hardlimits
    if (multips ne 1) then $
      p_sptv,plotset,id++,'image','Image Cubes (unmasked)',$
      imgxyz,saved,/image,hardlimits=hardlimits

    if (slice) then $
      p_sliceplots,id=id++,pxwrapper,oPixon,synthetic,midpoint=midpoint,$
      tag=plotset.tag,noscreen=plotset.noscreen,/bracket,saved=saved,$
      appendps=plotset.mega

    if (fullset) then begin

        p_sptv,plotset,id++,'Input & Output, no mask','inout',mset,saved

        p_sptv,plotset,id++,'Weight','wt',wt,saved
    
        p_sptv,plotset,id++,'Masks','masks',masks,saved
    
        p_sptv,plotset,id++,'Input & Output Cubes','imagemasked',$
          imgxyz,saved,/image
    
        p_sptv,plotset,id++,'Pixon map','pixonmap',pxnmap,saved
    
        p_sptv,plotset,id++,'FOV Coverage','fov',fovthree,saved,/image

    end

    if (flyby) then saved=[saved,p_genflyby(img,plotset.tag,noisy=noisy)]

    if (megafile eq 1 and multips ne 1) then begin

        ; close out our single massive postscript file
        device,/CLOSE_FILE
        set_plot,current_device

    end else if (multips eq 0) then begin

        if (noisy) then print,'Wrote these files:'
        for i=0,n_elements(saved)-1 do if (file_exist(saved[i])) then $
          if (noisy) then print,saved[i]
        print,'lpr '+strjoin(saved[1:*],' ') ; to help user print all ps files

    end


    ; now cleanup
    ;if (n_elements(img) ne 0) then delvar,img
    ;if (n_elements(mask3D) ne 0) then delvar,mask3D
    ;if (n_elements(sigma) ne 0) then delvar,sigma
    ;if (n_elements(qmap) ne 0) then delvar,qmap
    ;if (n_elements(wt) ne 0) then delvar,wt
    ;delvar,indata,render,s2n


END

;     header='FIGURE 5: ' + runname + ' Q History'
;     fname='5_'+runname2+'_q.ps'
;     if (pxwrapper.sc.iterations ne 0) then begin
;     if (noscreen eq 0 and subset eq 0) then $
;       p_plotstats,pxwrapper,solo='q',header=header,/log
;     if (subset eq 0) then $
;       p_plotstats,pxwrapper,solo='q',filename=fname,/ps,$
;       header=header,/log
;     if (subset eq 0) then saved=[saved,fname]

