pro br260_brf1951

;1. Assign phd,thd,lad,cth,sth 
    pzgrid,360,180,phd,thd,lad,cth,sth
;2. Calcualte ga & ha using nmax=9 and 'mdigh_1951.asc'
    pgh_mbrsc,'MDI.1951.fits',9,ga,ha,ghmf='mdigh_1951.asc',/ifl
  ; print,'ga(0:2,0:2):',ga(0:2,0:2)
;3. Calculate Br,Bt,Bp on a spherical surface of r <= rcp, say r=2.5
    csmph,phd,9,cmph,smph
    pdp,cth,9,P,dP
    psb_hc,ga,ha,2.5,sth,cmph,smph,P,dP,bbr,bbt,bbp
;4. Calculate gca,hca and produce 'ghfcsss_1951.asc'  based on ga,ha,nmax=9
    pghfcsss_gaha,ga,ha,gca,hca,'ghfcsss_1951.asc',nmax=9
  ; print,'gca(0:2,0:2): ',gca(0:2,0:2)
;5. Calculate bbr,bbt,bbp on a spherical surface of, say, r=2.61
    psb_csss,ga,ha,gca,hca,2.61,phd,thd,bbr,bbt,bbp
    writefits,'br_261_1951.fits',bbr
    writefits,'bt_261_1951.fits',bbt
    writefits,'bp_261_1951.fits',bbp

end