[neutron-mc] Error in component Source_simple.comp ?

[Frederik Zilly]

Hello Peter !

Thank you very much for your quick answer.
Indeed, I must have been blind not seeing the square root.
With this square root everything is alright as your jpg shows.

Sorry for this "false alarm".

Cheers,

Frederik.

Peter Willendrup wrote:

> Hello Frederik,
> 
> On Tue, 2006-02-14 at 16:12 +0100, Frederik Zilly wrote:
> 
>>Points near the center have greater probility than points with a higher 
>>radius.
> 
> 
> Actually, this is not the case. The square root actually gives zero
> probability for points exactly at the origin; more events must be
> distributed over a circular area element at the circle perimeter  than
> at the center to achieve a uniform distribution.
> 
> I wrote up a few lines of matlab code to show this in graphics (you need
> the attached hist2d code to reproduce the attached graphics):
> 
> chi=2*pi*rand(1e6,1);
> r=sqrt(rand(1e6,1)); % Source radius assumed one metre
> x=r.*cos(chi);
> y=r.*sin(chi);
> mXY=[x y];
> vXEdge = linspace(-1,1,128);
> vYEdge = linspace(-1,1,128);
> img=hist2d(mXY,vXEdge,vYEdge);
> imagesc(img);
> 
> As you can see in the jpg graphic, the result of the above code (which
> is equvalent to the c-code in Source_simple) is completely uniform
> across the source area.
> 
> Hope this helps?
> 
> 
> Cheers,
> 
> Peter
> 
> 
> ------------------------------------------------------------------------
> 
> %function mHist = hist2d ([vY, vX], vYEdge, vXEdge)
> %2 Dimensional Histogram
> %Counts number of points in the bins defined by vYEdge, vXEdge.
> %size(vX) == size(vY) == [n,1]
> %size(mHist) == [length(vYEdge) -1, length(vXEdge) -1]
> %
> %EXAMPLE
> %   mYX = rand(100,2);
> %   vXEdge = linspace(0,1,10);
> %   vYEdge = linspace(0,1,20);
> %   mHist2d = hist2d(mYX,vYEdge,vXEdge);
> %
> %   nXBins = length(vXEdge);
> %   nYBins = length(vYEdge);
> %   vXLabel = 0.5*(vXEdge(1:(nXBins-1))+vXEdge(2:nXBins));
> %   vYLabel = 0.5*(vYEdge(1:(nYBins-1))+vYEdge(2:nYBins));
> %   pcolor(vXLabel, vYLabel,mHist2d); colorbar
> function mHist = hist2d (mX, vYEdge, vXEdge)
> [nRow, nCol] = size(mX);
> 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);
>     
>     [mIdxRow] = find (vRow > rRowLB & vRow <= rRowUB);
>     vColFound = vCol(mIdxRow);
>     
>     if (~isempty(vColFound))
>         
>         
>         vFound = histc (vColFound, vXEdge);
>         
>         nFound = (length(vFound)-1);
>         
>         if (nFound ~= nCol)
>             [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
> 
> 
> 
> 
> ------------------------------------------------------------------------
> 
> 
> ------------------------------------------------------------------------
> 
> _______________________________________________
> neutron-mc mailing list
> neutron-mc at risoe.dk
> http://mailman.risoe.dk/mailman/listinfo/neutron-mc