{
/*
This macro sort wafers onto ladders

Jerome Baudot - baudot@in2p3.fr - IReS Strasbourg - 2001/07/07
*/

  gROOT->Reset();

  // Parameters of the sorting

  const Int_t Ntest = 500; // number of wafers to test (DEBUGGING)

  const Int_t Nladders = 23; // number of ladders to build
  const Int_t NwafersPerLadder = 17; // number of wafers to make a ladder
  const Int_t maxQuality = 3; // badest quality acceptable

  const Int_t sortMethod =1; // 0=random choice, 
                             // 1=use criteria,
                             // 2=maximize operating domain with minimal Vbias
                             // 3= ??

  const Float_t weightDomain = 15.;
  const Float_t weightVbias = 1.;
  const Float_t weightItot = 1.;
  const Float_t weightDead = 1.;

  const Int_t domainMin =10; // minimal width of operating domain in volts

  // Parameters for the outputs display

  const Int_t printDebug0 = 0; // print data from mySQl
  const Int_t printDebug1 = 0; // print steps of sorting algortihms
  const Int_t printSorted = 0; // print list of sorted wafers
  const Int_t printLadder = 1; // print contents of ladder
  const Int_t printHisto = 1;  // print histos into .gif files


/* A priori, you did not need to change a single line below this point . */

  // Some variables

  Int_t entries, entry=0, i, j, k, choice, Iladder, leftWafers;
  Int_t tries, maxTries=1000;
  Int_t newDomain[2];
  Char_t fileName[50];
  TMath math;
  TRandom rand;
  Char_t name[50],title[50];


  // Containers for wafers and ladders

  struct wafer_t {
    //  char name[9];
    Int_t vbias;
    Int_t vmax;
    Float_t itotal;
    Int_t deadp;
    Int_t deadn;
    Int_t position;
  };

  char  waferName[500][9];
  Int_t waferQuality[500];
  Int_t waferVbias[500];
  Int_t waferVmax[500];
  Float_t waferItot[500];
  Int_t waferDeadp[500];
  Int_t waferDeadn[500];
  Int_t waferAssembled[500];
  Float_t waferCriteria[500];

  Int_t sortedVbias[500];
  Int_t sortedDomain[500];
  Int_t sortedItot[500];
  Int_t sortedDead[500];
  Int_t sorted[500];

  Int_t indexVbias[500];
  Int_t indexDomain[500];
  Int_t indexItot[500];
  Int_t indexDead[500];

  Int_t ladderWafers[Nladders][NwafersPerLadder];
  Int_t ladderDomain[Nladders][2];
  Int_t ladderItot[Nladders];
  Int_t ladderDeadp[Nladders];
  Int_t ladderDeadn[Nladders];
  Int_t ladderIncomplete[Nladders];


  // Create files

  sprintf( name, "waferSort%1d.root", sortMethod );
  TFile *histFile = new TFile( name,"RECREATE");


  // Create histograms

  sprintf( name, "ladDomain%1d", sortMethod );
  sprintf( title, "Ladders operation domain - methode %1d", sortMethod );
  TH2S *ladDomain = new TH2S( name, title, Nladders, 0, Nladders, 40, 0, 40);
  ladDomain->SetXTitle( "ladder #" );
  ladDomain->SetYTitle( "(V)" );
 
  sprintf( name, "ladItot%1d", sortMethod );
  sprintf( title, "Ladders total leakage current - methode %1d", sortMethod );
  TH2F *ladItot = new TH2F( name, title, Nladders, 0, Nladders, 40, 0., 40.);
  ladItot->SetXTitle( "ladder #" );
  ladItot->SetYTitle( "Current (uA)" );

  sprintf( name, "ladDeadp%1d", sortMethod );
  sprintf( title, "Ladders dead strips P - methode %1d", sortMethod );
  TH2S *ladDeadp = new TH2S(name, title, Nladders, 0, Nladders, 100, 0, 100);
  ladDeadp->SetXTitle( "ladder #" );
  ladDeadp->SetYTitle( "number of dead strips" );

  sprintf( name, "ladDeadn%1d", sortMethod );
  sprintf( title, "Ladders dead strips N - methode %1d", sortMethod );
  TH2S *ladDeadn = new TH2S( name, title, Nladders, 0, Nladders, 100, 0, 100);
  ladDeadn->SetXTitle( "ladder #" );
  ladDeadn->SetYTitle( "number of dead strips" );

  sprintf( name, "ladVdVm%1d", sortMethod );
  sprintf( title, "Ladders operation domain - methode %1d", sortMethod );
  TH2S *ladVdVm = new TH2S( name, title, Nladders, 0, Nladders, 80, 10, 90);
  ladVdVm->SetXTitle( "ladder #" );
  ladVdVm->SetYTitle( "Vbias - Vmax (V)" );

  TH2S *wafVdVm[Nladders];
  for( i=0; i<Nladders; i++ ) {
    sprintf( name, "waferVdVm%2d", i);
    sprintf( title, "Ladder %2d, operation domain", i);
    wafVdVm[i] = new TH2S(name, title, NwafersPerLadder, 0, NwafersPerLadder, 80, 10, 90);
    wafVdVm[i]->SetXTitle( "wafer #" );
    wafVdVm[i]->SetYTitle( "Vbias - Vmax (V)" );
  }

  // histos for unassociated wafers
  sprintf( name, "leftVbias%1d", sortMethod );
  sprintf( title, "Vbias for unassociated wafers - methode %1d", sortMethod );
  TH1S *leftVbias = new TH1S(name, title, 80, 10, 90);
  leftVbias->SetXTitle( "depletion voltage (V)" );
  sprintf( name, "leftVmax%1d", sortMethod );
  sprintf( title, "Vmax for unassociated wafers - methode %1d", sortMethod );
  TH1S *leftVmax = new TH1S(name, title, 80, 10, 90);
  leftVmax->SetXTitle( "breakdown voltage (V)" );
  sprintf( name, "leftDomain%1d", sortMethod );
  sprintf( title, "Domain for unassociated wafers - methode %1d", sortMethod );
  TH1S *leftDomain = new TH1S(name, title, 80, 10, 90);
  leftDomain->SetXTitle( "Vmax-Vbias (V)" );
  sprintf( name, "leftItot%1d", sortMethod );
  sprintf( title, "Total leakage current for unassociated wafers - methode %1d", sortMethod );
  TH1F *leftItot = new TH1F(name, title, 50, 0., 5.);
  leftItot->SetXTitle( "current (uA)" );
  sprintf( name, "leftDeadp%1d", sortMethod );
  sprintf( title, "Dead strips P for unassociated wafers - methode %1d", sortMethod );
  TH1S *leftDeadp = new TH1S(name, title, 100, 0, 100);
  leftDeadp->SetXTitle( "number of dead strips" );
  sprintf( name, "leftDeadn%1d", sortMethod );
  sprintf( title, "Dead strips N for unassociated wafers - methode %1d", sortMethod );
  TH1S *leftDeadn = new TH1S(name, title,100 , 0, 100);
  leftDeadn->SetXTitle( "number of dead strips" );
  sprintf( name, "leftDomainVbias%1d", sortMethod );
  sprintf( title, "Domain vs Vbias for unassociated wafers - methode %1d", sortMethod );
  TH2S *leftDomainVbias = new TH2S(name, title, 80, 10, 90, 80, 10, 90);
  leftDomainVbias->SetXTitle( "depletion voltage (V)" );
  leftDomainVbias->SetYTitle( "Vmax-Vbias (V) (V)" );


  // Init connection to the mySQL database

  const char* host = "mysql://wwwstar-sbg.in2p3.fr/ssd";
  //const char* host = "mysql://wulkan.if.pw.edu.pl/ssd";
  //const char* host = "mysql://venus.if.pw.edu.pl/ssd";
  const char* user = "jeromeb";
  const char* pass = "jeromeb";

  TSQLServer *db = TSQLServer::Connect(host,user,pass);
    
  printf("Successfull conection to mySQL (%s) host/db (%s) \n",db->ServerInfo(), host);


  // Get the info for all wafer from wafer_test table

  TSQLRow *row, *row2;
  TSQLResult *result, *result2;

  char query2[500];    

  const char *line2 = "SELECT obj_quality, obj_assembled FROM object_header WHERE obj_name='%s'";

  const char *query ="SELECT obj_name, wafer_depletion_v, wafer_breakdown_v, wafer_it_at_stab, wafer_p_dead, wafer_n_dead FROM wafer_test";

  result = db->Query(query);

  while( ( row=result->Next() ) && ( entry<Ntest ) ) {

    // Query also quality from object_header
    sprintf(query2, line2, row->GetField(0));
    result2 = db->Query(query2);
    row2=result2->Next();

    // Select only wafers with quality<= minQuality
    if( atoi( row2->GetField(0) )<=maxQuality ) {

    if( printDebug0 ) printf("name=%s quality=%d assembled=%s vbias=%d vmax=%d itot=%6.3f deadp=%d deadn=%d entry=%d\n", row->GetField(0), atoi( row2->GetField(0) ), row2->GetField(1), atoi( row->GetField(1) ), atoi( row->GetField(2) ), atof( row->GetField(3) ), atoi( row->GetField(4) ), atoi( row->GetField(5) ), entry);
    
//     aWafer.name = row->GetField(0);
//     aWafer.vbias = atoi( row->GetField(1) );
//     aWafer.vmax = atoi( row->GetField(2) );
//     aWafer.itot = atof( row->GetField(3) );
//     aWafer.deadp = atoi( row->GetField(4) );
//     aWafer.deadn = atoi( row->GetField(5) );
//     printf("name=%s vbias=%d vmax=%d itot=%6.3f deadp=%d deadn=%d entry=%d\n", aWafer.name, aWafer.vbias, aWafer.vmax , aWafer.itot , aWafer.deadp , aWafer.deadn , entry++);
      sprintf( waferName[entry], row->GetField(0));
      waferQuality[entry] = atoi( row2->GetField(0) );
//    waferAssembled[entry] = row2->GetField(1)=='yes'?1:0;
      waferAssembled[entry] = 0;
      waferVbias[entry] = atoi( row->GetField(1) );
      waferVmax[entry] = atoi( row->GetField(2) );
      waferItot[entry] = atof( row->GetField(3) );
      waferDeadp[entry] = atoi( row->GetField(4) );
      waferDeadn[entry] = atoi( row->GetField(5) );
      if( printDebug0) printf("name=%s quality=%d assembled=%s vbias=%d vmax=%d itot=%6.3f deadp=%d deadn=%d entry=%d\n", waferName[entry], waferQuality[entry], waferAssembled[entry], waferVbias[entry], waferVmax[entry], waferItot[entry], waferDeadp[entry], waferDeadn[entry], entry);
      sorted[entry] = sortedVbias[entry] = sortedDomain[entry] = sortedItot[entry] = sortedDead[entry] = entry++;

    } // end select on quality

      delete row;
      delete row2;
      delete result2;

  } // end loop on db entries

  entries = leftWafers = entry; // number of available wafers
  printf("\n** %d wafers read from the database **\n", entries);

  delete result;


  // Sort wafers with respect to their Vbias (the lower, the better)

  for( i=0; i<entries; i++ ) {
    for( j=i+1; j<entries; j++ ) {
      if( waferVbias[sortedVbias[i]] > waferVbias[sortedVbias[j]] ) {
	k = sortedVbias[i];
	sortedVbias[i] = sortedVbias[j];
	sortedVbias[j] = k;
      }
    }
  }
  if( printSorted ) printf("\n**** VBias-sorted list of wafers:\n");
  for( entry=0; entry<entries; entry++ ) {
    indexVbias[sortedVbias[entry]] = entry;
    if( printSorted ) printf(" entry %3d: %s Vbias=%2d V\n", entry, waferName[sortedVbias[entry]], waferVbias[sortedVbias[entry]]);
  }


  // Sort wafers with respect to their operation domain Vmax-Vbias

  for( i=0; i<entries; i++ ) {
    for( j=i+1; j<entries; j++ ) {
      if( (waferVmax[sortedDomain[i]]-waferVbias[sortedDomain[i]]) < (waferVmax[sortedDomain[j]]-waferVbias[sortedDomain[j]]) ) {
	k = sortedDomain[i];
	sortedDomain[i] = sortedDomain[j];
	sortedDomain[j] = k;
      }
    }
  }
  if( printSorted ) printf("\n**** VDomain-sorted list of wafers:\n");
  for( entry=0; entry<entries; entry++ ) {
    indexDomain[sortedDomain[entry]] = entry;
    if( printSorted ) printf(" entry %3d: %s domain=%2d V\n", entry, waferName[sortedDomain[entry]], waferVmax[sortedDomain[entry]]-waferVbias[sortedDomain[entry]]);
  }


  // Sort wafers with respect to their total leackage current Itot

  for( i=0; i<entries; i++ ) {
    for( j=i+1; j<entries; j++ ) {
      if( waferItot[sortedItot[i]] > waferItot[sortedItot[j]] ) {
	k = sortedItot[i];
	sortedItot[i] = sortedItot[j];
	sortedItot[j] = k;
      }
    }
  }
  if( printSorted ) printf("\n**** Itot-sorted list of wafers:\n");
  for( entry=0; entry<entries; entry++ ) {
    indexItot[sortedItot[entry]] = entry;
    if( printSorted ) printf(" entry %3d: %s Itot=%6.3f uA\n", entry, waferName[sortedItot[entry]], waferItot[sortedItot[entry]]);
  }


  // Sort wafers with respect to their total number of dead strips

  for( i=0; i<entries; i++ ) {
    for( j=i+1; j<entries; j++ ) {
      if( (waferDeadp[sortedDead[i]]+waferDeadn[sortedDead[i]]) > (waferDeadp[sortedDead[j]]+waferDeadn[sortedDead[j]]) ) {
	k = sortedDead[i];
	sortedDead[i] = sortedDead[j];
	sortedDead[j] = k;
      }
    }
  }
  if( printSorted ) printf("\n**** Dead-sorted list of wafers:\n");
  for( entry=0; entry<entries; entry++ ) {
    indexDead[sortedDead[entry]] = entry;
    if( printSorted ) printf(" entry %3d: %s dead=%4d\n", entry, waferName[sortedDead[entry]], waferDeadp[sortedDead[entry]]+waferDeadn[sortedDead[entry]]);
  }


  // Sort wafers with respect to a criteria
  // a*Domain + b*Itot + c*Dead where a,b,c are respective weights

  // compute criteria first
  for( i=0; i<entries; i++ ) {
    waferCriteria[i] = weightVbias*(Float_t)indexVbias[i] + weightDomain*(Float_t)indexDomain[i] + weightItot*(Float_t)indexItot[i] + weightDead*(Float_t)indexDead[i];
  }

  for( i=0; i<entries; i++ ) {
    for( j=i+1; j<entries; j++ ) {
      if( waferCriteria[sorted[i]] > waferCriteria[sorted[j]] ) {
	k = sorted[i];
	sorted[i] = sorted[j];
	sorted[j] = k;
      }
    }
  }
  if( printSorted ) printf("\n**** Criteria-sorted list of wafers:\n");
  for( entry=0; entry<entries; entry++ ) {
    if( printSorted ) printf(" %3d: %s criteria=%7.2f, Vbias=%2dV Domain=%2dV Itot=%6.3fuA Dead=%4d\n", entry, waferName[sorted[entry]], waferCriteria[sorted[entry]], waferVbias[sorted[entry]], waferVmax[sorted[entry]]-waferVbias[sorted[entry]], waferItot[sorted[entry]], waferDeadp[sorted[entry]]+waferDeadn[sorted[entry]]);
  }


  // Start grouping wafers on ladders
  printf("\n** Trying to associate wafers onto %2d ladders of %2d wafers...\n", Nladders, NwafersPerLadder);

  // Initialize flag
  for( Iladder =0; Iladder<Nladders; Iladder++ ) {
    ladderIncomplete[Iladder] = 0;
  }



//********************************************************
  //*** ZERO ALGORITHM ***
  if( sortMethod == 0 ) {

    printf("\n *** Starting zero algorithm ***\n\n");
    printf("Random choice assuring a minimal operation domain.\n");

    // loop on ladders
    for( Iladder=0; Iladder<Nladders; Iladder++ ) {
      if( printDebug1 ) printf("\nAssembling ladder %2d\n\n", Iladder);
      ladderDomain[Iladder][0] = 0;
      ladderDomain[Iladder][1] = 100;

      for( i=0; i<NwafersPerLadder; i++ ) {
	if( printDebug1 ) printf(" Looking for wafer %d (%d left)\n", i, leftWafers);

	tries = 0;
	j = (int)rand->Rndm()*entries;
	while( ( tries++ < maxTries ) && ( waferAssembled[j] || waferVbias[j]>ladderDomain[Iladder][1] || ( math.Min( ladderDomain[Iladder][1], waferVmax[j]) - math.Max( ladderDomain[Iladder][0], waferVbias[j]) < domainMin ) ) ) { 
	  j = (int)(rand->Rndm()*entries);
	  if( printDebug1 ) printf( "   trying %3d...\n", j);
	}

	if( tries < maxTries ) {
	  ladderWafers[Iladder][i] = j;
	  ladderDomain[Iladder][0] = math.Max( ladderDomain[Iladder][0], waferVbias[j] );
	  ladderDomain[Iladder][1] = math.Min( ladderDomain[Iladder][1], waferVmax[j] );
	  ladderItot[Iladder] += waferItot[j];
	  ladderDeadp[Iladder] += waferDeadp[j];
	  ladderDeadn[Iladder] += waferDeadn[j];
	  waferAssembled[j] = 1;
	  leftWafers--;
	  if( printDebug1 ) printf("  Found wafer %2d = (j=%3d) %s, domain %2d-%2d\n", i, j, waferName[j], waferVbias[j], waferVmax[j]);
	  
	} // end test if one wafer has been found
	// else a wafer is missing
	else {
	  printf("MISSING WAFERS to make wafer %2d of ladder %2d (%d wafers left)!!!\n", i, Iladder, leftWafers);
	  ladderIncomplete[Iladder] = i;
	}

      } // end  loop on wafers    

    } // end loop on ladders

  }  // end of method zero



//********************************************************
  //*** FIRST ALGORITHM ***
  else if( sortMethod == 1 ) {
    printf("\n *** Starting first algorithm ***\n\n");
    printf("Parameters are:\n    maximum quality accepted=%1d\n    weights for domain=%4.1f, Vbias=%4.1f, Itot=%4.1f, Dead=%2d\n    minimal operating domain = %2d V\n\n", maxQuality, weightDomain, weightVbias, weightItot, weightDead, domainMin);

    // loop on the number of wafers per ladder
    for( i=0; i<NwafersPerLadder; i++ ) {
      if( printDebug1 ) printf("Looking for wafer %d (%d left)\n", i, leftWafers);
      // loop on ladders
      for( Iladder =0; Iladder<Nladders; Iladder++ ) {
	// select the next best wafer wrt criteria
	// if first wafer on ladder
	if( i==0 ) {
	  j = 0;
	  while( j<entries && waferAssembled[sorted[j]] ) { j++; }
	  ladderWafers[Iladder][i] = sorted[j];
	  ladderDomain[Iladder][0] = waferVbias[sorted[j]];
	  ladderDomain[Iladder][1] = waferVmax[sorted[j]];
	  ladderItot[Iladder] = waferItot[sorted[j]];
	  ladderDeadp[Iladder] = waferDeadp[sorted[j]];
	  ladderDeadn[Iladder] = waferDeadn[sorted[j]];
	  waferAssembled[sorted[j]] = 1;
	  leftWafers--;
	  if( printDebug1 ) printf("  Found for ladder %2d, first wafer = (j=%3d) %s\n", Iladder, j, waferName[sorted[j]]);
	}
	// else check it matches the operating range
	else {
	  j=0;
	  while( j<entries && ( waferAssembled[sorted[j]] || ( math.Min( waferVmax[sorted[j]], ladderDomain[Iladder][1]) - math.Max(waferVbias[sorted[j]], ladderDomain[Iladder][0]) )<domainMin ) ) { 
 	    if( printDebug0 ) {
 	      printf("ladder domain %2d-%2d, (j=%3d,ass=%1d) Vbias=%2d Vmax=%2d, max=%2d, min=%2d, diff=%2d (/%2d)\n", ladderDomain[Iladder][0], ladderDomain[Iladder][1], j, waferAssembled[sorted[j]], waferVbias[sorted[j]], waferVmax[sorted[j]], math.Max(waferVbias[sorted[j]], ladderDomain[Iladder][0]), math.Min( waferVmax[sorted[j]], ladderDomain[Iladder][1]), math.Min( waferVmax[sorted[j]], ladderDomain[Iladder][1]) - math.Max(waferVbias[sorted[j]], ladderDomain[Iladder][0]), domainMin); 
 	    }
	    j++; 
	  }
	  if( j<entries ) {
	    ladderWafers[Iladder][i] = sorted[j];
	    if( ladderDomain[Iladder][0] < waferVbias[sorted[j]] ) {
	      ladderDomain[Iladder][0] = waferVbias[sorted[j]];
	    }
	    if( ladderDomain[Iladder][1] > waferVmax[sorted[j]] ) {
	      ladderDomain[Iladder][1] = waferVmax[sorted[j]];
	    }
	    ladderItot[Iladder] += waferItot[sorted[j]];
	    ladderDeadp[Iladder] += waferDeadp[sorted[j]];
	    ladderDeadn[Iladder] += waferDeadn[sorted[j]];
	    waferAssembled[sorted[j]] = 1;
	    leftWafers--;
	    if( printDebug1 ) printf("  Found ladder %2d, wafer %2d = (j=%3d) %s\n", Iladder, i, j, waferName[sorted[j]]);
	  }
	  else {
	    printf("MISSING WAFERS to make wafer %2d of ladder %2d (%d wafers left)!!!\n", i, Iladder, leftWafers);
	    ladderIncomplete[Iladder] = i;
	  }
	}
      } // end loop on ladders
    } // end loop on the number of wafers per ladder
    

  } // end of first method



  

//********************************************************
  //*** SECOND METHOD ***
  else if( sortMethod == 2 ) {

    printf("\n *** Starting second algorithm ***\n\n");
    printf("Try to maximize operating domain for lowest Vbias: minimal operating domain = %2d V\n\n\n", domainMin);

    // loop on ladders
    for( Iladder=0; Iladder<Nladders; Iladder++ ) {
      if( printDebug1 ) printf("\nAssembling ladder %2d\n\n", Iladder);

      // Associate first ladder as the one with smallest Vbias available
      if( printDebug1 ) printf(" Looking for first wafer (%d left)\n", leftWafers);
      j = 0;
      while( j<entries && ( waferAssembled[sortedVbias[j]] || ( waferVmax[sortedVbias[j]]-waferVbias[sortedVbias[j]] ) < domainMin  ) ) { j++; }
      ladderWafers[Iladder][0]  = sortedVbias[j];
      ladderDomain[Iladder][0] = waferVbias[sortedVbias[j]];
      ladderDomain[Iladder][1] = waferVmax[sortedVbias[j]];
      ladderItot[Iladder]      = waferItot[sortedVbias[j]];
      ladderDeadp[Iladder]     = waferDeadp[sortedVbias[j]];
      ladderDeadn[Iladder]     = waferDeadn[sortedVbias[j]];
      waferAssembled[sortedVbias[j]] = 1;
      leftWafers--;
      if( printDebug1 ) printf("  Found first wafer = (j=%3d) %s, domain %2d-%2d\n", j, waferName[sortedVbias[j]], waferVbias[sortedVbias[j]], waferVmax[sortedVbias[j]]);

      // loop on following wafers to make the ladder
      for( i=1; i<NwafersPerLadder; i++ ) {
	if( printDebug1 ) printf(" Looking for wafer %d (%d left)\n", i, leftWafers);

	// look for the wafer which maximize the operating domain
	// starting from the one with the lowest Vbias
	j = 0;
	while( j<entries && ( waferAssembled[sortedVbias[j]] || waferVbias[sortedVbias[j]]>ladderDomain[Iladder][1] || ( math.Min( ladderDomain[Iladder][1], waferVmax[sortedVbias[j]]) - math.Max( ladderDomain[Iladder][0], waferVbias[sortedVbias[j]]) < domainMin ) ) ) { j++; }

	// Test if one wafer has been found
	if( j<entries ) {
	  choice = j;
	  newDomain[0] = math.Max( ladderDomain[Iladder][0], waferVbias[sortedVbias[choice]]);
	  newDomain[1] = math.Min( ladderDomain[Iladder][1], waferVmax[sortedVbias[choice]]);
	  if( printDebug0 ) printf(" -- found first choice=%2d %s, Vbias=%2d, Vmax=%2d\n", choice, waferName[sortedVbias[choice]], waferVbias[sortedVbias[choice]], waferVmax[sortedVbias[choice]]);
	  // try to find a wafer which provide a larger operating domain
	  // in the rest of the array
	  for( k=j+1; k<entries; k++ ) {
	    if( printDebug0 ) printf(" - testing %s, assembled?=%1d, Vbias=%2d, Vmax=%2d, new domain %2d-%2d, current domain %2d-%2d\n", waferName[sortedVbias[k]], waferAssembled[sortedVbias[k]], waferVbias[sortedVbias[k]], waferVmax[sortedVbias[k]], newDomain[0], newDomain[1], ladderDomain[Iladder][0], ladderDomain[Iladder][1]);
	    if(    waferAssembled[sortedVbias[k]]==0
                && waferVbias[sortedVbias[k]] <= ladderDomain[Iladder][1]
	        && ( math.Min( ladderDomain[Iladder][1], waferVmax[sortedVbias[k]]) - math.Max( ladderDomain[Iladder][0], waferVbias[sortedVbias[k]]) ) > ( newDomain[1] - newDomain[0] ) 
                && ( math.Min( ladderDomain[Iladder][1], waferVmax[sortedVbias[k]]) - math.Max( ladderDomain[Iladder][0], waferVbias[sortedVbias[k]]) ) >= domainMin
              ) {
	      choice = k;
	      if( printDebug0 ) printf(" -- found next choice=%2d %s, Vbias=%2d, Vmax=%2d\n", choice, waferName[sortedVbias[choice]], waferVbias[sortedVbias[k]], waferVmax[sortedVbias[k]]);
	      newDomain[0] = math.Max( ladderDomain[Iladder][0], waferVbias[sortedVbias[choice]]);
	      newDomain[1] = math.Min( ladderDomain[Iladder][1], waferVmax[sortedVbias[choice]]);
	    }
	  } // end loop on the rest of the array
	
	  ladderWafers[Iladder][i] = sortedVbias[choice];
	  if( ladderDomain[Iladder][0] < waferVbias[sortedVbias[choice]] ) {
	    ladderDomain[Iladder][0] = waferVbias[sortedVbias[choice]];
	  }
	  if( ladderDomain[Iladder][1] > waferVmax[sortedVbias[choice]] ) {
	    ladderDomain[Iladder][1] = waferVmax[sortedVbias[choice]];
	  }
	  ladderItot[Iladder] += waferItot[sortedVbias[choice]];
	  ladderDeadp[Iladder] += waferDeadp[sortedVbias[choice]];
	  ladderDeadn[Iladder] += waferDeadn[sortedVbias[choice]];
	  waferAssembled[sortedVbias[choice]] = 1;
	  leftWafers--;
	  if( printDebug1 ) printf("  Found wafer %2d = (j=%3d) %s, domain %2d-%2d\n", i, choice, waferName[sortedVbias[choice]], waferVbias[sortedVbias[choice]], waferVmax[sortedVbias[choice]]);

	} // end test if one wafer has been found
	// else a wafer is missing
	else {
	  printf("MISSING WAFERS to make wafer %2d of ladder %2d (%d wafers left)!!!\n", i, Iladder, leftWafers);
	  ladderIncomplete[Iladder] = i;
	}

      } // end loop on following wafers to make the ladder
    } // end loop on ladders
    
  } // end of second method



//********************************************************
  //*** THIRD ALGORITHM ***
  else if( sortMethod == 3 ) {

    printf("\n *** Starting third algorithm ***\n\n");
    printf("Not yet implemented, giving up!\n");

    leftWafers = entries; // re-initialize number of available wafers

    // set first wafers of each ladder

    for( Iladder =0; Iladder<Nladders; Iladder++ ) {
      // select the next best wafer wrt the domain
      if( printDebug1 ) printf("Looking for first wafer\n");
      j = 0;
      while( j<entries && waferAssembled[sortedDomain[j]] ) { j++; }
      ladderWafers[Iladder][0] = sortedDomain[j];
      ladderDomain[Iladder][0] = waferVbias[sortedDomain[j]];
      ladderDomain[Iladder][1] = waferVmax[sortedDomain[j]];
      ladderItot[Iladder] = waferItot[sortedDomain[j]];
      ladderDeadp[Iladder] = waferDeadp[sortedDomain[j]];
      ladderDeadn[Iladder] = waferDeadn[sortedDomain[j]];
      waferAssembled[sortedDomain[j]] = 1;
      leftWafers--;
      if( printDebug1 ) printf("  Found ladder %d, first wafer = (j=%d) %s\n", Iladder, j, waferName[sortedDomain[j]]);
    }

    // attribute successive wafers in order to maximize the operating domain
  
  } // end third method



  //****************************************************
  // Display the SSD constituents and characteristics

  if( printLadder ) {
    printf("\n\n******* THE SSD LOOKS LIKE:\n\n");
    for( i=0; i<Nladders; i++ ) {
      printf(" LADDER %2d: ", i);
      for( j=0; j<(ladderIncomplete[i]==0?NwafersPerLadder:ladderIncomplete[i]); j++ ) {
	printf("%s ", waferName[ladderWafers[i][j]]);
      }
      if( ladderIncomplete[i] ) {
	printf(" *INCOMPLETE* %2d wafers out of %2d\n\n", ladderIncomplete[i], NwafersPerLadder );
      }
      else {
	printf("\n\n");
      }
    }
    printf("\n");
  }
    
    printf("  LADDERS CHARACTERISTICS\n");
    for( i=0; i<Nladders; i++ ) {
      if( ladderIncomplete[i] ) {
	printf(" LADDER %2d: %2d-%2d V, %7.3f uA, P-side=%5d N-side=%5d **INCOMPLETE %2d<**\n", i, ladderDomain[i][0], ladderDomain[i][1], ladderItot[i], ladderDeadp[i], ladderDeadn[i], ladderIncomplete[i]);
      }
      else {
	printf(" LADDER %2d: %2d-%2d V, %7.3f uA, P-side=%5d N-side=%5d\n", i, ladderDomain[i][0], ladderDomain[i][1], ladderItot[i], ladderDeadp[i], ladderDeadn[i]);
      }
      // Fill histograms
      ladDomain->Fill(i, ladderDomain[i][1]-ladderDomain[i][0]);
      ladItot->Fill(i, ladderItot[i]);
      ladDeadp->Fill(i, ladderDeadp[i]);
      ladDeadn->Fill(i, ladderDeadn[i]);
      for( k=ladderDomain[i][0]; k<=ladderDomain[i][1];k++ ) {
	ladVdVm->Fill(i, k);
      } 
      for( j=0; j<(ladderIncomplete[i]==0?NwafersPerLadder:ladderIncomplete[i]); j++ ) {
	for( k=waferVbias[ladderWafers[i][j]]; k<=waferVmax[ladderWafers[i][j]]; k++ ) {
	  wafVdVm[i]->Fill(j, k);
	}
      }
    }
    printf("\n");


  // Fill histos for unassociated wafers
  for( i=0; i<entries; i++ ) {
    if( !waferAssembled[i] ) {
      leftVbias->Fill(waferVbias[i]);
      leftVmax->Fill(waferVmax[i]);
      leftItot->Fill(waferItot[i]);
      leftDomain->Fill(waferVmax[i]-waferVbias[i]);
      leftDeadp->Fill(waferDeadp[i]);
      leftDeadn->Fill(waferDeadn[i]);
      leftDomainVbias->Fill(waferVbias[i], waferVmax[i]-waferVbias[i]);
    }
  }

  printf("** %3d wafers were not associated (out of %3d) **\n", leftWafers, entries);
  

  // Write histos to file

  histFile->Write();
//   histFile->Close();



  //****************************************************
  // Draw histos

  if( printHisto ) {

    // Set graphics options

    gROOT->SetStyle("Plain");
    gStyle->SetOptTitle(1);
    gStyle->SetOptStat(0);
    
    sprintf( name, "methode%1d", sortMethod );
    sprintf( title, "* Methode %1d *", sortMethod );
    TCanvas*  cA = new TCanvas(name , title, 19,66,721,547);
    
    ladVdVm->Draw("col");
    cA->Modified();
    sprintf( name, "ladVdVm%1d.gif", sortMethod );
    cA->Print( name );

    ladDomain->Draw("col");
    cA->Modified();
    sprintf( name, "ladDomain%1d.gif", sortMethod );
    cA->Print( name );

    ladItot->Draw("col");
    cA->Modified();
    sprintf( name, "ladItot%1d.gif", sortMethod );
    cA->Print( name );

    ladDeadp->Draw("col");
    cA->Modified();
    sprintf( name, "ladDeadp%1d.gif", sortMethod );
    cA->Print( name );

    ladDeadn->Draw("col");
    cA->Modified();
    sprintf( name, "ladDeadn%1d.gif", sortMethod );
    cA->Print( name );

//     cA->Close();
//     cA = new TCanvas(name , title, 19,66,721,923);
//     j = (int)Nladders/3;
    for( i=0; i<Nladders; i++ ) {
      wafVdVm[i]->Draw("col");
      cA->Modified();
      if( i<=9 ) {
	sprintf( name, "wafVdVm0%1d_%1d.gif", i, sortMethod );
      }
      else {
	sprintf( name, "wafVdVm%2d_%1d.gif", i, sortMethod );
      }
      cA->Print( name );
    }

    gStyle->SetLineColor(2);
    gStyle->SetLineWidth(2);
    gStyle->SetFillColor(2);

    leftVbias->Draw();
    sprintf( name, "leftVbias%1d.gif", sortMethod );
    cA->Print( name );

    leftVmax->Draw();
    sprintf( name, "leftVmax%1d.gif", sortMethod );
    cA->Print( name );

    leftDomain->Draw();
    sprintf( name, "leftDomain%1d.gif", sortMethod );
    cA->Print( name );

    leftItot->Draw();
    sprintf( name, "leftItot%1d.gif", sortMethod );
    cA->Print( name );

    leftDeadp->Draw();
    sprintf( name, "leftDeadp%1d.gif", sortMethod );
    cA->Print( name );

    leftDeadn->Draw();
    sprintf( name, "leftDeadn%1d.gif", sortMethod );
    cA->Print( name );

    leftDomainVbias->Draw("col");
    sprintf( name, "lefttDomainVbias%1d.gif", sortMethod );
    cA->Print( name );

  }

  
  // End

  delete db;

}