March 13, 2006 for Nagoya group and UCSD group
* pzgrid,xsz,ysz,phd,thd,lad,cth,sth
* pgh_mbrsc,mbrsc,nmax,gm,hm,ghmf=,/ifl		; n=0,nmax
    pzgrid,xsz,ysz,phd,thd,lad,cth,sth
*   pdp,cth,nmax,P,dP        <=*pdpj,cthj,nmax,pj,dpj
*   csmph,phd,nmax,cmph,smph <=*csmphi,nmax,phri,cmphi,smphi
* ppb_hc,ga,ha,r,phri,cthj,sthj,Br,Bt,Bp,apar=
*   hchf,nmax,apar,r,Hr,dHr,d2Hr ???
    csmphi,nmax,phri,cmphi,smphi
*   gfni,ga,ha,cmphi,smphi,Gfi,dGfi
    pdpj,cthj,nmax,Pj,dPj
*   brtp,nmax,apar,r,sthj,Hr,dHr,Pj,dPj,Gfi,dGfi,Br,Bt,Bp
* psb_hc,ga,ha,r,sth,cmph,smph,P,dP,bbr,bbt,bbp
    hchf,nmax,apar,r,Hr,dHr,d2Hr
*   Gfn,ga,ha,cmph,smph,Gf,dGf <= gfni,ga,ha,cmphi,smphi,Gfi,dGfi
    brtp,nmax,apar,r,sthj,Hr,dHr,Pj,dPj,Gfi,dGfi,Br,Bt,Bp
* pghfcsss_gaha,ga,ha,gca,hca,ghcf
    pzgrid,iph,jth,phd,thd,lad,cth,sth
    csmph,phd,nmax,cmph,smph
    pdp,cth,nmax,P,dP
    g=100*ga
    h=100*ha                         ; from Gauss to mT
    psb_hc,ga,ha,rcp,sth,cmph,smph,P,dP,brcp,btcp,bpcp
    ghv_csss,brcp,btcp,bpcp,sth,cmph,smph,P,dP,apar,rcp,rss,gca,hca,ghcf
* ppb_csss,ga,ha,gca,hca,5.,phri,thrj,br,bt,bp
*   pb0_csss,gc,hc,rr,phri,cth,sth,apar,Rcp,Rss,Br,Bt,Bp
*     cssshf,Nmax,apar,Rcp,Rss,r,Hr,dHr,d2Hr
      csmphi,nmax,phri,cmphi,smphi
      gfni,ga,ha,cmphi,smphi,Gfi,dGfi
      pdpj,cthj,Nmax,Pj,dPj
*     brtp_csss,Nmax,r,sthj,Hr,dHr,Pj,dPj,Gfi,dGfi,Br,Bt,Bp
    ppb_hc,ga,ha,r,phri,cthj,sthj,Br,Bt,Bp
* psb_csss,ga,ha,gca,hca,5,phd,thd,bbr,bbt,bbp,apar=,rcp=,rss=
    ppb_csss,ga,ha,gca,hca,5.,phri,thrj,br,bt,bp
* iflij_csss,ga,ha,gca,hca,phri,thrj,fte,np,xw,yw,zw,rw,tw,pw