;+
; NAME: pa_mike_deadtime.pro
;
; PURPOSE: Apply the deadtime correction to P-Alice data
;
; CATEGORY: Data processing
;
; CALLING SEQUENCE: deadtime_correction = pa_deadtime(cntrate)
;
; INPUTS: 
;   cntrate -- the measured countrate (Hz) of the data
;
; OUTPUTS:
;   Returns the factor by which to multiply the data to yield the
;   count rate that would have been detected given infinitely fast
;   electronics. 
;
; COMMON BLOCKS: None
;
; SIDE EFFECTS: None
;
; NOTES: 
;   The deadtime correction time constant is 16.2 microseconds. At input
;   rates below 50 kHz, the detector electronics is non-paralyzable
;   (fixed deadtime per event that is not re-triggerable). To calculate
;   the detector output count rate, the following formula is used:
;
;   r_out = r_in / (exp[r_in * tau_p] + r_in * tau_n)
;   
;
; MODIFICATION HISTORY:
;  V 1  2007-02-19  Andrew J. Steffl
;  V 2  2008-01-23  AJS corrected algorithm to treat detector
;                   electronics as non-paralyzable.
;  V 3  2008-01-25  AJS Added check if count rate is above 50 kHz
;  v 4  2009-02-06  AJS Added hybirid deadtime model
;-


FUNCTION ra_mike_deadtime, cntrate

ncntrate = n_elements(cntrate)
IF ncntrate EQ 0 THEN BEGIN
   message, '%RA_MIKE_DEADTIME: CNTRATE not defined!'
   return, 1.
ENDIF

tau_n = 15.7e-6
tau_p = 11.3e-6

; For this set of time constants, the detected counts rolls over at 88.509 kHz
; input rate. Assume that whatever Rosetta is looking at, it has a lower
; input count rate than this.
r_in = findgen(88510) + 1.

r_out = r_in / (exp(r_in * tau_p) + r_in * tau_n)

deadtime_correction = interpol(r_in/r_out, r_out, cntrate)

return, deadtime_correction
end