gt:         double logTreeWeight = getLogTreeWeight();
context: /*
   * AbstractObservationProcess.java
   *
   * Copyright (c) 2002-2014 Alexei Drummond, Andrew Rambaut and Marc Suchard
   *
   * This file is part of BEAST.
   * See the NOTICE file distributed with this work for additional
   * information regarding copyright ownership and licensing.
   *
   * BEAST is free software; you can redistribute it and/or modify
   * it under the terms of the GNU Lesser General Public License as
   * published by the Free Software Foundation; either version 2
   * of the License, or (at your option) any later version.
   *
   *  BEAST is distributed in the hope that it will be useful,
   *  but WITHOUT ANY WARRANTY; without even the implied warranty of
   *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
   *  GNU Lesser General Public License for more details.
   *
   * You should have received a copy of the GNU Lesser General Public
   * License along with BEAST; if not, write to the
   * Free Software Foundation, Inc., 51 Franklin St, Fifth Floor,
   * Boston, MA  02110-1301  USA
   */
  
  package dr.evomodel.MSSD;
  
  import dr.evolution.alignment.AscertainedSitePatterns;
  import dr.evolution.alignment.PatternList;
  import dr.evolution.datatype.MutationDeathType;
  import dr.evolution.tree.NodeRef;
  import dr.evomodel.branchratemodel.BranchRateModel;
  import dr.evomodel.branchratemodel.DefaultBranchRateModel;
  import dr.evomodel.sitemodel.SiteRateModel;
  import dr.evomodel.tree.TreeModel;
  import dr.evomodel.treelikelihood.LikelihoodCore;
  import dr.evomodel.treelikelihood.LikelihoodPartialsProvider;
  import dr.evomodel.treelikelihood.ScaleFactorsHelper;
  import dr.inference.model.AbstractModel;
  import dr.inference.model.Model;
  import dr.inference.model.Parameter;
  import dr.inference.model.Variable;
  import dr.math.GammaFunction;
  
  /**
   * Package: AbstractObservationProcess
   * Description:
   * <p/>
   * <p/>
   * Created by
   * Alexander V. Alekseyenko (alexander.alekseyenko@gmail.com)
   * Date: Feb 19, 2008
   * Time: 12:41:01 PM
   */
  abstract public class AbstractObservationProcess extends AbstractModel {
      protected boolean[] nodePatternInclusion;
      protected boolean[] storedNodePatternInclusion;
      protected double[] cumLike;
      protected double[] nodePartials;
      protected double[] nodeLikelihoods;
      protected int nodeCount;
      protected int patternCount;
      protected int stateCount;
      protected TreeModel treeModel;
      protected PatternList patterns;
      protected double[] patternWeights;
      protected Parameter mu;
      protected Parameter lam;
  
      // update control variables
      protected boolean weightKnown;
      protected double logTreeWeight;
      protected double storedLogTreeWeight;
      private double gammaNorm;
      private double totalPatterns;
      protected MutationDeathType dataType;
      protected int deathState;
      protected SiteRateModel siteModel;
      private double logN;
      protected boolean nodePatternInclusionKnown = false;
      BranchRateModel branchRateModel;
  
      public AbstractObservationProcess(String Name, TreeModel treeModel, PatternList patterns, SiteRateModel siteModel,
                                        BranchRateModel branchRateModel, Parameter mu, Parameter lam) {
          super(Name);
          this.treeModel = treeModel;
          this.patterns = patterns;
          this.mu = mu;
          this.lam = lam;
          this.siteModel = siteModel;
          if (branchRateModel != null) {
              this.branchRateModel = branchRateModel;
          } else {
              this.branchRateModel = new DefaultBranchRateModel();
          }
          addModel(treeModel);
          addModel(siteModel);
          addModel(this.branchRateModel);
          addVariable(mu);
          addVariable(lam);
  
          nodeCount = treeModel.getNodeCount();
          stateCount = patterns.getDataType().getStateCount();
          this.patterns = patterns;
          patternCount = patterns.getPatternCount();
          patternWeights = patterns.getPatternWeights();
          totalPatterns = 0;
          for (int i = 0; i < patternCount; ++i) {
              totalPatterns += patternWeights[i];
          }
          logN = Math.log(totalPatterns);
  
          gammaNorm = -GammaFunction.lnGamma(totalPatterns + 1);
  
          dataType = (MutationDeathType) patterns.getDataType();
          this.deathState = dataType.DEATHSTATE;
          setNodePatternInclusion();
          cumLike = new double[patternCount];
          nodeLikelihoods = new double[patternCount];
          weightKnown = false;
      }
  
  //    public Parameter getMuParameter() {
  //        return mu;
  //    }
  //
  //    public Parameter getLamParameter() {
  //        return lam;
  //    }
  
      private double calculateSiteLogLikelihood(int site, double[] partials, double[] frequencies) {
          int v = site * stateCount;
          double sum = 0.0;
          for (int i = 0; i < stateCount; i++) {
              sum += frequencies[i] * partials[v + i];
          }
          return Math.log(sum);
      }
  
  
      private void calculateNodePatternLikelihood(int nodeIndex,
                                                  double[] freqs,
                                                  LikelihoodCore likelihoodCore,
                                                  double averageRate,
                                                  double[] cumLike) {
          // get partials for node nodeIndex
          likelihoodCore.getPartials(nodeIndex, nodePartials); // MAS
              /*
                  multiply the partials by equilibrium probs
                      this part could be optimized by first summing
                      and then multiplying by equilibrium probs
              */
          double prob = Math.log(getNodeSurvivalProbability(nodeIndex, averageRate));
  
          for (int j = 0; j < patternCount; ++j) {
              if (nodePatternInclusion[nodeIndex * patternCount + j]) {
                  cumLike[j] += Math.exp(calculateSiteLogLikelihood(j, nodePartials, freqs) + prob);
              }
          }
      }
  
      private double accumulateCorrectedLikelihoods(double[] cumLike, double ascertainmentCorrection,
                                                    double[] patterWeights) {
          double logL = 0;
          for (int j = 0; j < patternCount; ++j) {
              logL += Math.log(cumLike[j] / ascertainmentCorrection) * patternWeights[j];
          }
          return logL;
      }
  
      public final double nodePatternLikelihood(double[] freqs, LikelihoodPartialsProvider likelihoodCore,
                                                ScaleFactorsHelper scaleFactorsHelper) {
          int i, j;
          double logL = gammaNorm;
  
          double birthRate = lam.getParameterValue(0);
          double logProb;
          if (!nodePatternInclusionKnown)
              setNodePatternInclusion();
          if (nodePartials == null) {
              nodePartials = new double[patternCount * stateCount];
          }
  
          double averageRate = getAverageRate();
  
          for (j = 0; j < patternCount; ++j) cumLike[j] = 0;
  
          for (i = 0; i < nodeCount; ++i) {
              // get partials for node i
              likelihoodCore.getPartials(i, nodePartials);
              scaleFactorsHelper.rescalePartials(i, nodePartials);
              /*
                  multiply the partials by equilibrium probs
                      this part could be optimized by first summing
                      and then multiplying by equilibrium probs
              */
  //            likelihoodCore.calculateLogLikelihoods(nodePartials, freqs, nodeLikelihoods);   // MAS Removed
              logProb = Math.log(getNodeSurvivalProbability(i, averageRate));
  
              for (j = 0; j < patternCount; ++j) {
                  if (nodePatternInclusion[i * patternCount + j]) {
  //                    cumLike[j] += Math.exp(nodeLikelihoods[j] + logProb);  // MAS Replaced with line below
                      cumLike[j] += Math.exp(calculateSiteLogLikelihood(j, nodePartials, freqs)
                              + logProb);
                  }
              }
          }
  
          double ascertainmentCorrection = getAscertainmentCorrection(cumLike);
  //        System.err.println("AscertainmentCorrection: "+ascertainmentCorrection);
  
          for (j = 0; j < patternCount; ++j) {
              logL += Math.log(cumLike[j] / ascertainmentCorrection) * patternWeights[j];
          }
  
          double deathRate = mu.getParameterValue(0);
  
  