diff --git a/w-stacking-fftw.c b/w-stacking-fftw.c
index 26e8018c14dac343cbb20ffefcc188bb585532a1..fc8eb047c6702cda089ad120fd175ca3edb68b0b 100644
--- a/w-stacking-fftw.c
+++ b/w-stacking-fftw.c
@@ -39,33 +39,32 @@ void Push(struct sectorlist** headRef, long data) {
// Main Code
int main(int argc, char * argv[])
{
+ // Main MPI parameters
int rank;
int size;
+ // Define main filenames
FILE * pFile;
FILE * pFile1;
FILE * pFilereal;
FILE * pFileimg;
+ // Global filename to be composed
char filename[1000];
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/gridding/old/data/gauss2_t201806301100_SBL180.binMS/";
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/gridding/newgauss2noconj_t201806301100_SBL180.binMS/";
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/gridding/gauss1_t201806301100_SBL180.binMS/";
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/gridding/newgauss4_t201806301100_SBL180.binMS/";
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/gridding/gauss4_t201806301100_SBL180.binMS/";
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/gridding/hba-8hrs_t201806301100_SBH255i-test.binMS/";
- //
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/gridding/hba-8hrs_gauss4new.binMS/";
- //char datapath[900] = "/m100_scratch/userexternal/cgheller/Lofar/Observations/L798046_SB244_uv.uncorr_130B27932t_146MHz.pre-cal.binMS/";
+
+ // MS paths
char datapath[900];
char datapath_multi[NFILES][900];
+ // Bin MS files
char ufile[30] = "ucoord.bin";
- char vfile[30] = "vcoord.bin";
- char wfile[30] = "wcoord.bin";
- char weightsfile[30] = "weights.bin";
- char visrealfile[30] = "visibilities_real.bin";
- char visimgfile[30] = "visibilities_img.bin";
- char metafile[30] = "meta.txt";
+ char vfile[30] = "vcoord.bin";
+ char wfile[30] = "wcoord.bin";
+ char weightsfile[30] = "weights.bin";
+ char visrealfile[30] = "visibilities_real.bin";
+ char visimgfile[30] = "visibilities_img.bin";
+ char metafile[30] = "meta.txt";
+
+ // Mesh related files
char outfile[30] = "grid.txt";
char outfile1[30] = "coords.txt";
char outfile2[30] = "grid_real.bin";
@@ -78,6 +77,7 @@ int main(int argc, char * argv[])
char srank[4];
char timingfile[30] = "timings.dat";
+ // Visibilities related variables
double * uu;
double * vv;
double * ww;
@@ -100,58 +100,74 @@ int main(int argc, char * argv[])
double wmax,wmax0;
double resolution;
- // MESH SIZE
+ // Mesh related parameters: global size
int grid_size_x = 2048;
int grid_size_y = 2048;
- int local_grid_size_x;// = 8;
- int local_grid_size_y;// = 8;
+ // Split Mesh size (auto-calculated)
+ int local_grid_size_x;
+ int local_grid_size_y;
int xaxis;
int yaxis;
+
+ // Number of planes in the w direction
int num_w_planes = 8;
- // DAV: the corresponding KernelLen is calculated within the wstack function. It can be anyway hardcoded for optimization
+ // Size of the convoutional kernel support
int w_support = 7;
- int num_threads;// = 4;
+
+ // Number of OpenMP threads (input parameter)
+ int num_threads;
+
+ // Resolution
double dx = 1.0/(double)grid_size_x;
double dw = 1.0/(double)num_w_planes;
+ // Half support size
double w_supporth = (double)((w_support-1)/2)*dx;
+ // Internal profiling parameters
clock_t start, end, start0, startk, endk;
double setup_time, process_time, mpi_time, fftw_time, tot_time, kernel_time, reduce_time, compose_time, phase_time;
double setup_time1, process_time1, mpi_time1, fftw_time1, tot_time1, kernel_time1, reduce_time1, compose_time1, phase_time1;
double writetime, writetime1;
-
struct timespec begin, finish, begin0, begink, finishk;
double elapsed;
- long nsectors;
- /* GT get nymber of threads exit if not given */
- if(argc == 1) {
- fprintf(stderr, "Usage: %s number_of_OMP_Threads \n", argv[0]);
- exit(1);
- }
- // Set the number of OpenMP threads
- num_threads = atoi(argv[1]);
-
- if ( num_threads == 0 )
- {
- fprintf(stderr, "Wrong parameter: %s\n\n", argv[1]);
- fprintf(stderr, "Usage: %s number_of_OMP_Threads \n", argv[0]);
- exit(1);
- }
+ // Number of sectors in which the mesh is divided along the v direction
+ // IF nsectors < 0, nsectors = NUMBER OF MPI RANKS
+ long nsectors = NUM_OF_SECTORS;
+
+
+ /* Get number of threads, exit if not given */
+ if (argc == 1) {
+ fprintf(stderr, "Usage: %s number_of_OMP_Threads \n", argv[0]);
+ exit(1);
+ }
+ // Set the number of OpenMP threads
+ num_threads = atoi(argv[1]);
+
+ // Number of threads <= 0 not acceptable
+ if (num_threads <= 0) {
+ fprintf(stderr, "Wrong parameter: %s\n\n", argv[1]);
+ fprintf(stderr, "Usage: %s number_of_OMP_Threads \n", argv[0]);
+ exit(1);
+ }
+
+ // Internal profiling
clock_gettime(CLOCK_MONOTONIC, &begin0);
start0 = clock();
- // Intialize MPI environment
#ifdef USE_MPI
+ // Intialize MPI environment
MPI_Init(&argc,&argv);
MPI_Comm_rank(MPI_COMM_WORLD, &rank);
MPI_Comm_size(MPI_COMM_WORLD, &size);
if(rank == 0)printf("Running with %d MPI tasks\n",size);
#ifdef USE_FFTW
+ // Initialize FFTW environment
fftw_mpi_init();
#endif
#else
+ // If MPI is not used the two parameters are set for consistency
rank = 0;
size = 1;
#endif
@@ -173,7 +189,6 @@ if(rank == 0){
// set the local size of the image
local_grid_size_x = grid_size_x;
- nsectors = NUM_OF_SECTORS;
if (nsectors < 0) nsectors = size;
local_grid_size_y = grid_size_y/nsectors;
//nsectors = size;
@@ -185,7 +200,7 @@ if(rank == 0){
clock_gettime(CLOCK_MONOTONIC, &begin);
start = clock();
- // INPUT FILES (only the first ndatasets entries are used)
+ // INPUT FILES (only the first ndatasets entries are used)
int ndatasets = 1;
//strcpy(datapath_multi[0],"data/newgauss2noconj_t201806301100_SBL180.binMS/");
//strcpy(datapath_multi[0],"/m100_scratch/userexternal/cgheller/gridding/newgauss4_t201806301100_SBL180.binMS/");