;+
;
; Project   : renderer wrappers
;                   
; Name      : dsumidl
;               
; Purpose   : Thomson scattering rendering done completely in IDL
;               
; Explanation: part of Pixon and stand-alone renderer
;               
; Use       : see documentation in code
;    
; Inputs    : (see code)
;               
; Outputs   : (see code)
;
; Keywords  : (see code)
;               
; Common    : 
;               
; Restrictions: none
;               
; Side effects: 
;               
; Category    : display, rendering, physics
;               
; Prev. Hist. : None.
;
; Written     : Sandy Antunes, NRC, March-April 2009
;               
;-            
; 'no pixel interpolation' and 'pixel interpolation'
;  from pxndata__dsf_default.pro
; (later rewritten into disum.c)
;
; Includes experimental interpolation routine 'voxelcloudadd',
; which only gets called if you provide a 'deltavoxel' interpolation
; factor.

PRO dsumidl,interpolate,N,i,j,la,image,fizzix,dat,deltavoxel

    if (n_elements(deltavoxel) ne 0) then begin
        ; cloud interpolation, great for uneven pixel-voxel matchups
        for a=0,N-1 do for b=0,N-1 do begin             
            i0=round(i[a,b])
            j0=round(j[a,b])

            signal=image[a,b,la]*fizzix[a,b]
            voxelcloudadd,deltavoxel,signal,i0,j0,N,dat
        end

    end else if (interpolate ne 0) then begin
        ; interpolate v=1

        for a=0,N-1 do for b=0,N-1 do begin

            i0=round(i[a,b])
            si=i[a,b]-i0
            di=abs(si)
            si=round(si/di)
                                
            j0=round(j[a,b])
            sj=j[a,b]-j0
            dj=abs(sj)
            sj=round(sj/dj)

            f = image[a,b,la]*fizzix[a,b]

            ii=i0
            jj=j0
            if((ii ge 0) and (ii lt N) and (jj ge 0) and (jj lt N))then $ 
              dat[ii,jj]=dat[ii,jj]+f*(1-di)*(1-dj)
                                
            ii=i0+si
            jj=j0
            if((ii ge 0) and (ii lt N) and (jj ge 0) and (jj lt N))then $ 
              dat[ii,jj]=dat[ii,jj]+f*di*(1-dj)
                                
            ii=i0
            jj=j0+sj
            if((ii ge 0) and (ii lt N) and (jj ge 0) and (jj lt N))then $ 
              dat[ii,jj]=dat[ii,jj]+f*(1-di)*dj
                                
            ii=i0+si
            jj=j0+sj
            if((ii ge 0) and (ii lt N) and (jj ge 0) and (jj lt N))then $ 
              dat[ii,jj]=dat[ii,jj]+f*di*dj

        end


    end else begin
        ; no interpolation
        for a=0,N-1 do for b=0,N-1 do begin             
            i0=round(i[a,b])
            j0=round(j[a,b])

            if((i0 ge 0) and (i0 lt N) and            $
               (j0 ge 0) and (j0 lt N))then           $
              dat[i0,j0]=dat[i0,j0]+                 $
              image[a,b,la]*fizzix[a,b]
        end

    end


END