mydate = systime()

dispipB=0 & dispipB = bytarr(512,512) 
dispiB=0  & dispiB = bytarr(512,512) 
dispiM=0  & dispiM = bytarr(512,512) 
ymn=10.   & ymx=10000.
Bcol=255  &  pBcol=241  & Mcol=242


if( n_elements(initR) eq 0) then begin ;{
    initR=1
    device,window_state=winst
    if(winst(3) eq 0)then window,3,xsize=600,ysize=350,xpos=0,ypos=30
    if(winst(0) eq 0)then window,0,xsize=512*3,ysize=512

    pBcomp=intarr(512,512)
    filenm = 'pB.comp.ROT'          & print,'Trying to open: ',filenm
    close,1 & openr,1,filenm & readu,1,pBcomp & close,1

    Bcomp =intarr(512,512)
    filenm = 'B.comp.ROT'           & print,'Trying to open: ',filenm
    close,1 & openr,1,filenm & readu,1,Bcomp  & close,1

    mk3   =intarr(512,512)
    filenm = '91d192.17_28.ch0.rpB' & print,'Trying to open: ',filenm
    close,1 & openr,1,filenm & readu,1,mk3  & close,1

    dolog='y'
    printdata='y'

    _a=intarr(256)  &  _b=intarr(256)  &  _c=intarr(256)
    openr,11,'difcol.rgb' & readu,11,_a & readu,11,_b & readu,11,_c & close,11
	_a(  1)=  1 & _b(  1)=  1 & _c(  1)=255
	_a( 15)=  1 & _b( 15)=255 & _c( 15)=255
	_a( 16)=255 & _b( 16)=  1 & _c( 16)=255
	_a(220)=255 & _b(220)=201 & _c(220)=  1
	_a(240)=225 & _b(240)=225 & _c(240)= 33
	_a(241)=255 & _b(241)=  1 & _c(241)=  1
	_a(242)=  1 & _b(242)=255 & _c(242)=  1
	_a(243)=225 & _b(243)=225 & _c(243)= 33
	tvlct,_a,_b,_c

    xcen = 255.5  & ycen = 249.5   & pixrad = 104.44  & rnder = 0.505
    acirc = !PI * 2.0 / 360.0
    dpx = intarr(156)
    dpy = intarr(156)

    xcenM = 252.5 & ycenM = 272.5
    dpxM = intarr(156)
    dpyM = intarr(156)

    fakex   = (findgen(156) + pixrad)/pixrad

    apltrad = (findgen(156) + pixrad)/pixrad

endif;}

print,' Enter the angle you want to plot the radial scan --> '
read,pltang

dispipB = bytscl(pBcomp,min=ymn,max=5000,top=240)
dispiB  = bytscl(Bcomp,min=ymn,max=5000,top=240)
dispiM  = bytscl(mk3,min=ymn,max=5000,top=240)



;for pltang=0.,360.,15.  do begin ;{

    dp  =  pltang * acirc  +  acirc * 90.0
    dpM =  pltang * acirc  +  acirc * 90.0

    dpx = fix ( cos(dp)  * pixrad * apltrad + xcen + rnder )
    dpy = fix ( sin(dp)  * pixrad * apltrad + ycen + rnder )
    dpxM= fix ( cos(dpM) * pixrad * apltrad + xcenM+ rnder )
    dpyM= fix ( sin(dpM) * pixrad * apltrad + ycenM+ rnder )

    xdat = intarr(2) & xdat = [ 1. , 2.5 ]
    ydat = intarr(2) & ydat = [ 10 , ymx ]  
    mytitle = $
    'Radial scan of B and pB*1.0^10 BBsun at '+$
	string(format='(f5.1)',pltang) +' degrees' + '      ' + mydate
    myxtitle='Radius (Rsun)'


    mycol=246
    pltarrB  = Bcomp(dpx,dpy)
    pltarrpB = pBcomp(dpx,dpy)
    pltarrM  = mk3(dpxM,dpyM)

    wset,3

    if ( dolog eq 'y' ) then $
	plot_io,xdat,ydat,title=mytitle,/nodata,xtitle=myxtitle
    if ( dolog eq 'n' ) then $
	plot,xdat,ydat,title=mytitle,/nodata,xtitle=myxtitle

	oplot,fakex,pltarrpB,col=pBcol,lin=3.0,thi=2
	oplot,fakex,pltarrB,col=Bcol,lin=0.0,thi=2
	oplot,fakex,pltarrM,col=Mcol,lin=5.0,thi=2
	oplot,[2.25,2.35],[ymx-0.35*ymx,ymx-0.35*ymx],$
	    col=Bcol,lin=0.0,thi=2
	oplot,[2.25,2.35],[ymx-0.55*ymx,ymx-0.55*ymx],$
	    col=pBcol,lin=3.0,thi=2
	oplot,[2.25,2.35],[ymx-0.75*ymx,ymx-0.75*ymx],$
	    col=Mcol,lin=5.0,thi=2

	xyouts,2.,ymx-0.4*ymx,' B composite',col=Bcol,charsize=1.1
	xyouts,2.,ymx-0.6*ymx,'pB composite',col=pBcol,charsize=1.1
	xyouts,2.1,ymx-0.75*ymx,'mk3',col=Mcol,charsize=1.2


    wset,0
	dispiB (dpx,dpy)   = 246    & tv,dispiB
	dispipB(dpx,dpy)   = 246    & tv,dispipB,512  ,0
	dispiM (dpxM,dpyM) = 246    & tv,dispiM, 512*2,0

    if ( printdata eq 'y' ) then begin ;{
	filenm=string(format='(f4.2)',pltang)+'pBdat'
	filenm=strtrim(filenm,2)
	openw,1,filenm,/append
	printf,1,mytitle
	printf,1,' '
	printf,1,'      Radius      composite B       composite pB      mk3 '
	for ii=0,156-1,1 do begin ;{
	    printf,1,ii,'    ',pltarrB(ii),'      ',pltarrpB(ii),'      ',pltarrM(ii)
	endfor ;}
	printf,1,' '
	printf,1,' '
	close,1
    endif ;}

;endfor ;}

end
