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

fakex=0 & fakex=findgen(360)

if( n_elements(initT) eq 0) then begin ;{
    initT=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='n'

    _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
    dp = findgen(360) * acirc  +  acirc * 90.0
    dpx = intarr(360)
    dpy = intarr(360)

    xcenM = 252.5 & ycenM = 272.5
    dpM = findgen(360) * acirc  +  acirc * 90.0 
    dpxM = intarr(360)
    dpyM = intarr(360)

endif;}

print,' Enter the radius you want to plot the azimuthal scan --> '
read,pltrad

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 pltrad=1.10,1.18,0.2  do begin ;{

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

    xdat = intarr(2) & xdat = [ 0 , 360 ]
    ydat = intarr(2) & ydat = [ 10 , ymx ]  
    mytitle = $
    'Theta scan of B and pB*1.0^10 BBsun at '+$
	string(format='(f4.2)',pltrad) +' Rsun' + '      ' + mydate
    myxtitle='Angle (degrees)'


    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,[360.,380.],[ymx-0.35*ymx,ymx-0.35*ymx],$
	    col=Bcol,lin=0.0,thi=2
	oplot,[360.,380.],[ymx-0.55*ymx,ymx-0.55*ymx],$
	    col=pBcol,lin=3.0,thi=2
	oplot,[360.,380.],[ymx-0.75*ymx,ymx-0.75*ymx],$
	    col=Mcol,lin=5.0,thi=2

	xyouts,300.,ymx-0.4*ymx,' B composite',col=Bcol,charsize=1.1
	xyouts,300.,ymx-0.6*ymx,'pB composite',col=pBcol,charsize=1.1
	xyouts,300.,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)',pltrad)+'pBdat'
	filenm=strtrim(filenm,2)
	openw,1,filenm,/append
	printf,1,mytitle
	printf,1,' '
	printf,1,'      Angle      composite B       composite pB      mk3 '
	for ii=0,360,1 do begin ;{
	    printf,1,ii,'    ',pltarrB(ii),'      ',pltarrpB(ii),'      ',pltarrM(ii)
	endfor ;}
	printf,1,' '
	printf,1,' '
	close,1
    endif ;}

;endfor ;}

end