Re: weighted 2D histogram

[lynx . abraxas]

On 19/12/09 16:12:45, Søren Hauberg wrote:
> lør, 19 12 2009 kl. 13:11 +0100, skrev address@hidden:
>  >  Now  I  need  a weighted version of this such that each cell of the  2D-
histogram
> > is not just the amount of points in that cells area but  the  sum  of  all
weights
> > that come with the points that ly in that cell.
>
> Essentially, the code you sent is merely a wrapper around 'histc'. It
> seems that 'histc' does not support weights, so I don't think the
> changes will be trivial. I think the 'histc' API could be extended to
> support optional weights, but that would need some discussion.


Thanks Søren.
I  had a look at histc and changed it to histc_w used that in my hist2d_w. The
code of both is below.
In the end it was quite simple so I now wonder, have I  done  anything  wrong?
Have I forgotten something??? I'm unsure about my vW calculation.


Regards,
Lynx

_______________________________

function freq = histc_w (X, e, W)

freq = zeros(length(e),1);

for i=1:length(X)
  for k=1:length(e)-1
    if(e(k) <= X(i) && X(i) < e(k+1))
      freq(k)=freq(k) + W(i);
    end
  end

  if ( X(i)== e(end))
    freq(end) = freq(end) + W(i);
  end
end





function mHist = hist2d_w (mX, vYEdge, vXEdge, mW)#change
nCol = size(mX, 2);
if nCol < 2
    error ('mX has less than two columns')
end

nRow = length (vYEdge)-1;
nCol = length (vXEdge)-1;

vRow = mX(:,1);
vCol = mX(:,2);

mHist = zeros(nRow,nCol);

for iRow = 1:nRow
    rRowLB = vYEdge(iRow);
    rRowUB = vYEdge(iRow+1);
    
    vColFound = vCol((vRow > rRowLB) & (vRow <= rRowUB));
    vW = mW((vRow > rRowLB) & (vRow <= rRowUB)); #change
    
    if (~isempty(vColFound))
        
        
        vFound = histc_w (vColFound, vXEdge, vW);#change
        nFound = (length(vFound)-1);
        
        if (nFound ~= nCol)
            disp([nFound nCol])
            error ('hist2d error: Size Error')
        end
        
        [nRowFound, nColFound] = size (vFound);
        
        nRowFound = nRowFound - 1;
        nColFound = nColFound - 1;
        
        if nRowFound == nCol
            mHist(iRow, :)= vFound(1:nFound)';
        elseif nColFound == nCol
            mHist(iRow, :)= vFound(1:nFound);
        else
            error ('hist2d error: Size Error')
        end
    end
    
end