 			if (writeIndex+1<10) ss << "0"; 
 			ss << writeIndex+1 << "psi.fits";
 			writeFitsImg(ss.str().c_str(), dimensions, 
-				pixelsPerImage, psiImages[writeIndex]);
+				pixelsPerImage, phiImages[writeIndex]);
 			ss.str("");
 			ss.clear();
 		}
 	}
 	std::cout << "Done.\n";
 
-	/* Write psiImagess file */
+	/* Write results to file */
 	// std::cout needs to be rerouted to output to console after this...
 	std::freopen("results.txt", "w", stdout);
 	std::cout << "# t0_x t0_y theta_par theta_perp v_x v_y likelihood est_flux\n";

	CUDA_CHECK_RETURN(cudaMemcpy(deviceImgTimes, imageTimes,
			sizeof(float)*imageCount, cudaMemcpyHostToDevice));

	// Copy over psi images one at a time
	for (int i=0; i<imageCount; ++i)
	{
		CUDA_CHECK_RETURN(cudaMemcpy(deviceImages+pixelsPerImage*i, psiImages[i],
			sizeof(float)*pixelsPerImage, cudaMemcpyHostToDevice));
	}

	dim3 blocks(dimensions[0]/32+1,dimensions[1]/32+1);
	dim3 threads(32,32);

	// Launch Search
	searchImages<<<blocks, threads>>> (dimensions[0], dimensions[1], imageCount, deviceImages,
				trajCount, deviceTests, deviceSearchResults, backgroundLevel, padding);

	// Read back psiImagess
	CUDA_CHECK_RETURN(cudaMemcpy(bestTrajects, deviceSearchResults,
				sizeof(trajectory)*pixelsPerImage, cudaMemcpyDeviceToHost));

	CUDA_CHECK_RETURN(cudaFree(deviceTests));
	CUDA_CHECK_RETURN(cudaFree(deviceImgTimes));
	CUDA_CHECK_RETURN(cudaFree(deviceSearchResults));
	CUDA_CHECK_RETURN(cudaFree(deviceImages));

	
	// Sort psiImagess by likelihood
	qsort(bestTrajects, pixelsPerImage, sizeof(trajectory), compareTrajectory);
	if (debug)
	{

	delete[] psiImages;
	delete[] phiImages;
	delete[] imageTimes;
	delete[] interleavedPsiPhi;		

	delete[] angles;
	delete[] velocities;
	delete[] trajectoriesToSearch;	
