public static double erlangBprobabilityNIR(double totalTraffic, double totalResources) {

		double prob = 0.0D; // blocking probability

		// numerator
		double lognumer = totalResources * Math.log(totalTraffic) - totalTraffic;
		// denominator (incomplete Gamma Function)
		double oneplustr = 1.0D + totalResources;
		double crigf = Stat.complementaryRegularisedGammaFunction(oneplustr, totalTraffic);
		if (crigf != crigf || crigf == 0.0) {
			int low = (int) Math.floor(totalResources) - 2;
			if (low < 0) low = 0;
			int ni = 6;
			int[] tr = new int[ni];
			double[] trd = new double[ni];
			double[] pp = new double[ni];
			for (int i = 0; i < ni; i++) {
				tr[i] = low + i;
				trd[i] = tr[i];
				pp[i] = Stat.erlangBprobability(totalTraffic, tr[i]);
			}
			CubicSpline cs = new CubicSpline(trd, pp);
			prob = cs.interpolate(totalResources);
		} else {
			double logdenom = Math.log(crigf) + Stat.logGammaFunction(oneplustr);
			prob = Math.exp(lognumer - logdenom);
		}
		return prob;
	}