private static double[][] normalizeDensity(double[][] unNormDensity){
		
		// The number of gridlines across the density domain
		int numGridLines = getNumGridlines();
		
		int iter = 0;                        // How many normalization cycles have been performed
		double threshold = Math.pow(10, -8); // The maximum acceptable error
		double densityDomainSize = Settings.getMaximumRange() - Settings.getMinimumRange();
		double[][] normDens = unNormDensity;
		
		while (iter < 1000) {
			
			// Zero negative points
			for (int i1 = 0; i1 < numGridLines; i1++) {
				
				for (int i2 = 0; i2 < numGridLines; i2++) {
					
					if (normDens[i1][i2] < 0.0) {
						normDens[i1][i2] = 0.0;
					}
				}
			}
			
			// Sum over probability density over interval
			double integralSum = 0.0;
			for (int i1 = 0; i1 < numGridLines; i1++) {
				for (int i2 = 0; i2 < numGridLines; i2++) {
					integralSum += normDens[i1][i2] * Math.pow(Settings.discretization, 2);
				}
			}
			
			// Return if error is under threshold
			if (Math.abs(integralSum - 1) < threshold) {
				return (normDens);
			}
			
			// Modify density so that it integrates to 1
			double normalizeConstant = (integralSum - 1) / Math.pow(densityDomainSize, 2);
			for (int i1 = 0; i1 < numGridLines; i1++) {
				for (int i2 = 0; i2 < numGridLines; i2++) {
					normDens[i1][i2] -= normalizeConstant;
				}
			}
			
			iter++;
		}
		
		// Settle for the current approximation
		return normDens;
	} //end normalizeDensity