;+
;PROCEDURE: thm_part_time_interpolate
;PURPOSE:  Interpolate particle data to match the time grid of another distribution, or to an arbitary time grid
;
;INPUTS:
;  source: A particle dist_data structure to be modified by interpolation to match target
;  target: A particle dist_data structure whose times will be matched, or an array of target times for interpolatuion.
;  
;OUTPUTS:
;   Replaces source with a time interpolated dist_data structure
;
;KEYWORDS: 
;  error: Set to 1 on error, zero otherwise
;  
; NOTES:
;   Any target times that occur between modes will contain samples filled with NANs.  Effective interpolation is very very tricky between modes.(Read: It will probably never happen)
;   Accepts any keywords normally accepted by IDL interpol
;   
;   This has a done of TBDs, we need to come back and fix them when time is available.
; SEE ALSO:
;   thm_part_dist_array, thm_part_smooth, thm_part_subtract,thm_part_omni_convert
;
;  $LastChangedBy: pcruce $
;  $LastChangedDate: 2013-03-07 09:47:09 -0800 (Thu, 07 Mar 2013) $
;  $LastChangedRevision: 11743 $
;  $URL: svn+ssh://thmsvn@ambrosia.ssl.berkeley.edu/repos/thmsoc/tags/tdas_8_00/idl/themis/spacecraft/particles/thm_part_time_interpolate.pro $
;-


pro thm_part_time_interpolate,source,target,error=error,_extra=ex

  compile_opt idl2

  error = 1
  
  ;input validation TBD
  
  ;handling of time array inputs TBD 
  
  ;TBD interpolate all fields
  
  ;concatenate time list from all target modes
  for i = 0,n_elements(target)-1 do begin
    target_mid_times = array_concat(((*target[i]).time+(*target[i]).end_time)/2d,target_mid_times)        
    target_delta_times = array_concat(((*target[i]).end_time-(*target[i]).time),target_delta_times)    
  endfor
  
  
  for i = 0,n_elements(source)-1 do begin
  
    ;handle target times that occur between modes.  Since we can't interpolate these reliably in the general case, we're treating them as missing data
    if i eq 0 then begin ;This handles the case where we're checking for all targets before the first source mode
      idx = where(target_mid_times lt min((*source[0]).time),c)
      if c gt 0 then begin
        out_dist = ptr_new(replicate((*source[0])[0],c))
        (*out_dist[0]).time = target_mid_times[idx]-target_delta_times[idx]/2d ;using the target time width to resize the time bin after interpolation only works because the spacecraft runs in snapshot mode 
        (*out_dist[0]).end_time = target_mid_times[idx]+target_delta_times[idx]/2d 
        (*out_dist[0]).data = !VALUES.D_NAN
      endif
    endif else begin ;this handles the case where we're checking for all targets between modes(we can't interpolate between modes because the data may not match)
      idx = where(target_mid_times lt min((*source[i]).time) and target_mid_times gt max((*source[i-1]).end_time),c)
      if c gt 0 then begin
        out_dist = array_concat(ptr_new(replicate((*source[i])[0],c)),out_dist)
        (*out_dist[n_elements(out_dist)-1]).time = target_mid_times[idx]-target_delta_times[idx]/2d ;using the target time width to resize the time bin after interpolation only works because the spacecraft runs in snapshot mode 
        (*out_dist[n_elements(out_dist)-1]).end_time = target_mid_times[idx]+target_delta_times[idx]/2d
        (*out_dist[n_elements(out_dist)-1]).data = !VALUES.D_NAN
      endif
    endelse
    
    ;This handles the case where we're checking for targets within modes
    idx = where (target_mid_times ge min((*source[i]).time) and target_mid_times le max((*source[i]).end_time),c)
    if c gt 0 then begin 
      out_data = replicate((*source[i])[0],c)
      
      ;unneed, time assignment is handled by the interpolator above
      ;out_data.time = target_mid_times[idx]-target_delta_times[idx]/2d
      ;out_data.end_time = target_mid_times[idx]+target_delta_times[idx]/2d
  
      tag_names=tag_names(out_data[0])
  
      for j=0,n_elements(tag_names)-1 do begin
        if is_num(out_data[0].(j)) then begin
          for k = 0,n_elements(out_data[0].(j))-1 do begin
            out_data.(j)[k] = interpol((*source[i]).(j)[k],((*source[i]).time+(*source[i]).end_time)/2d,target_mid_times[idx],_extra=ex) ;The assumption with this interpolation is that it will be valid for any data that are changing over time within a mode.  For others, it should not change the data.(because it is constant)
          endfor
        endif else if n_elements(out_data[0].(j)) eq 1 then begin
          out_data.(j) = (*source[i])[0].(j)
        endif else begin
          message,'Panic!'
        endelse
      endfor
      
      out_dist = array_concat(ptr_new(out_data),out_dist)
      
    endif
  endfor
  
  ;This handles the case where we're checking for all targets after the last mode
  idx = where(target_mid_times gt max((*source[n_elements(source)-1]).end_time),c)
  if c gt 0 then begin
    out_dist = array_concat(ptr_new(replicate((*source[n_elements(source)-1])[0],c)),out_dist)
    (*out_dist[n_elements(out_dist)-1]).time = target_mid_times[idx]-target_delta_times[idx]/2d ;using the target time width to resize the time bin after interpolation only works because the spacecraft runs in snapshot mode 
    (*out_dist[n_elements(out_dist)-1]).end_time = target_mid_times[idx]+target_delta_times[idx]/2d
    (*out_dist[n_elements(out_dist)-1]).data = !VALUES.D_NAN
  endif
  
  source=out_dist 
  heap_gc ; remove pointers whose references were just deleted
  
  error = 0

end