diff --git a/Be/BeAnalysis.cpp b/Be/BeAnalysis.cpp
index 2e1e7c6923b0285b62ad9327dc5eafd04ac71c37..06e8b5e58129afabbcfc29aec340748e92a60e48 100644
--- a/Be/BeAnalysis.cpp
+++ b/Be/BeAnalysis.cpp
@@ -46,6 +46,13 @@ BeAnalysis::BeAnalysis() : che(0), spectra(1) {
 	epsilonY[0] = 1;
 	epsilonY[2] = 1;
 
+	cosThetaY = new Bool_t[noIntervals];
+	for (Int_t i = 0; i < noIntervals; i++) {
+		cosThetaY[i] = 0;
+	}
+	cosThetaY[0] = 1;
+	cosThetaY[2] = 1;
+	cosThetaY[4] = 1;
 
 	kVerbose = 1;
 	kRangeProportionIn = 2.0;
@@ -676,6 +683,91 @@ void BeAnalysis::EpsilonY() {
 }
 
 void BeAnalysis::CosThetaYk() {
+
+	TCanvas *cThetaY[noIntervals];
+		for (Int_t i = 0; i<noIntervals; i++) {	//energy intervals (canvases)
+
+		TH1F *hscoskY[noIntervals][6];
+		TH1F *hecoskY[noIntervals][6];
+		TH1F *hicoskY[noIntervals][6];
+
+		if (cosThetaY[i]) continue;
+		//			{
+
+		cThetaY[i] = new TCanvas();
+
+		canvasTitle.Form("cos thetak in %s Y-system;\t(%d,%d) degrees", (const char*)cBeE[i], kMinAngle, kMaxAngle);
+		cThetaY[i]->SetTitle(canvasTitle.Data());
+		cThetaY[i]->Divide(2, 3);
+
+		che->SetLineColor(kBlack);
+		for (Int_t j = 0; j < 6; j++) {	//different files
+			cThetaY[i]->cd(j+1);
+			chs[j]->SetLineColor(kGray+1);
+			chs[j]->SetFillColor(kGray+1);
+			hsName.Form("hscoskY%d_%d", i, j);
+			drawCommand.Form("fCosThetaYk>>%s(50,-1,1)", hsName.Data());
+			chs[j]->Draw(drawCommand.Data(), cQ && crBeE[i] && crAngles && sRatio[j], "", sEvents[j]);
+			hscoskY[i][j] = (TH1F*)gPad->FindObject(hsName.Data());
+			heName.Form("hecoskY%d_%d", i, j);
+			drawCommand.Form("fCosThetaYk>>%s(50,-1,1)", heName.Data());
+			che->Draw(drawCommand.Data(), cQ && cBeE[i] && cAngles, "E same", eEvents[j]);
+//			TH1F *hecoskY[i][j] = (TH1F*)gPad->FindObject(heName.Data());
+			hecoskY[i][j] = (TH1F*)gPad->FindObject(heName.Data());
+
+			hscoskY[i][j]->Draw();
+			//				hscoskY[i][j]->GetYaxis()->SetRangeUser(0, coskYRange[i][j]);
+			hecoskY[i][j]->Draw("E same");
+			if (kAutoRange) {
+				Float_t leftMaxMC = kRangeProportion*hecoskY[i][j]->GetMaximum();
+				Float_t leftMaxE = kRangeProportion*hscoskY[i][j]->GetMaximum();
+				hscoskY[i][j]->GetYaxis()->SetRangeUser(0, leftMaxE);
+				leftMaxMC > leftMaxE ? hscoskY[i][j]->GetYaxis()->SetRangeUser(0, leftMaxMC) : hscoskY[i][j]->GetYaxis()->SetRangeUser(0, leftMaxE);
+			}
+			else {
+				hscoskY[i][j]->GetYaxis()->SetRangeUser(0, eTRange[i][j]);
+			}
+
+
+
+			cThetaY[i]->Update();
+
+			hiName.Form("hicoskY%d_%d", i, j);
+			drawCommand.Form("sCosThetaYk>>%s", hiName.Data());
+			ti[j]->Draw(drawCommand.Data(), ciBeE[i] && sRatio[j], "same");
+			hicoskY[i][j] = (TH1F*)gPad->FindObject(hiName.Data());
+			Float_t rightmax = 1.1*hicoskY[i][j]->GetMaximum();
+			Float_t scale = cThetaY[i]->GetPad(j+1)->GetUymax()/rightmax;
+			hicoskY[i][j]->SetLineColor(kRed);
+			hicoskY[i][j]->Scale(scale);
+
+			//draw an axis on the right side
+			TGaxis *axis = new TGaxis(gPad->GetUxmax(),gPad->GetUymin(),
+					gPad->GetUxmax(), gPad->GetUymax(),0,rightmax,510,"+L");
+			axis->SetLineColor(kRed);
+			axis->SetLabelColor(kRed);
+			axis->Draw();
+
+			cThetaY[i]->Update();
+
+			if (kVerbose) {
+				Info("sfAngInt_cosThetaYk.cxx", "cut%d; case%d: %3.1f (exp)/ %3.1f (sim) = %3.3f",
+						i, j, hecoskY[i][j]->Integral(0,hecoskY[i][j]->GetNbinsX()),
+						hscoskY[i][j]->Integral(0,hscoskY[i][j]->GetNbinsX()),
+						hecoskY[i][j]->Integral(0,hecoskY[i][j]->GetNbinsX())/hscoskY[i][j]->Integral(0,hscoskY[i][j]->GetNbinsX()) );
+			}
+
+
+		}//for j
+		/*if (savePictures) {
+					canvasName.Form("%sfig%dcosThetakY:%s%s", ppath.Data(), i, configuration.Data(), pictFormat.Data());
+					cThetaY[i]->SaveAs(canvasName.Data());
+					cThetaY[i]->cd();
+					cThetaY[i]->Close();
+				}*/
+		//			}//if
+		}//for i
+
 }
 
 void BeAnalysis::ThetaP1T() {
diff --git a/Be/BeAnalysis.h b/Be/BeAnalysis.h
index c9caa17568b11d719e564d7746d1352a9fe7ea55..0f8519156fdf4c8e1595e2e2f17e77982b90d668 100644
--- a/Be/BeAnalysis.h
+++ b/Be/BeAnalysis.h
@@ -69,6 +69,7 @@ private:
 	Bool_t *epsilonT;	//!
 	Bool_t *cosThetaT;	//!
 	Bool_t *epsilonY;	//!
+	Bool_t *cosThetaY;	//!
 
 	Bool_t kAutoRange;
 	Float_t kRangeProportion;
diff --git a/macros/BeCorrPRC/lib_test.cxx b/macros/BeCorrPRC/lib_test.cxx
index 7a8f4015beb5fd97f1eeb526ccaacc94a89016bc..62f04aac5bdfd27c7ab21730fa7010d28d084a43 100755
--- a/macros/BeCorrPRC/lib_test.cxx
+++ b/macros/BeCorrPRC/lib_test.cxx
@@ -16,9 +16,10 @@ void lib_test() {
 	ana.ExpEventsECuts();
 	ana.SimEventsECuts();
 	ana.EpsilonTRange();
-//	ana.EpsilonT();
 
-	ana.CosThetaTk();
+//	ana.EpsilonT();
+//	ana.CosThetaTk();
 
 	ana.EpsilonY();
+	ana.CosThetaYk();
 }