// Value of CheckMC tells what type of events this is:
// 0 = on peak data, 1 = off peak, 2 = Jpsi K+ MC, 3 = Other B-->Jpsi MC
// 4 = Other BB MC, 5 = continuum MC

// C. Hearty, 21-Jun-2006. This version for rootTuples produced 
// with BtaTupleMaker

Int_t histAnalysis::CheckMC(Int_t mcLund[200], Int_t mothIdx[200], Int_t dauLen[200], Int_t dauIdx[200], Int_t mcLen, Float_t sqrts) {

  // The numbers are the absolute values of the lund IDs; see PDT/pdt.table
  Int_t idb0(511),idbch(521), idjpsi(443), idk(321), idgam(22), 
    idmu(13), ide(11), idtau(15), idd(1), idu(2), ids(3), idc(4);

  Int_t b0flag(0), bchflag(0), jpsiflag(0), signalflag(0);
    
  // Return value
  Int_t etm = -1;

  //..Real data
  if (mcLen == 0) {
    if( sqrts > 10.56 ) {
      etm = 0;
    } else {
      etm = 1;
    }
    return etm;
  }

  //..Some kind of MC. The second particle in the list will tell us a lot
  Int_t lundID1 = abs(mcLund[1]);

  //..Continuum
  if (lundID1 == idc || lundID1 == idu || lundID1 == ids || lundID1 == idd
      || lundID1 == idtau) {
    etm = 5;
    return etm;
  }

  //..Ought to be a B event. If it isn't, just give up.
  if (lundID1 == idb0) b0flag = 1;
  if (lundID1 == idbch) bchflag = 1;
  if( b0flag == 0 && bchflag == 0 ) return etm;

  //..An inclusive Jpsi event is a B event that includes Jpsi --> ee or mumu
  Int_t Jpsitoll[200];
  for (Int_t k=0; k<mcLen; k++) {
    if( mcLund[k] == idjpsi ) {

      //..This is a Jpsi. Check if any of the daughters are leptons.
      Jpsitoll[k] = 0;
      for (Int_t kd=dauIdx[k]; kd<dauIdx[k]+dauLen[k]; kd++) {
	if( mcLund[kd] == ide || mcLund[kd] == idmu ) { 
	  jpsiflag = 1;
	  Jpsitoll[k] = 1;
	}
      }

      //..If this is an inclusive Jpsi event, check that the first daughter is
      //  a Jpsi and the 2nd a K+. If there are any other daughters, they
      //  are required to be photons.
      if( Jpsitoll[k]==1 && abs(mcLund[mothIdx[k]])==idbch ) {
	Int_t ibmoth = mothIdx[k];
	Int_t numK(0), numGam(0);
	for (Int_t kd=dauIdx[ibmoth]; kd<dauIdx[ibmoth]+dauLen[ibmoth]; kd++) {
	  if( abs(mcLund[kd])==idk ) numK++;
	  if( mcLund[kd]==idgam ) numGam++;
	}
	if( numK==1 && numGam==(dauLen[ibmoth]-2) ) signalflag=1;
      }
    }
  }

  if (bchflag) etm = 4;
  if (b0flag) etm = 4;
  if (jpsiflag) etm = 3;
  if (signalflag) etm = 2;

  // Dump the MC table if requested
  Int_t dumpTable(0);
  if( dumpTable ) {
    cout << "       k" << "  mcLund" << " mothIdx" << "  dauLen" << "  dauIdx" << endl;
    for (Int_t k=0; k<mcLen; k++) {
      cout << setw(8) << k  << setw(8) << mcLund[k]  << setw(8) << mothIdx[k]  << setw(8) << dauLen[k]  << setw(8) << dauIdx[k]  << endl;
    }
  }

  return etm;
}