;+
; NAME:
;	DISCP_RESP
;
; PURPOSE:
;	Return the count rate response of the DISCSP to a model spectrum
;	based on angle and energy boundaries.
;
; CATEGORY:
;	BATSE, SPEX
;
; CALLING SEQUENCE:
;	DISCSP_RESP, Edges, Theta, Apar, Resp, SCALE=SCALE
;
; INPUTS:
;	Edges: Lower and upper edge in keV for 1 channel.
;		if more than 1 channel, should be an array [lo,hi] x Number_channels
;	Theta: Angle between source and detector axis in degreees.
;	Apar:  Parameters needed for spectral model function.
;	
; KEYWORD PARAMETERS:
;	MODEL: String name of spectral model, default, "F_BPOW"
;	SCALE: Scale factor for resolution above PHA values, passed onto
;	       BATSE_SPEC_DRM
;
; OUTPUTS:
;   	Resp: Count rate in interval.
;
; COMMON BLOCKS:
;	DISCSP_RESP:	edges_out, adrm, e_in, atheta, adrm, scale_save
;
; PROCEDURE:
;	Builds a grid of response matrices, computes count rates, then interpolates.
;
; RESTRICTIONS:
;	Edges should be between 3. and 1000. keV
; EXAMPLE:
;	IDL> print,apar                           
;   	0.00250560      2.12000     0.000000     0.000000
;	IDL> discsp_resp, [5.,300],20.,apar, resp &print,resp
;	      50.6533
;
; MODIFICATION HISTORY:
;	ras, 22-dec-1995
;	ras, 23-sep-1996 keyword inheritance to model function call
;	ras, 10-oct-1996 added multi-channel capability
;	ras, 17-oct-1996, added protection against narrow output edges
;               RAS, 7-Jul-1997, changed PERM_DATA to SSWDB_BATSE
;-

pro discsp_resp, edges, theta, apar, resp, model=model, scale=scale, _extra=_extra, $
	ein_lim=ein_lim

;on_error,2
parcheck, edges, 1, [4,5,2,3], indgen(2)+1,'Energy Edges'
parcheck, theta, 2, [4,5,2,3], indgen(2),'Zenith Angle in Degrees'
parcheck, apar, 3, indgen(4)+2, indgen(2),'Model Parameters'

sedges = size(edges)
if n_elements(edges) gt 2 then begin
	resp = fltarr(sedges(2))
	for i=0,sedges(2)-1 do begin
	discsp_resp, edges(*,i), theta, apar, respi, model=model, scale=scale, _extra=_extra, $
	ein_lim=ein_lim
	resp(i)=respi
	endfor
	return
endif

common discsp_resp, edges_out, adrm, e_in, atheta, scale_save
if not exist(edges_out) then restore, concat_dir((chklog('SSWDB_BATSE'))(0), 'discsp_resp.sav')

checkvar, model, 'f_bpow'
if exist(scale) then begin
	checkvar, scale, fcheck(scale_save,1.)
	scale_save = scale
	edges_out=[3+findgen(14)*.5,10*exp(findgen(080)*.03), 110*exp(findgen(33)*.07)]
	edge_products, edges_out,edges_2=edges_out, wid=wedge, mean=emdge

	batse_spec_drm,edges_out=edges_out, edges_in=e_in,theta=0.,$
	  det=0,drm=drm,sigma=sigma,e_in=eph, scale=scale 

	sdrm=size(drm)

	atheta=[0., 15., 30., 40.+findgen(7)*7.]
	adrm = fltarr(sdrm(1),sdrm(2),n_elements(atheta))

	adrm(0,0,0)=drm
	for i=1,n_elements(atheta)-1 do begin 
	  batse_spec_drm,edges_out=edges_out, edges_in=e_in,theta=atheta(i), $
		det=0,drm=drm,sigma=sigma,e_in=eph 
	  adrm(0,0,i)=drm 
	endfor

endif
edge_products,e_in,gm=em, width=win
edge_products,edges_out, width=wedge
wlim = where( em lt fcheck(ein_lim,max(em)) )
if keyword_set(_extra) then $
	modl=call_function( model, em(wlim), apar, _extra=_extra ) else $
	modl=call_function( model, em(wlim), apar)

;Identify the two bounding values of atheta
wtheta= where( atheta gt theta, ntheta)
ibnd = [-1,0]+wtheta(0)

;Identify the two energy intervals containing the input energies
iedge0 = (where( edges(0) ge edges_out(0,*) and edges(0) le edges_out(1,*)))(0)
iedge1 = ( where( edges(1) ge edges_out(0,*) and edges(1) le edges_out(1,*)))(0)

n = iedge1-1 - (iedge0+1) + 1
cpb = fltarr( n_elements(wedge), 2) ;what are the counts per output bin?

for j=0,1 do $
	cpb(iedge0,j) = adrm( iedge0:iedge1,wlim,ibnd(j) )# (modl*win)(wlim) * wedge(iedge0:iedge1)

;Integrate the response over the range for both values of theta
if n ge 1 then begin
	w = indgen(n)+iedge0+1
	resp = total( cpb(w,*),1)
endif else resp = cpb(0,*)*0.0 ;no complete bins to sum over
;Add in the fractional bins
for j=0,1 do resp(j) = resp(j) + ( cpb(iedge0,j)/wedge(iedge0))*(edges_out(1,iedge0)-edges(0)) +$
			         ( cpb(iedge1,j)/wedge(iedge1))*(edges(1)-edges_out(0,iedge1))

resp = resp * 127.
resp = interpol( resp, atheta(ibnd), theta)
	
end