#include <iostream.h>
//
// chips_per_board should be 2 for SL1, 3 for SL2, 4 for SL3, 5 for SL4 and 12 for linker boards
// xflag should be 0 when bins are in chips, 1 when bins are in COT cells, 2 when bins are in COT wires   
// Sensitvity allows you to change the level at which cold/hot spots are reported!

int xft_hotcold(TH1* h1, Int_t chips_per_board=1, Int_t xflag=0, Double_t sensitivity =10.)
{	
	//Some constants: useful for going from histogram bin to board
	Int_t boards_per_2pi = 24;
	Int_t board_offset_wrt_0pi = 6;
	Int_t wires_per_cell = 12;
	Int_t cells_per_chip =  4;
	Int_t wires_per_chip = wires_per_cell * cells_per_chip;	

	Int_t chip_offset = 1; //Chip 0 is finder chip 1B, 2B, 4B
	if(chips_per_board == 4)  chip_offset = 2; //Chip 0 is finder chip 3C
	if(chips_per_board == 12) chip_offset = 0; //Chip 0 is linker chip 0

	//Print the histogram title and fit the occupancy histogram with a flat line
	Char_t *title = h1->GetTitle();
	cout << " " << endl;
	cout << "**********************************************************************************************************" << endl;
	cout << "*** " << title << " ***" << endl;

	h1->Fit("pol0","Q");
	Double_t chi2  = h1->GetFunction("pol0")->GetChisquare();
	Double_t ndf   = h1->GetFunction("pol0")->GetNDF();
	Double_t mean  = h1->GetFunction("pol0")->GetParameter(0);
	Double_t error = h1->GetFunction("pol0")->GetParError(0);
	Double_t sigma = sqrt(fabs(mean));

	cout << "----------------------------------------------------------------------------------------------------------" << endl;
	cout << "Mean occupancy "<< mean << " +- "<< error;
	cout << "   Chi2/dof "<< chi2 << " / " << ndf << " = "<< chi2/ndf << endl;

	cout << "List bins outside ( " << (mean-sensitivity*sigma) << ", " << (mean + sensitivity*sigma) << " )"
		  << ", ie "<< sensitivity << " * " << sigma << " distant" << endl;
	cout << "----------------------------------------------------------------------------------------------------------" << endl;
	
	Int_t nx = h1->GetNbinsX();
	for (Int_t i=0; i < nx; i++){ // Note well that i=0 is underflow, i=nx+1 is overflow
		Double_t occupancy = h1->GetBinContent(i+1);
		Double_t delta = occupancy-mean; 
		if(fabs(delta) > sensitivity * sigma){
			if(delta > 0.){
				printf("\n HOT!   Bin %4d  Occupancy %8d  sigma away %6.1f",i,occupancy,delta/sigma);
 			}else if (delta < 0.){
				printf("\n COLD!  Bin %4d  Occupancy %8d  sigma away %6.1f",i,occupancy,delta/sigma);
			}

			//Calculate board, chip (cell, wire) from histogram bin number
			Int_t absolute_cell = i*cells_per_chip;
			Int_t chipBin = i;
			Int_t wire    = 0;
			Int_t cell    = 0;
			
			if(xflag == 1){ 							// Bins are COT cells
				absolute_cell = i;
				cell    = i%cells_per_chip;
				chipBin = i/cells_per_chip;
			}else if (xflag == 2){  				// Bins are COT wires	
				absolute_cell =  i/wires_per_cell;
				wire    =  i%wires_per_cell;
				cell    = (i/wires_per_cell) % cells_per_chip;
				chipBin =  i/wires_per_chip; 
			}
			Int_t chip  =   (chipBin + chip_offset) % chips_per_board;
			Int_t board = (((chipBin + chip_offset) / chips_per_board) + board_offset_wrt_0pi) % boards_per_2pi;
			Double_t degrees = absolute_cell * 360. / (cells_per_chip * chips_per_board * boards_per_2pi);
			Double_t radians = degrees * 3.1416 /180.;
			
			printf("   Board %2d-%2d",board,chip);
			if(xflag>0) printf("-%1d",cell);
			if(xflag>1) printf(" wire %2d",wire);
			if(chips_per_board != 12) printf("   Cell %4d",absolute_cell);
			printf(" (%5.1f deg, %5.3f rad)", degrees, radians);
		}
	}	
	return 0;

}
int xft_hotcoldhits(TH1* h1, TH1* h2, TH1* h3, Int_t chips_per_board=1, Int_t xflag=0, Double_t sensitivity =10.)
{	
	//Some constants: useful for going from histogram bin to board
	Int_t boards_per_2pi = 24;
	Int_t board_offset_wrt_0pi = 6;
	Int_t wires_per_cell = 12;
	Int_t cells_per_chip =  4;
	Int_t wires_per_chip = wires_per_cell * cells_per_chip;	

	Int_t chip_offset = 1; //Chip 0 is finder chip 1B, 2B, 4B
	if(chips_per_board == 4)  chip_offset = 2; //Chip 0 is finder chip 3C
	if(chips_per_board == 12) chip_offset = 0; //Chip 0 is linker chip 0

//Print the histogram title and fit the occupancy histogram with a flat line
	Char_t *title = h1->GetTitle();
	cout << " " << endl;
	cout << "**********************************************************************************************************" << endl;
	cout << "*** " << title << " ***" << endl;

	h1->Fit("pol0","Q");
	Double_t chi2   = h1->GetFunction("pol0")->GetChisquare();
	Double_t ndf    = h1->GetFunction("pol0")->GetNDF();
	Double_t mean1  = h1->GetFunction("pol0")->GetParameter(0);
	Double_t error  = h1->GetFunction("pol0")->GetParError(0);
	Double_t sigma1 = sqrt(fabs(mean1));

	cout << "----------------------------------------------------------------------------------------------------------" << endl;
	cout << "Mean occupancy "<< mean1 << " +- "<< error;
	cout << "   Chi2/dof "<< chi2 << " / " << ndf << " = "<< chi2/ndf << endl;

	cout << "List bins outside ( " << (mean1-sensitivity*sigma1) << ", " << (mean1 + sensitivity*sigma1) << " )"
		  << ", ie "<< sensitivity << " * " << sigma1 << " distant" << endl;
	cout << "----------------------------------------------------------------------------------------------------------" << endl;

	//Print the histogram title and fit the occupancy histogram with a flat line
	title = h2->GetTitle();
	cout << "*** " << title << " ***" << endl;

	h2->Fit("pol0","Q");
	chi2  			  = h2->GetFunction("pol0")->GetChisquare();
	ndf   			  = h2->GetFunction("pol0")->GetNDF();
	Double_t mean2   = h2->GetFunction("pol0")->GetParameter(0);
	error 			  = h2->GetFunction("pol0")->GetParError(0);
	Double_t sigma2  = sqrt(fabs(mean2));

	cout << "----------------------------------------------------------------------------------------------------------" << endl;
	cout << "Mean occupancy "<< mean2 << " +- "<< error;
	cout << "   Chi2/dof "<< chi2 << " / " << ndf << " = "<< chi2/ndf << endl;

	cout << "List bins outside ( " << (mean2-sensitivity*sigma2) << ", " << (mean2 + sensitivity*sigma2) << " )"
		  << ", ie "<< sensitivity << " * " << sigma2 << " distant" << endl;
	cout << "----------------------------------------------------------------------------------------------------------" << endl;
	
	//Print the histogram title and fit the occupancy histogram with a flat line
	title = h3->GetTitle();
	cout << "*** " << title << " ***" << endl;

	h3->Fit("pol0","Q");
	chi2  			 = h3->GetFunction("pol0")->GetChisquare();
	ndf   			 = h3->GetFunction("pol0")->GetNDF();
	Double_t mean3  = h3->GetFunction("pol0")->GetParameter(0);
	error 			 = h3->GetFunction("pol0")->GetParError(0);
	Double_t sigma3 = sqrt(fabs(mean3));

	cout << "----------------------------------------------------------------------------------------------------------" << endl;
	cout << "Mean occupancy "<< mean3 << " +- "<< error;
	cout << "   Chi2/dof "<< chi2 << " / " << ndf << " = "<< chi2/ndf << endl;

	cout << "List bins outside ( " << (mean3-sensitivity*sigma3) << ", " << (mean3 + sensitivity*sigma3) << " )"
		  << ", ie "<< sensitivity << " * " << sigma3 << " distant" << endl;
	cout << "----------------------------------------------------------------------------------------------------------" << endl;


	Int_t nx = h1->GetNbinsX();
	for (Int_t i=0; i < nx; i++){ // Note well that i=0 is underflow, i=nx+1 is overflow
		Double_t occupancy1 = h1->GetBinContent(i+1);
		Double_t delta1 = occupancy1-mean1; 
		
		Double_t occupancy2 = h2->GetBinContent(i+1);
		Double_t delta2 = occupancy2-mean2;
		 
		Double_t occupancy3 = h3->GetBinContent(i+1);
		Double_t delta3 = occupancy3-mean3; 
		
		
		if((fabs(delta1) > sensitivity * sigma1) || 
			(fabs(delta2) > sensitivity * sigma2) ||
			(fabs(delta3) > sensitivity * sigma3)){
		
			printf("\n Bin %4d",i);
			printf(" %8d %6.1f",occupancy1,delta1/sigma1);
			if(fabs(delta1) > sensitivity * sigma1){
				if(delta1 > 0.){
					printf(" P HOT  ");
				}else{
					printf(" P COLD ");
				}
			}else{printf("        ");}
			
			printf(" %8d %6.1f",occupancy2,delta2/sigma2);
			if(fabs(delta2) > sensitivity * sigma2){
				if(delta2 > 0.){
					printf(" D HOT  ");
				}else{
					printf(" D COLD ");
				}
			}else{printf("        ");}
			printf(" %8d %6.1f",occupancy3,delta3/sigma3);
			if(fabs(delta3) > sensitivity * sigma3){
				if(delta3 > 0.){
					printf(" PD HOT  ");
				}else{
					printf(" PD COLD ");
				}
			}else{printf("         ");}



			//Calculate board, chip (cell, wire) from histogram bin number
			Int_t absolute_cell = i*cells_per_chip;
			Int_t chipBin = i;
			Int_t wire    = 0;
			Int_t cell    = 0;
			
			if(xflag == 1){ 							// Bins are COT cells
				absolute_cell = i;
				cell    = i%cells_per_chip;
				chipBin = i/cells_per_chip;
			}else if (xflag == 2){  				// Bins are COT wires	
				absolute_cell =  i/wires_per_cell;
				wire    =  i%wires_per_cell;
				cell    = (i/wires_per_cell) % cells_per_chip;
				chipBin =  i/wires_per_chip; 
			}
			Int_t chip  =   (chipBin + chip_offset) % chips_per_board;
			Int_t board = (((chipBin + chip_offset) / chips_per_board) + board_offset_wrt_0pi) % boards_per_2pi;
			Double_t degrees = absolute_cell * 360. / (cells_per_chip * chips_per_board * boards_per_2pi);
			Double_t radians = degrees * 3.1416 /180.;
			
			printf("   Board %2d-%2d",board,chip);
			if(xflag>0) printf("-%1d",cell);
			if(xflag>1) printf(" wire %2d",wire);
			if(chips_per_board != 12) printf("   Cell %4d",absolute_cell);
			printf(" (%5.1f deg, %5.3f rad)", degrees, radians);
		}
	}	
	return 0;

}