{

	// Plot the conditional probability P(p|K,N)
	// N is the number of trials
  // p is the success probability
  // K is the number of success

	gROOT->LoadMacro("distributionEfficiencyFromOccurence.C");
	
	Int_t N = 400;
	Int_t nK = 7; // number of cases to plot
	
	Int_t Kmin = 0;
	Int_t Kmax = N;
	Int_t K; // number of successes

	// Function for the conditional probability
	TF1 fp( "fp", distributionEfficiencyFromOccurence, 0., 1., 2);
	fp.SetLineWidth(1.5);
	fp.SetParameter(0, N);

	// Histogram for the conditional probability
	TH1F *hp;
	
	gStyle->SetOptStat(0);
	TCanvas c("c","Conditional probability of efficiency", 600, 600);
	c.SetGrid(1,1);
	
	Char_t title[200];
	sprintf( title, "Conditional probability of p with K, N=%d", N);
	TH1F h("h", title, 100, -0.2, 1.2);
	h.SetXTitle("p");
	h.SetMinimum(0);
	fp.SetParameter( 1, Kmin);
	h.SetMaximum( fp.Eval(0)/2. );
	h.Draw();
	
	Double_t x[1], p[2];
	p[0] = N;
	Char_t text[20], name[10];
	TLatex l;
	l.SetTextSize(.022);
	for( Int_t iK=0; iK<nK; iK++) {
		K = Kmin + iK*(Kmax-Kmin)/(nK-1);
		p[1] = K;
		cout << endl << "K= " << K << endl;
		sprintf( text, "K=%d", K);
		sprintf( name, "hp%d", K);
		hp = new TH1F( name, text, 100, 0., 1.);
		hp->SetXTitle("p");
		//hp->SetLineColor(iK+1);
		for( Short_t iBin=1; iBin<=hp->GetNbinsX(); iBin++) {
			x[0] = (Double_t)iBin/hp->GetNbinsX();
			hp->SetBinContent( iBin, distributionEfficiencyFromOccurence(x,p)/hp->GetNbinsX() );
			//cout << x[0] << " = " << distributionEfficiencyFromOccurence(x,p)/hp->GetNbinsX() << endl;
		}
		if( iK==0) hp->SetMaximum(.4);
		hp->Draw( iK?"Lsame":"L");
		//fp.SetParameter(1, K);
		//fp.SetRange( TMath::Max(0., (Double_t)K/N-.2), TMath::Min(1., (Double_t)K/N+.2) );
		//cout << "Plotting " << K << "Range " << TMath::Max(0., (Double_t)K/N-.2) << " " << TMath::Min(1., (Double_t)K/N+.2) << endl;
		//fp.DrawClone("same");
		l.DrawLatex( (Float_t)K/N, hp->GetMaximum(), text);
		//if( K==0 || K==N ) 		l.DrawLatex( (Float_t)K/N, fp.Eval((Float_t)K/N)/2., text);
	}
	cout << "end" << endl;
		
}