Commit 99500eb4 authored by Emanuele De Rubeis's avatar Emanuele De Rubeis
Browse files

Full working parallel fits writing with CFITSIO

parent 736c2f2a

w-stacking-fftw.c

+45 −6
Original line number Diff line number Diff line
@@ -986,13 +986,52 @@ if(rank == 0){ 
        long * fpixel = (long *) malloc(sizeof(long)*naxis);

        long * lpixel = (long *) malloc(sizeof(long)*naxis);



        #ifdef USE_MPI

        MPI_Barrier(MPI_COMM_WORLD);

        #endif

	#ifdef PARALLEL_FITS

        #ifdef FITSIO



        fpixel[0] = 1;

        fpixel[1] = rank*yaxis+1;

        lpixel[0] = xaxis;

        lpixel[1] = (rank+1)*yaxis;



        status = 0;

        fits_open_image(&fptreal, testfitsreal, READWRITE, &status);

        fits_write_subset(fptreal, TDOUBLE, fpixel, lpixel, image_real, &status);

        fits_close_file(fptreal, &status);



        status = 0;

        fits_open_image(&fptrimg, testfitsimag, READWRITE, &status);

        fits_write_subset(fptrimg, TDOUBLE, fpixel, lpixel, image_imag, &status);

        fits_close_file(fptrimg, &status);



        #endif



        pFilereal = fopen (fftfile2,"ab");

        pFileimg = fopen (fftfile3,"ab");



	long global_index = rank*(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

        #else

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

	{

	    #ifdef USE_MPI

	    MPI_Barrier(MPI_COMM_WORLD);

            #endif

	    if(isector == rank)

            //if(rank == 0)

	    {

               #ifdef FITSIO
@@ -1004,8 +1043,6 @@ if(rank == 0){ 
               fpixel[1] = isector*yaxis+1;

               lpixel[0] = xaxis;

               lpixel[1] = (isector+1)*yaxis;

               //printf("fpixel %d, %d\n", fpixel[0], fpixel[1]);

	       //printf("lpixel %d, %d\n", lpixel[0], lpixel[1]);

               

               status = 0;

	       fits_open_image(&fptreal, testfitsreal, READWRITE, &status);
@@ -1033,10 +1070,12 @@ if(rank == 0){ 
               fclose(pFileimg);

	    }

	}

	#endif

	#ifdef USE_MPI

	MPI_Barrier(MPI_COMM_WORLD);

        #endif





#endif //WRITE_IMAGE