matthiasm@8: function result = logsum(logv)
matthiasm@8: 
matthiasm@8: len = length(logv);
matthiasm@8: if (len<2);
matthiasm@8:   error('Subroutine logsum cannot sum less than 2 terms.');
matthiasm@8: end;
matthiasm@8: 
matthiasm@8: % First two terms
matthiasm@8: if (logv(2)<logv(1)),
matthiasm@8:   result = logv(1) + log( 1 + exp( logv(2)-logv(1) ) );
matthiasm@8: else,
matthiasm@8:   result = logv(2) + log( 1 + exp( logv(1)-logv(2) ) );
matthiasm@8: end;
matthiasm@8: 
matthiasm@8: % Remaining terms
matthiasm@8: for (i=3:len),
matthiasm@8:   term = logv(i);
matthiasm@8:   if (result<term),
matthiasm@8:     result = term   + log( 1 + exp( result-term ) );
matthiasm@8:   else,
matthiasm@8:     result = result + log( 1 + exp( term-result ) );
matthiasm@8:   end;    
matthiasm@8: end