;
; PURPOSE:   Compute GENERAL G & H based on MDI, WSO or KPK Br SC data,
;              bra (Br array) and Nmax 
; WRITTEN:   29Mar2001  Xuepu
; MODIFIED:  14FEB2003  Xuepu add CMPT=cmpt,pp,cmph,smph
;
pro ghf_bra,bra,Nmax,gm,hm,CMPT=cmpt,pp,cmph,smph,GHF=ghf,scrn=scrn

    if N_params() LT 1 then begin
       print,'pro ghf_bra,bra,Nmax,gm,hm,/cmpt,pp,cmph,smph,GHMF=,scrn='
       print,'called by JOB/ghf_obs.pro'
       return
    endif
    sz=SIZE(bra)
    xsz=sz(1) & ysz=sz(2)
    if not keyword_set(cmpt) then begin
       zgrid,xsz,ysz,phd,thd,lad,cth,sth
       pdp,cth,Nmax,pp,dP
       csmph,phd,Nmax,cmph,smph
    endif

    if keyword_set(ghf) then begin
       get_lun,unit
       openw,unit,ghf
       printf,unit,'CR'+scrn
       printf,unit,Nmax
    endif
    gm=fltarr(Nmax+1,Nmax+1)
    hm=gm
    ij=xsz*ysz
    for n=0,Nmax do begin         ; n=0 monopole 
        coe=(2*n+1.0)/ij
        for m=0,n do begin
            gm(n,m)=0
            hm(n,m)=0
            for i=0,xsz-1 do begin
                for j=0,ysz-1 do begin
                    gm(n,m)=gm(n,m)+bra(i,j)*pp(n,m,j)*cmph(m,i)
                    hm(n,m)=hm(n,m)+bra(i,j)*pp(n,m,j)*smph(m,i)
                endfor
            endfor
            gm(n,m)=coe*gm(n,m)
            hm(n,m)=coe*hm(n,m)
            if keyword_set(ghf) then printf,unit,n,m,gm(n,m),hm(n,m)
        endfor
    endfor
    if keyword_set(ghf) then free_lun,unit

end