Commit b0ce573b authored by Claudio Gheller's avatar Claudio Gheller
Browse files

Makefile and gitignore extended to inverse-imaging, inverse-imaging updated

parent 36a65593

.gitignore

+1 −0
Original line number Diff line number Diff line
@@ -7,3 +7,4 @@ w-stackingCfftw 
w-stackingfftw

w-stackingCfftw_serial

w-stackingfftw_serial

inverse-imaging

Makefile

+1 −0
Original line number Diff line number Diff line
@@ -72,6 +72,7 @@ serial_cuda: 


mpi: $(COBJ)

	$(MPICC) $(OPTIMIZE) -o w-stackingCfftw   $^  $(CFLAGS) $(LIBS)

	$(MPICC) $(OPTIMIZE) -o inverse-imaging inverse-imaging.c w-stacking.c $(CFLAGS) $(LIBS)



mpi_cuda:

	$(NVCC)   $(NVFLAGS) -c w-stacking.cu phase_correction.cu $(NVLIB)

inverse-imaging.c

+28 −62
Original line number Diff line number Diff line
@@ -79,8 +79,8 @@ int main(int argc, char * argv[]) 
	char timingfile[FILENAMELENGTH] = "timings.dat";



        // Image related files

	char imagepath[900] = "./"

	imagename[FILENAMELENGTH] = "image_name.bin"

	char imagepath[900] = "./";

	char imagename[FILENAMELENGTH] = "image_name.bin";



	// Visibilities related variables

	double * uu;
@@ -115,7 +115,7 @@ int main(int argc, char * argv[]) 
	int yaxis;



	// Number of planes in the w direction

	int num_w_planes = 8;

	int num_w_planes = 1;



	// Size of the convoutional kernel support

	int w_support = 7;
@@ -278,6 +278,8 @@ if(rank == 0){ 
	visreal = (float*) calloc(Nvis,sizeof(float));

	visimg = (float*) calloc(Nvis,sizeof(float));



        // reading baselines



  if(rank == 0)printf("READING DATA\n");

	// Read data

	strcpy(filename,datapath);
@@ -316,9 +318,29 @@ if(rank == 0){ 
	setup_time = ((double) (end - start)) / CLOCKS_PER_SEC;

	setup_time1 = (finish.tv_sec - begin.tv_sec);

	setup_time1 += (finish.tv_nsec - begin.tv_nsec) / 1000000000.0;

        // baselines read 



        // define image variable

	double* image_real = (double*) calloc(xaxis*yaxis,sizeof(double));



        // reading image

	strcpy(filename,imagepath);

        strcat(filename,imagename);

        printf("Reading Image %s\n",filename);



	// all MPI tasks read together: parallel filesystem required

        pFilereal = fopen (filename,"rb");	

        long global_index = rank*(xaxis*yaxis)*sizeof(double);

        fseek(pFilereal, global_index, SEEK_SET);

	fwrite(image_real, xaxis*yaxis, sizeof(double), pFilereal);

        fclose(pFilereal);

	// image read

	

  if(rank == 0)printf("GRIDDING DATA\n");

  if(rank == 0)printf("FFT TRANSFORMING (from Real to Complex Fourier space)\n");

         





  if(rank == 0)printf("DEGRIDDING DATA\n");



	// Create histograms and linked lists

        clock_gettime(CLOCK_MONOTONIC, &begin);
@@ -919,62 +941,6 @@ if(rank == 0){ 


	fftw_free(fftwgrid);



	// Phase correction

        clock_gettime(CLOCK_MONOTONIC, &begin);

        start = clock();

	if(rank == 0)printf("PHASE CORRECTION\n");

        double* image_real = (double*) calloc(xaxis*yaxis,sizeof(double));	

        double* image_imag = (double*) calloc(xaxis*yaxis,sizeof(double));	



        phase_correction(gridss,image_real,image_imag,xaxis,yaxis,num_w_planes,grid_size_x,grid_size_y,resolution,wmin,wmax,num_threads);



        end = clock();

        clock_gettime(CLOCK_MONOTONIC, &finish);

        phase_time = ((double) (end - start)) / CLOCKS_PER_SEC;

        phase_time1 = (finish.tv_sec - begin.tv_sec);

        phase_time1 += (finish.tv_nsec - begin.tv_nsec) / 1000000000.0;

#ifdef WRITE_IMAGE



        if(rank == 0)

        {

            pFilereal = fopen (fftfile2,"wb");

            pFileimg = fopen (fftfile3,"wb");

            fclose(pFilereal);

            fclose(pFileimg);

        }

	#ifdef USE_MPI

	MPI_Barrier(MPI_COMM_WORLD);

        #endif

        if(rank == 0)printf("WRITING IMAGE\n");

	for (int isector=0; isector<size; isector++)

	{

	    #ifdef USE_MPI

	    MPI_Barrier(MPI_COMM_WORLD);

            #endif

	    if(isector == rank)

            {

	       printf("%d writing\n",isector);

               pFilereal = fopen (fftfile2,"ab");

               pFileimg = fopen (fftfile3,"ab");



	       long global_index = isector*(xaxis*yaxis)*sizeof(double);



               fseek(pFilereal, global_index, SEEK_SET);

               fwrite(image_real, xaxis*yaxis, sizeof(double), pFilereal);

               fseek(pFileimg, global_index, SEEK_SET);

               fwrite(image_imag, xaxis*yaxis, sizeof(double), pFileimg);



               fclose(pFilereal);

               fclose(pFileimg);

	    }

	}

	#ifdef USE_MPI

	MPI_Barrier(MPI_COMM_WORLD);

        #endif



#endif //WRITE_IMAGE





#endif //USE_FFTW



	end = clock();