// Spc.C  PHOBOS track finder mockup for MDC1
//   Adapted from rsFromRaw.C by A. Olszewski, K. Wozniak       DMcL 9/14/98
//   For use with phat.C    The global definitions are in Spc.h
//       See help on phat.C for details on how to convert macros
//*************************************************************************
//  Spc.C function wrapper:
void Spc(Int_t phase = 0, Int_t dflag = 0)
{
  if(ReportLevel >= 3)
     printf("Spc entered, phase = %d, dflag = %d\n",phase,dflag);
  if(phase == kIni1) {
  tman  = new TPhTemplateManager();
  tsel  = new TPhTrackSelector();
  tman->SetVerbose(1);
  tman->ReadTemplateList(gSystem->ExpandPathName(TemplFile));
  tman->SetVerbose(0);
  tsel->SetVerbose(0);
  nevent0 = 0;
  nevent1 = 0;
  nevent2 = 0;
  nevent3 = 0;
  }  // End of initialization phase wrapper 
//*************************************************************************

  //  Event loop: I've removed debug lines; these can be inserted again
  //  from rsFromRaw.C if necessary or made switchable by a flag test

  if(phase == kRun3) { 
        if( !(spectP = det->GetSubDetector("SpectP")) ) 
           printf("rsFromRaw: Positive spectrometer arm not found, shall not be processed.\n");
        if( !(spectN = det->GetSubDetector("SpectN")) )
           printf("rsFromRaw: Negative spectrometer arm not found, shall not be processed.\n");

        spectPIter = spectP->GetSensorIter();
        spectNIter = spectN->GetSensorIter();
 if (ReportLevel >= 3) printf("spect*Iters set\n");

        spectPIter->Reset();
        while( ssens = (TPhSpectSensor *) spectPIter->Next() ) {
	  if( ssens->MakeHitList("dummy","HitsFromFECn0sHitArrays") ) {
              ssens->ClearHitArray();

	      hitlist = NULL;
              if( hitlist = ssens->GetHitList() ) {
		NextHit = new TIter(hitlist);
		while( hit = (TPhHit *) NextHit() ) 
                { 
                  x = hit->GetX();
                  z = hit->GetZ();
//    ===============================================
//    selection criterion for hits for reconstruction
//    ===============================================
                  if(x>5.0 && x<80.0 && z>0.0 && z<80.0 && x<0.9*z )   
                  {
                     hitlist->Remove(hit);
  //                   printf("Remove x,z:  %11.4f  %11.4f\n", x, z);
                  }
                }
	      }
	      delete NextHit;
		  
	  }
	}
        spectNIter->Reset();
        while( ssens = (TPhSpectSensor *) spectNIter->Next() ) {
	  if( ssens->MakeHitList("dummy","HitsFromFECn0sHitArrays") ) {
              ssens->ClearHitArray();

	      hitlist = NULL;
              if( hitlist = ssens->GetHitList() ) {
		NextHit = new TIter(hitlist);
		while( hit = (TPhHit *) NextHit() ) 
		  if( hit->GetZ() > 2000.0 ) hitlist->Remove(hit);
	      }
	      delete NextHit;
		  
	  }
	}
        if (ReportLevel >= 2) printf("Start Reconstruction\n");///////////

//===============
// Reconstruction
//===============

//        CPU time limit 60 s
    templtested = tman->CreateTrackCandidates(60.0);
    ncandidates = tman->GetCandidateList()->GetEntries();

    tsel->SelectTracks(tman->GetCandidateList());
	
    tracklist = tsel->GetTrackList();
    ntrfound = tracklist->GetEntries();
    if (ReportLevel >= 3) printf("End Reconstruction\n");///////////

//===============
// Results:
//  templtested - fraction of templates tested up to CPU time limit
//  ncandidates - number of track candidates
//  ntrfound    - number of tracks after selection 
//  tracklist   - list of tracks found
// summarized:
//  nevent0     - number of events with 0 tracks found
//  nevent1     - number of events with 1 track found
//   ..
//  nevent3     - number of events with 3 or more tracks found
//===============
        if(ntrfound == 0)  nevent0++;
        if(ntrfound == 1)  nevent1++;
        if(ntrfound == 2)  nevent2++;
        if(ntrfound >= 3)  nevent3++;

	printf("Found %d candidates   and %d tracks in the fraction of templates %f \n", ncandidates, ntrfound, templtested);
	for(i=0;i<tracklist->GetEntries();i++)
	{
	    track = (TPhTrack *) tracklist->At(i);
	    //	    mom->Fill(track->GetP());
	    printf("Found track %g (%g , %g , %d)\n",
              track->GetP(),track->GetPx(),track->GetPz(),track->GetCharge());
	}
        tman->GetCandidateList()->Delete();
  }  // End of run phase wrapper 
//*************************************************************************

  if(phase == kDisf) {
//==========================================
// display histogram and make final cleaning
//==========================================
    //    delete det;   // This causes seg viol in interactive mode..???
    printf(" Events with  0  tracks found: %8d\n", nevent0);
    printf(" Events with  1  tracks found: %8d\n", nevent1);
    printf(" Events with  2  tracks found: %8d\n", nevent2);
    printf(" Events with  >2 tracks found: %8d\n", nevent3);

    //  Put these lines in later, after defining canvas explicitly!
    //  TH1F mom("mom","tot. momentum of reconstructed ptls",200,0.350,0.550);
    //  mom.DrawCopy();
    //   c2->Modified();  // Need to name the canvas since there are others
    //   c2->cd();
  }  // end of final display wrapper
//*************************************************************************
}