;
; PURPOSE: Overplot NL, FPs and/or FLs onto Image
; WRTTEN:  09Feb2006	Xuepu
; MODIFIED: Feb. 2008	Xuepu Zhao
;
pro pimgnlfpfl_phofbrssfiflfoflf08,phof,brssf,iflf,oflf,sy,r1,r2,$
                                   XTIC=xtic,OPFP=opfp,OPFL=opfl


    if n_params( ) lt 1 then begin
       print,'pimgnlfpfl_phofbrssfiflfoflf08,phof,brssf,iflf,oflf,sy,r1,r2,$'
       print,'                               xtic=,/opfp,/opfl'
       print,"sy='sla','lat' or 'nml'"
       return
    end
    outfile='pimgnlfpfl_phofbrssfiflfoflf08.ps'
    psout,012,outfile
    bla=readfits(phof,hd,/silent)
    bbla=BYTSCL(bla,max=200,min=-200)
    if keyword_set(xtic) then $
       scimg_bba,bbla,ttl,px,py,sy=sy,csz=1.2,xtic=xtic else $
       scimg_bba,bbla,ttl,px,py,sy=sy,csz=1.2
    if keyword_set(opfl) then opfl_oflf,oflf,r1,r2,sy,/clr
    if keyword_set(opfp) then begin
       opfp_iflf,iflf,0.4,sy,/clr
       opnl_brf,brssf,sy,px,py,/clr
    endif
    print,'Outfile: ',outfile
    psdone

end
;
; PURPOSE: produce an Carrington image that can be overlad by CONTOURs & lines
; MODIFIED: 24March2005 Xuepu Zhao for adding keyword sy
;
PRO scimg_bba,bba,tt,px,py,SY=sy,CSZ=csz,XTIC=xtic

    if N_params() LT 1 then begin
       print,'PRO pimg0,bba,tt,px,py,sy=,csz='
       print,"sy = 'lat','sla','nml', csz=1.0"
       return
    endif
    if not keyword_set(csz) then csz=1.0
    if not keyword_set(sy) then sy='nml'
    interp=1
    CASE 1 OF                       ;set window used by contour or plot
      sy eq 'lat': begin
                   plot,[[0,0],[1,1]],/nodata,xstyle=5,ystyle=5,title=tt,$
                        xrange=[0,360],yrange=[-90,90],charsize=csz
      ;            if keyword_set(xtic) then sc_annot,1,xtic=xtic,/lat $
      ;                                 else sc_annot,1,/lat
                   end
      sy eq 'sla': begin
                   plot,[[0,0],[1,1]],/nodata,xstyle=5,ystyle=5,title=tt,$
                        xrange=[0,360],yrange=[-1,1],charsize=csz
      ;            if keyword_set(xtic) then sc_annot,1,xtic=xtic else $
      ;                                      sc_annot,1
                   end
      sy eq 'nml': begin
                   plot,[[0,0],[1,1]],/nodata,xstyle=5,ystyle=5,title=tt,$
                        xrange=[0,1],yrange=[0,1],charsize=csz
      ;            if keyword_set(xtic) then sc_annot,1,xtic=xtic,/norm $
      ;                                 else sc_annot,1,/norm
                   end
    ENDCASE
    px = !x.window * !d.x_vsize     ;Get size of window in device units
    py = !y.window * !d.y_vsize
    swx = px(1)-px(0)               ;Size in x in device units
    swy = py(1)-py(0)               ;Size in Y
    sz=size(bba)
    six = float(sz(1))              ;Image sizes
    siy = float(sz(2))
    if (!d.flags and 1) ne 0 then begin           ;Scalable pixels
       tv,bba,px(0),py(0),xsize = swx, ysize = swy, /device
    endif else begin                              ;Not scalable pixels
       bimgn = poly_2d((bba),[[0,0],[six/swx,0]], $
             [[0,siy/swy],[0,0]],interp,swx,swy)
       tv,bimgn,px(0),py(0)
    endelse
    CASE 1 OF                       ;set window used by contour or plot
      sy eq 'lat': begin
                   if keyword_set(xtic) then sc_annot,1,xtic=xtic,/lat $
                                        else sc_annot,1,/lat
                   end
      sy eq 'sla': begin
                   if keyword_set(xtic) then sc_annot,1,xtic=xtic else $
                                             sc_annot,1
                   end
      sy eq 'nml': begin
                   if keyword_set(xtic) then sc_annot,1,xtic=xtic,/norm $
                                        else sc_annot,1,/norm
                   end
    ENDCASE
    return

END
;
; WRITTEN:  30JAN2006   Xuepu
;
pro opnl_brf,brf,sy,px,py,CLR=clr

    if N_params() lt 1 then begin
       print,'pnl_brf,brf,sy,/clr'
       return
    endif
;
    bra=readfits(brf,hd,/silent) 
    sz=SIZE(bra)
    zgrid,sz(1),sz(2),phd,thd,lad
    CASE 1 OF
     sy eq 'sla': begin
       lad=SIN(lad*!DTOR)
       yr=[-1,1]
                  end
     sy eq 'lat': yr=[-90,90]
     else: 
    ENDCASE
    xy_pt,phd,lad,xx,yy
    lv=0.025
    if keyword_set(clr) then begin
       tek_color
       contour,bra,xx,yy,position=[px(0),py(0),px(1),py(1)],xstyle=5,ystyle=5,$
               level=[-lv,0,lv],c_colors=[2,0,3],c_thick=[1,5,1],/device,/noerase
    endif else begin
       contour,bra,xx,yy,position=[px(0),py(0),px(1),py(1)],xstyle=5,ystyle=5,$
               level=[-lv,0,lv],c_lines=[1,0,2],c_thick=[1,5,1],/device,/noerase
    endelse

end
;
; WRITTEN:  30JAN2006	Xuepu
;
pro opfp_iflf,iflf,psz,sy,CLR=clr

    phfp=fltarr(2)
    thfp=phfp
    get_lun,unit
    openr,unit,iflf
    WHILE ~ EOF(unit) DO BEGIN 
      readf,unit,np,fte  
      anp=ABS(np)
      for n=0,anp-1 do begin
          readf,unit,x0,y0,z0,r0,t0,p0                 
          if n eq anp-1 then begin
             phfp(0)=p0                    
             thfp(0)=90-t0
             CASE 1 OF
               sy eq 'sla': thfp(0)=SIN(thfp(0)*!DTOR)
               sy eq 'nml': begin
                  thfp(0)=SIN(thfp(0)*!DTOR)
                  thfp(0)=(thfp(0)+1)/2.0
                  phfp(0)=phfp(0)/360.0
                            end
               else:
             ENDCASE
          endif
      endfor
      if keyword_set(clr) then begin
         tek_color
         mkpsym8,psz
      ;  if np gt 0.0 then begin
      ;     oplot,phfp(0:0),thfp(0:0),psym=8,color=4
      ;  endif else begin
      ;     oplot,phfp(0:0),thfp(0:0),psym=8,color=2
      ;  endelse
         if np gt 0.0 then oplot,phfp(0:0),thfp(0:0),psym=8,color=3
         if np lt 0.0 then oplot,phfp(0:0),thfp(0:0),psym=8,color=2
      endif else begin
         if np ge 0.0 then begin
            mkpsym8,psz
            oplot,phfp(0:0),thfp(0:0),psym=8
         endif else begin
            mkpsym8,psz,/open
            oplot,phfp(0:0),thfp(0:0),psym=8
         endelse
      endelse
    ENDWHILE
    free_lun,unit

end
;
; PURPOSE:      Overplot foot-points of IFLs at rr
; NOTE:         bottom=1.001, dstep=0.02
; WRITTEN:      16Aug2000    refer to ~/ZIDLLIB/JOB/opfp_brfiflf.pro;
pro opfl_oflf,oflf,r1,r2,sy,CLR=clr

    if keyword_set(clr) then tek_color
    mkpsym8,0.2
    get_lun,unit
    openr,unit,oflf     
    while eof(unit) ne 1 do begin
       readf,unit,np        
       anp=abs(np)
       phfl=fltarr(anp)
       lafl=phfl & rrfl=phfl
       for n=0,anp-1 do begin            
           readf,unit,x0,y0,z0,r0,t0,p0
           phfl(n)=p0
           lafl(n)=90-t0
           rrfl(n)=r0
       endfor
       if sy eq 'sla' then lafl=SIN(lafl*!DTOR)
       if anp ge 3 then begin
          rmax=max(rrfl)
          if rmax gt r1 and rmax lt r2 then begin
             ind=where(rrfl eq rmax)
             in_rmax=ind(0)
             if in_rmax lt (anp-1) then begin
                if np gt 0 then begin
                   for k=0,in_rmax-1 do oplot,phfl(k:k),lafl(k:k),psym=8,color=3
                   for k=in_rmax,anp-1 do oplot,phfl(k:k),lafl(k:k),psym=8,color=2
                endif else begin
                   for k=0,in_rmax-1 do oplot,phfl(k:k),lafl(k:k),psym=8,color=2
                   for k=in_rmax,anp-1 do oplot,phfl(k:k),lafl(k:k),psym=8,color=3
                endelse
             endif else begin
                if np gt 0 then begin
                   for k=0,anp-1 do oplot,phfl(k:k),lafl(k:k),psym=8,color=3 
                endif else begin
                   for k=0,anp-1 do oplot,phfl(k:k),lafl(k:k),psym=8,color=2 
                endelse
             endelse
          endif    
       endif
    endwhile
    free_lun,unit

end