QDPIO::cout << "Before Layout" << endl;
int main(int argc, char **argv){  // Put the machine into a known state  QDP_initialize(&argc, &argv);  // Input parameter structure  Spectrum_input_t  input;  QDPIO::cout << "Before DATA" << endl;  // Instantiate xml reader for DATA  XMLReader xml_in("DATA");  QDPIO::cout << "Before read spectrum" << endl;  // Read data  read(xml_in, "/spectrum_w", input);  QDPIO::cout << "Before Layout" << endl;  // Specify lattice size, shape, etc.  Layout::setLattSize(input.param.nrow);  Layout::create();  /*   * Sanity checks   */  // Check that only one type of source smearing is specified  // Must match how the stored propagator was smeared  if (input.param.Pt_src && input.param.Sl_src) {    QDPIO::cerr << "Error if Pt_src and Sl_src are both set." << endl;    QDPIO::cerr << "Choose the one which matches the stored propagators." << endl;    QDP_abort(1);  }  // Figure out what to do about boundary conditions  // GTF HACK: only allow periodic boundary conditions  for (int i=0; i<Nd; ++i) {    if (input.param.boundary[i] != 1) {      QDPIO::cerr << "Only periodic boundary conditions supported." << endl;      QDPIO::cerr << "  boundary[" << i << "] = " << input.param.boundary[i] << endl;      QDP_abort(1);    }  }  for (int i=0; i<Nd; ++i) {    if (input.param.t_srce[i] < 0 || input.param.t_srce[i] >= input.param.nrow[i]) {      QDPIO::cerr << "Quark propagator source coordinate incorrect." << endl;      QDPIO::cerr << "t_srce[" << i << "] = " << input.param.t_srce[i] << endl;      QDP_abort(1);    }  }  QDPIO::cout << endl << "     Gauge group: SU(" << Nc << ")" << endl;  QDPIO::cout << "     volume: " << input.param.nrow[0];  for (int i=1; i<Nd; ++i) {    QDPIO::cout << " x " << input.param.nrow[i];  }  QDPIO::cout << endl;  // Read in the configuration along with relevant information.  multi1d<LatticeColorMatrix> u(Nd);  XMLReader gauge_xml;  switch (input.param.cfg_type)   {  case CFG_TYPE_SZIN :    readSzin(gauge_xml, u, input.cfg.cfg_file);    break;  default :    QDP_error_exit("Configuration type is unsupported.");  }  // Check if the gauge field configuration is unitarized  unitarityCheck(u);  // Instantiate XML writer for XMLDAT  XMLFileWriter xml_out("XMLDAT");  push(xml_out, "spectrum_w");  // Write out the input  write(xml_out, "Input", xml_in);  // Write out the config info  write(xml_out, "Config_info", gauge_xml);  push(xml_out, "Output_version");  write(xml_out, "out_version", 6);  pop(xml_out);  xml_out.flush();  // First calculate some gauge invariant observables just for info.  Double w_plaq, s_plaq, t_plaq, link;  multi1d<DComplex> pollp(Nd);  MesPlq(u, w_plaq, s_plaq, t_plaq, link);  for(int mu = 0; mu < Nd; ++mu)    polylp(u, pollp[mu], mu);  push(xml_out, "Observables");  Write(xml_out, w_plaq);  Write(xml_out, s_plaq);  Write(xml_out, t_plaq);  Write(xml_out, link);  Write(xml_out, pollp);  pop(xml_out);  xml_out.flush();  // Initialize the slow Fourier transform phases  SftMom phases(input.param.mom2_max, input.param.avg_equiv_mom, input.param.j_decay);  // Keep an array of all the xml output buffers  XMLArrayWriter xml_array(xml_out,input.param.numKappa);  push(xml_array, "Wilson_hadron_measurements");  // Flags  int t0      = input.param.t_srce[input.param.j_decay];  int bc_spec = input.param.boundary[input.param.j_decay];  // Now loop over the various fermion masses  for (int loop=0; loop < input.param.numKappa; ++loop)  {    // Read the quark propagator    LatticePropagator quark_propagator;    {      XMLReader prop_xml;      stringstream prop_file;      prop_file << "propagator_" << loop;      readSzinQprop(prop_xml, quark_propagator, prop_file.str());    }    push(xml_array);         // next array element - name auto-written    Write(xml_array, loop);    write(xml_array, "Kappa_mes", input.param.Kappa[loop]);    write(xml_array, "t_srce", input.param.t_srce);    // Do the mesons first    if (input.param.MesonP)     {      // Construct {Point|Shell}-Point mesons, if desired      if (input.param.Pt_snk)       {	if (input.param.Pt_src)	  mesons(quark_propagator, quark_propagator, phases, t0,		 xml_array, "Point_Point_Wilson_Mesons");        	if (input.param.Sl_src)	  mesons(quark_propagator, quark_propagator, phases, t0,		 xml_array, "Shell_Point_Wilson_Mesons");      } // end if (Pt_snk)      // Convolute the quark propagator with the sink smearing function.      // Make a copy of the quark propagator and then overwrite it with      // the convolution.       if (input.param.Sl_snk)       {	LatticePropagator quark_prop_smr;	quark_prop_smr = quark_propagator;	sink_smear2(u, quark_prop_smr, input.param.Wvf_kind, input.param.wvf_param[loop],		    input.param.WvfIntPar[loop], input.param.j_decay);	if (input.param.Pt_src)	  mesons(quark_prop_smr, quark_prop_smr, phases, t0,		 xml_array, "Point_Shell_Wilson_Mesons");	if (input.param.Sl_src)	  mesons(quark_prop_smr, quark_prop_smr, phases, t0,		 xml_array, "Shell_Shell_Wilson_Mesons");      } // end if (Sl_snk)    } // end if (MesonP)    // Do the currents next    if (input.param.CurrentP)     {      // Construct the rho vector-current and the pion axial current divergence      if (input.param.Pt_src)	curcor2(u, quark_propagator, quark_propagator, phases, 		t0, input.param.j_decay, 4,		xml_array, "Point_Point_Meson_Currents");              if (input.param.Sl_src)	curcor2(u, quark_propagator, quark_propagator, phases, 		t0, input.param.j_decay, 4,		xml_array, "Shell_Point_Meson_Currents");    } // end if (CurrentP)    // Do the baryons    if (input.param.BaryonP)     {      // Construct {Point|Shell}-Point mesons, if desired      if (input.param.Pt_snk)       {	if (input.param.Pt_src)	  baryon(quark_propagator, phases, 		 t0, bc_spec, input.param.time_rev, 		 xml_array, "Point_Point_Wilson_Baryons");        	if (input.param.Sl_src)	  baryon(quark_propagator, phases, 		 t0, bc_spec, input.param.time_rev, 		 xml_array, "Shell_Point_Wilson_Baryons");      } // end if (Pt_snk)      // Convolute the quark propagator with the sink smearing function.      // Make a copy of the quark propagator and then overwrite it with      // the convolution.       if (input.param.Sl_snk)       {	LatticePropagator quark_prop_smr;	quark_prop_smr = quark_propagator;	sink_smear2(u, quark_prop_smr, input.param.Wvf_kind, input.param.wvf_param[loop],		    input.param.WvfIntPar[loop], input.param.j_decay);	if (input.param.Pt_src)	  baryon(quark_propagator, phases, 		 t0, bc_spec, input.param.time_rev, 		 xml_array, "Point_Shell_Wilson_Baryons");        	if (input.param.Sl_src)	  baryon(quark_propagator, phases, 		 t0, bc_spec, input.param.time_rev, 			 xml_array, "Shell_Shell_Wilson_Baryons");      } // end if (Sl_snk)    } // end if (BaryonP)    pop(xml_array);  // array element  } // end for(loop)  pop(xml_array);  // Wilson_spectroscopy  pop(xml_out);  // spectrum_w  xml_out.close();  xml_in.close();  // Time to bolt  QDP_finalize();  exit(0);}