; docformat = 'rst'
;
; NAME: 
;    ssoSynthSFD
; PURPOSE:
;    Create histograms and synthetic Size-Frequency Distribution of SSOs
;+
; :Description:
;    Create histograms and synthetic Size-Frequency Distribution of SSOs
;
; :Categories:
;    Database, SSO
;
; :Params:
;    pop: in, required, type=fltarr
;        The absolute magnitude of each member of a population.
;
; :Returns: 
;
; :Keywords:
;    range: in, optional, type=fltarr(2), default=[2,24]
;        The absolute magnitude range over which to compute the histograms.
;    binSize: in, optional, type=float, default=0.5
;        The size of absolute magnitude bins for histogram.
;    reference: in, optional, type=fltarr(2), default=No default
;        If set, compute a synthetic population by extrapolating a
;        power law SFD, fitted over the reference range 
;
; :Examples:
;
; :Uses:
;
; :Author:
;    B.Carry (OCA)
;
; :History:
;    Change History::
;      Original Version written in January 2016, B. Carry (OCA)
;      2017 Feb. - B. Carry (OCA) - Switched the fitting from cumulative to differential distribution
;-
function ssoSynthSFD, pop, range=range, binSize=binSize, reference=reference
  COMPILE_OPT hidden, idl2
;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;
;---  TAG ---    I -- Initialization And Input Verification              -----------------------;
;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;

 ;--I.1-- Test Input Population  ---------------------------------------------------------------
  if not keyword_set(pop) then begin
    message, /ioError, 'syntax: sfd = ssoSynthSFD( pop, range=, binSize=, reference=])'
    return, -1
  endif

 ;--I.2-- Default Histogram Parameters  --------------------------------------------------------
  if not keyword_set(binSize) then binSize=0.5
  if not keyword_set(range)   then range=[2,24]
    nbBin = ( range[1]-range[0] )/binSize + 1



;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;
;---  TAG ---   II -- Compute Population Size-Frequency Distribution     -----------------------;
;-----------------------------------------------------------------------------------------------;
;-----------------------------------------------------------------------------------------------;

 ;--II.1-- Default Histogram Parameters  -------------------------------------------------------
  histoH = histogram( pop, binSize=binSize, min=range[0], max=range[1], loc=hArr )
  cumulH = total( histoH, /Cumulative )


 ;--II.2-- Not Synthetic Population to Create  -------------------------------------------------
  if not keyword_set(reference) then begin

    return, {h:hArr, histo:{obs:histoH}, cumul:{obs:cumulH}}


 ;--II.3-- Extrapolate if Requested  -----------------------------------------------------------
  endif else begin



;-oldd   ;--II.2.1-- Synthetic Cumulative Distribution - log N(<H) = C + aH
;-oldd    low=where(hArr le reference[0])
;-oldd    sel = where( hArr ge reference[0] and hArr le reference[1] )
;-oldd    lin = linfit( hArr[sel], alog10(cumulH[sel]) )
;-oldd      cumulM = 10^(lin[0] + hArr * lin[1])
;-oldd
;-oldd   ;--II.2.2-- Synthetic Distribution
;-oldd    histoM=cumulM*0.
;-oldd    for kB=nbBin-1,1,-1 do $
;-oldd       histoM[kB] = cumulM[kB] - cumulM[kB-1]
;-oldd
;-oldd   ;--II.2.3-- Fix Synthetic Distributions
;-oldd    cumulM[low]=cumulH[low]
;-oldd    histoM[low]=histoH[low]
;-oldd


   ;---new     
   ;--II.2.1-- Synthetic Cumulative Distribution - log N(=H) = C + aH
    low=where(hArr le reference[0])
    sel = where( hArr ge reference[0] and hArr le reference[1] )
    lin = linfit( hArr[sel], alog10(histoH[sel]) )
      histoM = 10^(lin[0] + hArr * lin[1])

   ;--II.2.2-- Synthetic Distribution
    cumulM = total( histoM, /Cumulative )

   ;--II.2.3-- Fix Synthetic Distributions
    cumulM[low]=cumulH[low]
    histoM[low]=histoH[low]

    return, {h:hArr, histo:{obs:histoH, model:histoM}, cumul:{obs:cumulH, model:cumulM, expo:lin}}


  endelse

end