 			}
 			
 			// For each iteration, run belief propagation
+			dim3 lvlGridSize(GINF_DIM(w / 2 + 1, blockSize.x), GINF_DIM(h, blockSize.y));
 			for (int i = 0; i < numItersPerLevel; i++) {
-				gpuBpKernel<<<gridSize, blockSize>>>(i, dGrid, w, h, dMsgs[t], dFrom[t], dData[t]);
+				gpuBpCollectKernel<<<gridSize, blockSize>>>(dGrid, w, h, dMsgs[t], dFrom[t], dData[t]);
+
+				// For each direction, send a message
+				for (int d = 0; d < GINF_NUM_DIR; d++) {
+					gpuBpKernel<<<lvlGridSize, blockSize>>>(i, dGrid, w, h, dMsgs[t], dFrom[t], dData[t], d);
+				}
 			}
 		}
 		
 		// Finally get the beliefs
+		gpuBpCollectKernel<<<gridSize, blockSize>>>(dGrid, grid->getWidth(), grid->getHeight(), dMsgs[0], dFrom[0], dData[0]);
 		gpuBpGetBeliefKernel<<<gridSize, blockSize>>>(dGrid, dMsgs[0], dFrom[0], dResult);
 		
 		// Transfer the results back to CPU

			}
			
			// Write the result
			result[dim.idx(x, y)] = minLabel;
		}
	}
	

	template void gpuDecodeHbp<int>(Grid<int>*, int, int, Matrix<int>*);
	
	template <typename T>
	void gpuDecodeBp(Grid<T> *grid, int numIters, Matrix<int> *result) {
		gpuDecodeHbp(grid, 1, numIters, result);
	}
	
	template void gpuDecodeBp<float>(Grid<float>*, int, Matrix<int>*);
	template void gpuDecodeBp<int>(Grid<int>*, int, Matrix<int>*);
}
