private static double[][] addWaveDensity(double[][] density) {
		
		int numGridLines = getNumGridlines();
			
		// Loop through for each resolution level
		for (double j = Settings.startLevel; j <= Settings.stopLevel; j++) {
			
			int waveInd = (int) Math.floor(j - Settings.startLevel);
			ArrayList<ArrayList<Double>> jPhis = wavePhisHere.get(waveInd);
			ArrayList<ArrayList<Double>> jPsis = wavePsisHere.get(waveInd);
			
			double i1 = Settings.getMinimumRange();
			// Calculate un-normalized density for each point in domain, looping across X1
			for (int x1Ind = 0; x1Ind < numGridLines; x1Ind++)
				{
				i1 += Settings.discretization;
				int[] i1RelevantIndices = findRelevantKIndices(i1, j);
				int k1Max = i1RelevantIndices[1];
				int k1Min = i1RelevantIndices[0];					
				
				// Cycle through relevant X1 translates for the line
				for (int k1Ind = k1Min; k1Ind < k1Max; k1Ind++) {
					
					double phi1Here = jPhis.get(x1Ind).get(k1Ind - k1Min);
					double psi1Here = jPsis.get(x1Ind).get(k1Ind - k1Min);
					
					// Loop across X2 dimension
					double i2 = Settings.getMinimumRange();
					for (int x2Ind = 0; x2Ind < numGridLines; x2Ind++)
					{
						i2 += Settings.discretization;
						int[] i2RelevantIndices = findRelevantKIndices(i2, j);
						int k2Max = i2RelevantIndices[1];
						int k2Min = i2RelevantIndices[0];
											
						// Cycle through relevant X2 translates
						for (int k2Ind = k2Min; k2Ind < k2Max; k2Ind++) {
							
							double phi2Here = jPhis.get(x2Ind).get(k2Ind - k2Min);
							double psi2Here = jPsis.get(x2Ind).get(k2Ind - k2Min);
							
							density[x1Ind][x2Ind] += Transform.psiPsiCoefficients.get(waveInd)[k1Ind][k2Ind] *
									psi1Here*psi2Here;
							density[x1Ind][x2Ind] += Transform.psiPhiCoefficients.get(waveInd)[k1Ind][k2Ind] *
									psi1Here*phi2Here;
							density[x1Ind][x2Ind] += Transform.phiPsiCoefficients.get(waveInd)[k1Ind][k2Ind] *
									phi1Here*psi2Here;
						} // End cycling relevant X2 translates
					} // End cycling across X2 gridlines
				} // End cycling relevant X1 translates
			} // End cycling across X1 gridlines
		}
		
		return density;
	}