public static double[] determineMarkThresholdsFromBinnedDataArray(int[][][] grid, boolean[] bpresent, 
								      double dpoissonthresh,double dfoldthresh,boolean bcontainsthresh)
    {
       double dcumthreshold = 1-dpoissonthresh;

       double[] sumtags =null;
       double[] thresholds = null;
       int nummarks= grid[0][0].length;
       String szcell;

       int ntotallocs = 0;

       sumtags = new double[nummarks];
       thresholds = new double[nummarks];

       //computes for those chromosomes considered present
       //the total number of locations in ntotallocs
       //and the total signal for the mark sumtags[nmark]
       for (int nchrom = 0; nchrom < grid.length; nchrom++)
       {
	   if (bpresent[nchrom])
	   {
	      int[][] grid_nchrom = grid[nchrom];
	      for (int nbin = 0; nbin < grid_nchrom.length; nbin++)
	      {
	         int[] grid_nchrom_nbin = grid[nchrom][nbin];
                 for (int nmark = 0; nmark < nummarks; nmark++)
                 {		    
	            sumtags[nmark] += grid_nchrom_nbin[nmark];
		 }	      
	      }
	      ntotallocs += grid_nchrom.length;
	   }
        }	   

       for (int nj = 0; nj < sumtags.length; nj++)
       {
	   //computes expected number of reads in a bin
          double dlambda = sumtags[nj]/(double) ntotallocs;
	  double dcum = 0;
	  short nthresh = 0;

	  double dlogfactorial = 0;
	  while (dcum <= dcumthreshold)
          {
	      //lambda^x times e^(-lambda) / (x!)
	      //taking the exponential of the log term
	      double dprob = Math.exp(Math.log(dlambda)*nthresh-dlambda-dlogfactorial);
	      dcum += dprob;
	      nthresh++;
	      dlogfactorial += Math.log(nthresh);
	  }

	  if (bcontainsthresh)
	  {
	      //want to ensure  dpoissonthresh is included with the cutoff
             nthresh--;
	  }

	  //threshold set to greater of poisson based and fold threshold
          thresholds[nj] = Math.max(Math.max((int) Math.ceil(dfoldthresh*dlambda),nthresh),1);
	  if (ChromHMM.BVERBOSE)
	  {
             System.out.println("Threshold\t"+nj+"\t"+thresholds[nj]+"\t"+sumtags[nj]+"\t"+((int) Math.ceil(dfoldthresh*dlambda))+"\t"+ntotallocs);
	  }
       }
       return thresholds;
    }