;+ ; NAME: ; INERTIATENS ; ; PURPOSE: ; Calculates the 2nd moment tensor (sometimes incorrectly referred to ; as the moment of inertia tensor) of a mass distribution specified by ; a list of particle positions. ; ; CATEGORY: ; Astro ; ; CALLING SEQUENCE: ; Result = INERTIATENS(Pos) ; ; INPUTS: ; Pos: An Nx3 array specifying the 3d positions of the N particles ; that make up the mass distribution. ; ; KEYWORD PARAMETERS: ; MASSES: An N-element vector of the mass of each point. If not ; specified, all masses are assumed to be unity. ; ; R2WEIGHT: If /R2WEIGHT is specified then particles are downweighted ; by a factor of 1/r^2 so that all particles have equal ; effect regardless of radius. ; ; OUTPUTS: ; The function returns a 3x3 symmetric array containing the 2nd moment ; of the mass distribution tensor, i.e. Result[i,j] is the sum over each ; particle k of MASSES[k] * Pos[k,i] * Pos[k,j]. ; ; EXAMPLE: ; Calculate the inertia tensor of 6 equal-mass points distributed on ; the vertices of a cube: ; ; xmasspos = [-1,-1,-1,-1,1,1,1,1] ; ymasspos = [-1,-1,1,1,-1,-1,1,1] ; zmasspos = [-1,1,-1,1,-1,1,-1,1] ; itens = INERTIATENS([[xmasspos],[ymasspos],[zmasspos]]) ; ; MODIFICATION HISTORY: ; Written by: Jeremy Bailin ; 10 June 2008 Public release in JBIU ; 22 July 2011 Bug fix for /R2WEIGHT to actually do something useful ;- function inertiatens, pos, masses=masses, r2weight=r2weight Itens = dblarr(3,3) npts = (size(pos,/dimen))[0] weights = replicate(1.0,npts) if keyword_set(r2weight) then begin posnot0 = where(total(pos^2,2) ne 0., nposnot0) if nposnot0 gt 0 then weights[posnot0] /= total(pos^2, 2) endif if n_elements(masses) ne 0 then weights *= masses for i=0,2 do begin for j=i,2 do begin Itens[i,j] = total( weights * pos[*,i] * pos[*,j] ) Itens[j,i] = Itens[i,j] endfor endfor return, Itens end