;+
; NAME:
;	brss_csss_scbf
; PURPOSE:
; CALLING SEQUENCE:
;	brss_csss_scbf,scbf,obs,ghfcsss,brfrss,IPH=iph,JTH=jth, NMAX=nmax,APAR=apar,RCP=rcp,RSS=rss
;
;	brss_csss_scbf
;		--- gh_mbrsc (non-WSO)
;		--- gh_wbrsc (WSO)  --- readwsosc
;							--- zgrid
;							--- pdp			--- pdpi
;							--- csmph
;		--- zgrid
;		--- pdp			--- pdpi
;		--- csmph
;		--- sb_hc
;		--- ghv_csss
;		--- sb_csss
;		--- writefits
;-

;+
; NAME:
;	readwsosc
; PURPOSE:
;	Read a WSO SC data file obtained from WSO Web
; CALLING SEQUENCE:
;	data = readwsosc(scfile,/cali,/br,/disp)
; INPUTS:
;	scfile			scalar; type: string
;						fully-qualified name of file containing
;						photospheric magnetic field data.
; OPTIONAL INPUT PARAMETERS:
;	/cali			multiplies output with 4.5-2.5*cos^2(theta)
;	/br				multiplies output with 1.0/sin(theta)
;	/disp			shows contour plot of data array
; CALLS:
;	zgrid, pimg0, sc_annot
; OUTPUTS:
;	data			array[72,30]; type: float
;						photospheric magnetic field data
;						longitudinal grid in 5 degree steps
;						latitudinal grid in equal steps of sine(lat)
;						Longitude and latitude increases monotonically with
;						array index.
; RESTRICTIONS:
;	XuePus original procedure assumed headers of two lines (prior to the
;	first record starting with CT). Instead we now use wso_read to read
;	to read the photospheric files (which does not have this restriction).
; PROCEDURE:
;	The standard WSO file format contains 73x30 data points, organized
;	in 73 groups of 4 records for 73 heliographic longitudes at 5 degree
;	intervals in monotonic decreasing order. I.e. the first and last group
;	of four records correspond to the same heliographic longitude but
;	one rotation later.
;
;	Each group of 4 records contains data on a latitudinal grid equally
;	spaced in heliographic latitude covering -14.5/15 to +14.5/15 in
;	monotonically decreasing order (i.e. going from north to south).
;
;	Only the first 72 groups of records are read. The last is discarded.
;	Both the longitudinal and latitudinal order in the ouput array is
;	reversed from the order in the file, i.e. both longitude and latitude
;	in the 72x30 output array increases monotonically with array index.
;-

;+
; NAME:
;	zgrid
; PURPOSE:
;	Calculate theta, phi angles and their functions
; CALLING SEQUENCE:
;	zgrid,jph,ith, phd,thd,lad,cth,sth,cmp,scs, /rup
; INPUTS:
;	jph				scalar; type: integer
;						# longitudes (72)
;	ith				scalar; type: integer
;						# latitudes (30)
;	cmp				scalar; type: float
;						offset angle in degrees
; OPTIONAL INPUT PARAMETERS:
;	/rup			reverses order of all output arrays
;					(phd,thd,lad,cth,sth,scs)
; OUTPUTS:
;	phd				array[jph]; type: float
;						heliographic longitude set up as equally
;						spaced grid in degrees from
;						0+0.5*360/jph, .., +360-0.5*360/jph
;
;	thd				array[ith]; type: float
;						theta, polar angle in degrees
;						theta=0 at N, theta=pi at S
;	lad				array[ith]; type: float
;						heliographic latitude in degrees
;						lad = 90-thd
;	cth				array[ith]; type: float
;						cos(theta); cosine of polar angle
;						cth is set up as an equally spaced array
;						from -1+1/ith, .., +1-1/ith
;	sth				array[ith]; type: float
;						sin(theta); sine of polar angle
;
;	scs				array[jph]; type: float
;						1.0/cos( (phd-cmp)*!dtor )
;						(only calculated if cmp is specified)
;-

;+
; NAME:
;	gh_wbrsc
; PURPOSE:
;	Compute GENERAL G & H based on WSO data
; CALLING SEQUENCE:
;	gh_wbrsc, wbrsc, Nmax, gw, hw, GHWF=ghwf, IFL=ifl, CALI=cali
; INPUTS:
;	wbrsc		/ifl SET:
;				scalar; type: float
;					name of file with photospheric magnetic field data
;					(/cali is used if specified; see below)
;				/ifl NOT set:
;				array[72,30];
;					photospheric magnetic field data
;					Note that the number of longitudes is 72, NOT 73
;					(see readwsosc). (/cali is not used in this case)
;	Nmax		scalar; type: integer
;					??
;
; OPTIONAL INPUT PARAMETERS:
;	ghwf=ghwf	scalar; type: string
;					specifies file where calculated gw and hw are saved
;
;	/ifl		if set then argument 'wbrsc' must be a filename
;					the file name is passed to 'readwsosc'
;	/cali		used only if argument 'wbrsc' is a filename (and /ifl is set)
;					passed to 'readwsosc'
;
; OUTPUTS:
;	gw			array[Nmax+1,Nmax+1]; type: float
;					??
;	hw			array[Nmax+1,Nmax+1]; type: float
;					??
; CALLS:
;	readwsosc, pdp, csmph
;-

;+
; NAME:
;	pdp
; PURPOSE:
;	Produce arrays for P and dP
; CALLING SEQUENCE:
;	pdp, cth, Nmax, P, dP
;
;	Nmax le 63 OK for the case of ith=30    04Dec1999
;-

;+
; NAME:
;	pdpi
; PURPOSE:
; CALLING SEQUENCE:
;	pdpi, cth, Nmax, ps, dps
; INPUTS:
;	cth
;	Nmax
; OUTPUTS:
;	ps
;	dps
;
;  pdpi:        Calculate Schmidt-type associated Legendre function using
;               P_m^m=(2-dm0)^0.5 ((2m-1)!!/2m!!)^0.5 (1-x^2)^{m/2}
;                     dm0=1 when m=0 and dm0=0 when m>0
;  P_{l+1}^m=((l+1)^2-m^2)^{-0.5}[(2l+1) x P_l^m - (l^2-m^2)^0.5 P_{l-1}^m]
;               dP_m^m=m (2-dm0)^0.5 ((2m-1)!!/2m!!)^0.5 x (1-x^2)^{(m-1)/2}
;  dP_{l+1}^m=((l+1)^2-m^2)^{-0.5} [(2l+1)[-(1-x^2)^0.5 P_l^m + x dP_l^m] -
;                                    (l^2-m^2)^0.5 dP_{l-1}^m]
;  output:      P(Nmax+1,Nmax+1), dP(Nmax+1,Nmax+1) for specified thetad
;               including (Nmax+1)(Nmax+2)/2 elements
;  wrote:       Feb. 28, 1996
;-