PRO PLOT_CH,chfile=chfile,brfile=brfile,PS=ps,OUTFILE=outfile,POL_CH=pol_ch,_EXTRA=e

;+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

; set some basic parameters
!p.font = 1
DEVICE, SET_FONT='Times', /TT_FONT
!p.background = 255
!p.charsize=2
IF KEYWORD_SET(PS) THEN DEV,OPTION=4,FILENAME=outfile+'.ps'

;+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

; read in the r1.hdf files for displaying the CH boundaries
filename = chfile
ch_raw = READ_MAS_2D(filename=filename,coord1=xch,coord2=ych)
ch_raw = transpose(ch_raw)
lonch = xch/!DTOR
latch = ych/!DTOR - 90.
nlon = n_elements(ch_raw(*,0))
nlat = n_elements(ch_raw(0,*))
; Read in Br and slice first dim
filename =brfile
br3d = READ_MAS_3D(filename=filename,R=r,THETA=ybr,PHI=xbr)
brphoto_orig = (REFORM(br3d(0,*,*)) + REFORM(br3d(1,*,*)))/2.
brphoto_orig = transpose(brphoto_orig)
lonbr = xbr/!DTOR
latbr = ybr/!DTOR-90.

lonbr_int = INTERPOL(INDGEN(N_ELEMENTS(lonbr)),lonbr,lonch)
latbr_int = INTERPOL(INDGEN(N_ELEMENTS(latbr)),latbr,latch)
brphoto_int = INTERPOLATE(brphoto_orig, lonbr_int, latbr_int, MISSING=0,/GRID)
;brphoto_int = SMOOTH(brphoto_int,100)

;+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

ch2 = ch_raw
c_ch = 255
c_ch_pos = 250 ; 250
c_ch_neg = 50 ; 50
rad_cut =2.0
ch2(where(ch_raw lt rad_cut)) = c_ch
for ii=0,nlon-1 do begin
for jj=0,nlat-1 do begin
    if ((ch_raw(ii,jj) ge rad_cut) and (brphoto_int(ii,jj) lt 0)) then $
        ch2(ii,jj) = c_ch_neg
    if ((ch_raw(ii,jj) ge rad_cut) and (brphoto_int(ii,jj) gt 0)) then $
        ch2(ii,jj) = c_ch_pos
endfor
endfor

;+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

loadct,39
set_viewport,.15,.9,.2,.85
CONTOUR,ch2,lonch,latch,/fill,xstyle=1,ystyle=1,$
    xrange=[0,360],yrange=[-90,90],xticks=6,yticks=6, $
    levels=[c_ch_neg,c_ch_pos,c_ch],$
    c_color=[c_ch_neg,c_ch_pos,c_ch], $
    XTITLE='Longitude (!Uo!N)',  $
    YTITLE='Latitude (!Uo!N)', $
    TITLE=title, $
	color=0,CHARSIZE=3
loadct,0
loadct,39

CONTOUR,ch2,lonch,latch,/NOERASE,xstyle=5,ystyle=5,$
    xrange=[0,360],yrange=[-90,90],xticks=6,yticks=6, $
    levels=[c_ch_neg,c_ch_pos,c_ch],color=0,THICK=2
pol_ch = ch2

;+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
IF KEYWORD_SET(PS) THEN DEV,OPTION=7

RETURN

END