ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/yangyong/combineICEB/getCalibConstants.cc

Comparing UserCode/yangyong/combineICEB/getCalibConstants.cc (file contents):
Revision 1.2 by yangyong, Tue Apr 10 19:41:08 2012 UTC vs.
Revision 1.3 by yangyong, Mon Aug 27 21:37:39 2012 UTC

#	Line 30 | Line 30 | void scaleMeanToUnit(float C[170][360]){
30
31
32
33	–	// void scaleMeanToUnitEndcap(float C[2][101][101]){
34	–
35	–
36	–	//   double mean = 0;
37	–	//   int n = 0;
38	–	//   for(int j=0; j< 170; j++){
39	–	//     for(int k=0; k< 360; k++){
40	–
41	–	//       if( C[j][k] >0){
42	–	//      mean += C[j][k];
43	–	//      n ++;
44	–	//       }
45	–	//     }
46	–	//   }
47	–	//   mean /= n;
48	–	//   cout<<"mean " << mean<<" "<<n<<endl;
49	–
50	–	//   for(int j=0; j< 170; j++){
51	–	//     for(int k=0; k< 360; k++){
52	–	//       if(C[j][k] >0){
53	–	//      C[j][k] /= mean;
54	–	//       }
55	–	//     }
56	–	//   }
57	–
58	–	// }
59	–
33
34
35
#	Line 389 | Line 362 | void readInterCalibConstEBSimple(const c
362
363
364		txtin1 >> eta >> phi >> cc;
392	–
393	–	int iphi = phi;
394	–	int ieta = eta;
395	–
396	–	// cout<< eta <<" "<< phi <<" "<< cc <<endl;
397	–
365		convxtalid(phi,eta);
366		eta += 85;
367		C[eta][phi] = cc;
#	Line 429 | Line 396 | void readInterCalibConstEBSimplev1(const
396		while(txtin1.good()){
397
398		txtin1 >> eta >> phi >> cc;
432	–
433	–	int iphi = phi;
434	–	int ieta = eta;
399
400		if( !(eta>=0 && eta <= 169 && phi >=0 && phi<=359)) {
401		cout<<"wrong eta/phi " << eta<<" "<<phi <<endl;
#	Line 461 | Line 425 | void readInterCalibConstEBSimplev1(const
425
426
427
464	–	void readInterCalibBarre_Phiv1new(){
465	–
466	–	char *file_input = new char[500];
467	–	sprintf(file_input,"phisym_EcalIntercalib-160M-corrected.txt");
468	–
469	–	sprintf(file_input,"EcalIntercalibConstants_Run2010A_corr.xml.txt");
470	–
471	–	cout<<"READING fro inter-calib file "<<file_input<<endl;
472	–	ifstream inputcc(file_input,ios::in);
473	–	if (inputcc.fail()){
474	–	cout<<"error open file barrel.. " << file_input<<endl;
475	–	exit(1);
476	–	}
477	–
478	–	float c1;
479	–	for(int ieta=-85; ieta <= 85; ieta++){
480	–	if( ieta==0) continue;
481	–	for(int iphi=1; iphi <= 360; iphi++){
482	–	inputcc >> c1;
483	–
484	–	int eta = ieta;
485	–	int phi = iphi;
486	–	convxtalid(phi,eta);
487	–
488	–	if( c1 <=0){
489	–	c1 =1;
490	–	}
491	–
492	–	if( c1 == 1 ){
493	–	c1 = - 1;
494	–	}
495	–	if( c1 <0.6 || c1 >2.5){
496	–	c1 = -1;
497	–	}
498	–
499	–
500	–	CphiCorr[eta+85][phi] = c1;
501	–	}
502	–	}
503	–
504	–	}
505	–
506	–
507	–	void readInterCalibBarre_Phi2011a(){
508	–
509	–	char *file_input = new char[500];
510	–
511	–
512	–	sprintf(file_input,"EcalIntercalib-2011A-corrected.xmlv1");
513	–
514	–	cout<<"READING fro inter-calib file "<<file_input<<endl;
515	–	ifstream inputcc(file_input,ios::in);
516	–	if (inputcc.fail()){
517	–	cout<<"error open file barrel.. " << file_input<<endl;
518	–	exit(1);
519	–	}
520	–
521	–	float c1;
522	–	for(int ieta=-85; ieta <= 85; ieta++){
523	–	if( ieta==0) continue;
524	–	for(int iphi=1; iphi <= 360; iphi++){
525	–	inputcc >> c1;
526	–
527	–	int eta = ieta;
528	–	int phi = iphi;
529	–	convxtalid(phi,eta);
530	–
531	–	if( c1 <=0){
532	–	c1 =1;
533	–	}
534	–
535	–	if( c1 == 1 ){
536	–	c1 = - 1;
537	–	}
538	–	if( c1 <0.6 || c1 >2.5){
539	–	c1 = -1;
540	–	}
541	–
542	–
543	–	CphiCorr[eta+85][phi] = c1;
544	–	}
545	–	}
546	–
547	–
548	–
549	–	}
550	–
551	–
552	–
553	–
554	–	void readInterCalibBarre_Phiv1(){
555	–
556	–	char *file_input = new char[500];
557	–	sprintf(file_input,"phisym_EcalIntercalib-160M-corrected.txt");
558	–
559	–
560	–	cout<<"READING fro inter-calib file "<<file_input<<endl;
561	–	ifstream inputcc(file_input,ios::in);
562	–	if (inputcc.fail()){
563	–	cout<<"error open file barrel.. " << file_input<<endl;
564	–	exit(1);
565	–	}
566	–
567	–	int eta;
568	–	int phi;
569	–	float c = 1;
570	–	float c1 = 1;
571	–
572	–	int n = 0;
573	–
574	–	eta  = -85;
575	–	phi = 0;
576	–
577	–	float mean = 0;
578	–
579	–	int cureta = -85;
580	–
581	–	int ngood = 0;
582	–
583	–
584	–
585	–	while (inputcc.good()){
586	–
587	–	inputcc >> c1;
588	–
589	–	if(c1==1) c1 = -1;
590	–
591	–
592	–	if(  n% 360 ==0 && n>0 ) {
593	–	eta++;
594	–	if(eta==0) eta = 1;
595	–	}
596	–
597	–	if( eta == cureta){
598	–	phi ++;
599	–	}else{
600	–	cureta = eta;
601	–	phi = 1;
602	–
603	–	}
604	–
605	–
606	–	n++;
607	–
608	–
609	–	int iSM = (phi-1)/20+1;
610	–	if( eta<0) iSM += 18;
611	–
612	–
613	–	int ieta = eta;
614	–	int iphi = phi;
615	–
616	–
617	–	convxtalid(iphi,ieta);
618	–
619	–	///txtout<<eta<<" "<<phi<<endl;
620	–
621	–	CphiCorr[ieta+85][iphi] = c1;
622	–
623	–
624	–	if(c1 >0 ){
625	–	mean += c1;
626	–	ngood ++;
627	–	}
628	–
629	–
630	–	if(n >= 360*170) break;
631	–
632	–	}
633	–
634	–
635	–	mean /=  ngood;
636	–	cout<<n<<" phi  constants read. "<< mean <<endl;
637	–
638	–
639	–	//   for(int j=0; j< 170; j++){
640	–	//     mean_ieta[j] /= nmean_ieta[j];
641	–	//     txtout<<"mean_ieta_phisym "<< j<<" "<< mean_ieta[j] <<" "<< nmean_ieta[j] <<endl;
642	–	//   }
643	–
644	–
645	–	//   for(int j=0; j< 170; j++){
646	–	//     for(int k=0; k< 360; k++){
647	–
648	–	//       if(  CphiCorr[j][k] >0 ){
649	–	//      CphiCorr[j][k] /=  mean_ieta[j] ;
650	–
651	–	//       }
652	–	//       Cphi[j][k] =   CphiCorr[j][k] *  iCalConst_data[j][k];
653	–	//     }
654	–	//   }
655	–
656	–	//   NormSMScaleToUnit(CphiCorr);
657	–
658	–	}
659	–
660	–
661	–	///read initi calibration constants for each crystal
662	–	void readInterCalibBarre_BSv1(){
663	–
664	–
665	–	char *file_input = new char[500];
666	–	sprintf(file_input,"IC_barrel_splash09_AllRuns.txt");
667	–
668	–
669	–	cout<<"READING fro inter-calib bsfile "<<file_input<<endl;
670	–	txtout<<"READING fro inter-calib bsfile "<<file_input<<endl;
671	–	ifstream inputcc(file_input,ios::in);
672	–	if (inputcc.fail()){
673	–	cout<<"error open file barrel.. " << file_input<<endl;
674	–	exit(1);
675	–	}
676	–
677	–	int eta;
678	–	int phi;
679	–	float c = 1;
680	–	float c1 = 1;
681	–
682	–	int n = 0;
683	–	int ngood = 0;
684	–
685	–
686	–	float mean = 0;
687	–
688	–
689	–
690	–
691	–	float mean_ieta[170] = {0};
692	–	int nmean_ieta[170] = {0};
693	–
694	–
695	–	while (inputcc.good()){
696	–
697	–	inputcc >> eta >> phi >>c1;
698	–
699	–	//  cout<<phi<<" "<<eta<<" "<<c<<endl;
700	–
701	–	convxtalid(phi,eta);
702	–
703	–	eta += 85;
704	–
705	–	if( c1 == 1){
706	–	c1 = -1;
707	–	}
708	–
709	–	if( c1 <0.6 || c1 >2.5){
710	–	c1 = -1;
711	–	}
712	–
713	–
714	–	///iCalConstPi0[eta][phi] = c1 ;
715	–
716	–	CBSCorr[eta][phi] = c1;
717	–
718	–	if(c1 >0 ){
719	–	mean += c1;
720	–	ngood ++;
721	–	}
722	–
723	–	if(c1 >0 ){
724	–	mean_ieta[eta] += c1;
725	–	nmean_ieta[eta] ++;
726	–	}
727	–
728	–	n++;
729	–
730	–	/// hh_imcEB->Fill(c);
731	–	if(n >= 360*170) break;
732	–
733	–	}
734	–
735	–
736	–
737	–	mean /= ngood;
738	–
739	–
740	–
741	–	cout<<n<<" bs  constants read. "<< mean << endl;
742	–
743	–
744	–
745	–	//   for(int j=0; j< 170; j++){
746	–	//     for(int k=0; k< 360; k++){
747	–
748	–	//       if(  CBSCorr[j][k] >0 ){
749	–	//      CBSCorr[j][k] /=  mean_ieta[j];
750	–	//       }
751	–
752	–	//       CBS[j][k] = CBSCorr[j][k] * iCalConst_data[j][k];
753	–
754	–
755	–	//     }
756	–	//   }
757	–	//   NormSMScaleToUnit(CBSCorr);
758	–
759	–
760	–
761	–	}
762	–
763	–
764	–
765	–
766	–
767	–	///read initi calibration constants for each crystal
768	–	void readInterCalibBarre_GR09_V8(){
769	–
770	–
771	–
772	–	char *file_input = new char[500];
773	–	//file_input ="/data/yangyong/data/p
774	–	file_input = "/afs/cern.ch/cms/cit/yongy/interCalib_GR09_R_V6A_barrel.txt";
775	–
776	–	///sprintf(file_input,"/uscms_data/d2/yongy/calib/backup/interCalib_STARTUP3X_V8P_barrel_mc.txt");  ///now the same used in data
777	–
778	–	//sprintf(file_input,"/uscms/home/yongy/work/interCalib_STARTUP3X_V8P_barrel.txt");
779	–
780	–	//  sprintf(file_input,"/uscms/home/yongy/work/interCalib_GR_R_36X_V11A_Barrel.txt");
781	–
782	–	cout<<"READING fro inter-calib file "<<file_input<<endl;
783	–	ifstream inputcc(file_input,ios::in);
784	–	if (inputcc.fail()){
785	–	cout<<"error open file barrel.. " << file_input<<endl;
786	–	exit(1);
787	–	}
788	–
789	–	int eta;
790	–	int phi;
791	–	float c = 1;
792	–	float c1 = 1;
793	–
794	–	int n = 0;
795	–
796	–	double mean = 0;
797	–
798	–
799	–	while (inputcc.good()){
800	–
801	–	inputcc >> eta >> phi >>c1;
802	–
803	–	//  cout<<phi<<" "<<eta<<" "<<c<<endl;
804	–
805	–	convxtalid(phi,eta);
806	–
807	–	eta += 85;
808	–
809	–	///iCalConst[eta][phi] = c;
810	–
811	–	//iCalConst_data[eta][phi] = c1 * 0.038012/ 0.038;  ///Now this is exactly the same as the inter-calibration used in data.
812	–
813	–	interCalib_preCalib[eta][phi] = c1 ; /// now OK iN MC startupv8p.
814	–
815	–
816	–	n++;
817	–	//  cout<<phi<<" "<<eta<<" "<<c<<endl;
818	–	mean += c1;
819	–
820	–	/// hh_imcEB->Fill(c);
821	–	if(n >= 360*170) break;
822	–
823	–	}
824	–
825	–
826	–	mean /= 61200;
827	–
828	–	cout<<n<<" constants read. "<< mean <<endl;
829	–
830	–
831	–	}
832	–
833	–
834	–
835	–
836	–
837	–
838	–
839	–
840	–	///read initi calibration constants for each crystal
841	–	void readInterCalibEndcap_GR09_V8(){
842	–	char *file_input = new char[500];
843	–	//file_input ="/data/yangyong/data/p
844	–	///file_input = "/afs/cern.ch/cms/cit/yongy/interCalib_GR09_R_V6A_barrel.txt";
845	–
846	–	file_input = "/afs/cern.ch/cms/cit/yongy/interCalib_GR09_H_V6OFF_endcap.txt";
847	–
848	–	cout<<"READING fro inter-calib file "<<file_input<<endl;
849	–	ifstream inputcc(file_input,ios::in);
850	–	if (inputcc.fail()){
851	–	cout<<"error open file barrel.. " << file_input<<endl;
852	–	exit(1);
853	–	}
854	–
855	–	int eta;
856	–	int phi;
857	–	int iz;
858	–	float c = 1;
859	–	int n = 0;
860	–
861	–
862	–	while (inputcc.good()){
863	–
864	–	inputcc >> iz >> eta >> phi  >> c;
865	–
866	–	int izz = iz < 0? 0:1;
867	–
868	–	interCalibEndcap_preCalib[izz][eta][phi] = c;
869	–	validRecHitEndCap[izz][eta][phi] = 1;
870	–
871	–	n ++;
872	–
873	–	}
874	–
875	–	cout<<n<<" constants read. "<<endl;
876	–
877	–	}
878	–
879	–
880	–	void getInterCalibEndcapv1(const char *inputfile,float C[2][101][101]){
881	–
882	–
883	–
884	–	ifstream inputcc(inputfile);
885	–	if (inputcc.fail()) {
886	–	cout<<"error "<< inputfile<<" can not be opened."<<endl;
887	–	exit(1);
888	–	}
889	–	int eta;
890	–	int phi;
891	–	int iz;
892	–	float c = 1;
893	–	int n = 0;
894	–
895	–	double mean = 0;
896	–
897	–	while (inputcc.good()){
898	–
899	–	inputcc >> iz >> eta >> phi  >> c;
900	–
901	–	int izz = iz < 0? 0:1;
902	–
903	–	C[izz][eta][phi] = c;
904	–	mean += c;
905	–
906	–	n ++;
907	–
908	–	if( n>= 14648){
909	–	break;
910	–	}
911	–
912	–	}
913	–
914	–	cout<<n<<" constants read. "<< mean /n <<endl;
915	–
916	–	}
917	–
918	–
919	–
920	–	void getInterCalibEndcap(const char *inputfile, float C[2][101][101]){
921	–
922	–
923	–	ifstream inputcc(inputfile);
924	–	if (inputcc.fail()) {
925	–	cout<<"error "<< inputfile<<" can not be opened."<<endl;
926	–	exit(1);
927	–	}
928	–
929	–	int ix,iy,iz;
930	–	float cc;
931	–
932	–	float ccMeaniring[50] ;
933	–	int nccMeaniring[50];
934	–
935	–
936	–	//  for(int j=0; j<2;j++){
937	–	for(int k=0; k< kEndcEtaRings; k++){
938	–	ccMeaniring[k] = 0;
939	–	nccMeaniring[k] = 0;
940	–	}
941	–	//  }
942	–
943	–	int nread = 0;
944	–	while( inputcc.good()){
945	–	inputcc >> iz >> ix >> iy >> cc;
946	–	C[iz][ix][iy] = cc;
947	–	int iring = iRingEndCap(2*iz-1,ix,iy); ///input -1/1,
948	–	if( cc > 0){
949	–	ccMeaniring[iring] += cc;
950	–	nccMeaniring[iring] += 1;
951	–	}
952	–
953	–	nread ++;
954	–	if( iz==1 && ix == 100 && iy == 60) break;
955	–	}
956	–
957	–	cout<<"nb of crystals endcap " << nread <<endl;
958	–
959	–
960	–	//for(int j=0; j<2;j++){
961	–	for(int k=0; k< kEndcEtaRings; k++){
962	–	ccMeaniring[k] /= nccMeaniring[k];
963	–	//  }
964	–	}
965	–
966	–	double meanall = 0;
967	–	int nall = 0;
968	–	for(int iz=0; iz<2; iz++){
969	–	for(int j=0; j<101; j++){
970	–	for(int k=0; k< 101; k++){
971	–	if( validRecHitEndCap[iz][j][k] <1) continue;
972	–	int iring = iRingEndCap(2*iz-1,j,k); ///input -1/1,
973	–	if( C[iz][j][k] >0) {
974	–	C[iz][j][k] /= ccMeaniring[iring];
975	–	meanall += C[iz][j][k];
976	–	nall ++;
977	–	}
978	–
979	–	}
980	–	}
981	–	}
982	–
983	–	meanall /= nall;
984	–	cout<<" corrfactoriZiXiY read " << meanall <<" "<< nall <<endl;
985	–
986	–
987	–	}
988	–
989	–
990	–
991	–	void NormCrystalDeadFlagEndcap_v1(float C[2][101][101],int ndeadflag[2][101][101]){
992	–
993	–	float mean_deadflag[20] = {0};
994	–	float nmean_deadflag[20] = {0};
995	–
996	–	int ndeadCrystals = 0;
997	–
998	–
999	–	for(int iz =0; iz < 2; iz++){
1000	–	for(int j=0; j< 101; j++){
1001	–	for(int k=0; k<101; k++){
1002	–
1003	–	if( validRecHitEndCap[iz][j][k] ==0) continue;
1004	–	if( C[iz][j][k] <0){
1005	–	continue;
1006	–	}
1007	–	int ndead = ndeadflag[iz][j][k];
1008	–	if( ndead >=0){
1009	–	mean_deadflag[ndead] += C[iz][j][k];
1010	–	nmean_deadflag[ndead] ++;
1011	–	}
1012	–
1013	–	}
1014	–	}
1015	–	}
1016	–	cout<<" NormCrystalDeadFlag ndeadCrystals: " << ndeadCrystals <<endl;
1017	–	for(int j=0; j<20; j++){
1018	–	if( nmean_deadflag[j] >0){
1019	–	mean_deadflag[j] /= nmean_deadflag[j];
1020	–	cout<<" mean_deadflag "<<j<<"  "<< mean_deadflag[j]<< " "<< nmean_deadflag[j] <<endl;
1021	–	}
1022	–	}
1023	–
1024	–
1025	–	for(int iz =0; iz < 2; iz++){
1026	–	for(int j=0; j< 101; j++){
1027	–	for(int k=0; k<101; k++){
1028	–
1029	–	if( validRecHitEndCap[iz][j][k] ==0) continue;
1030	–	if( C[iz][j][k] <0){
1031	–	continue;
1032	–	}
1033	–
1034	–	int ndead = ndeadflag[iz][j][k];
1035	–	if( ndead >=0){
1036	–	C[iz][j][k] /= mean_deadflag[ndead];
1037	–	}
1038	–
1039	–	}
1040	–	}
1041	–	}
1042	–
1043	–
1044	–
1045	–	}
1046	–
1047	–
1048	–
1049	–
1050	–
1051	–
1052	–	void NormCrystalDeadFlagEndcap(float C[2][101][101]){
1053	–
1054	–	float mean_deadflag[20] = {0};
1055	–	float nmean_deadflag[20] = {0};
1056	–
1057	–	int ndeadCrystals = 0;
1058	–
1059	–
1060	–	for(int iz =0; iz < 2; iz++){
1061	–	for(int j=0; j< 101; j++){
1062	–	for(int k=0; k<101; k++){
1063	–
1064	–	if( validRecHitEndCap[iz][j][k] ==0) continue;
1065	–	if( C[iz][j][k] <0){
1066	–	continue;
1067	–	}
1068	–	int ndeadflag = ndeadflag_endcap[iz][j][k];
1069	–	if( ndeadflag >=0){
1070	–	mean_deadflag[ndeadflag] += C[iz][j][k];
1071	–	nmean_deadflag[ndeadflag] ++;
1072	–	}
1073	–
1074	–	}
1075	–	}
1076	–	}
1077	–	cout<<" NormCrystalDeadFlag ndeadCrystals: " << ndeadCrystals <<endl;
1078	–	for(int j=0; j<20; j++){
1079	–	if( nmean_deadflag[j] >0){
1080	–	mean_deadflag[j] /= nmean_deadflag[j];
1081	–	cout<<" mean_deadflag "<<j<<"  "<< mean_deadflag[j]<< " "<< nmean_deadflag[j] <<endl;
1082	–	}
1083	–	}
1084	–
1085	–
1086	–	for(int iz =0; iz < 2; iz++){
1087	–	for(int j=0; j< 101; j++){
1088	–	for(int k=0; k<101; k++){
1089	–
1090	–	if( validRecHitEndCap[iz][j][k] ==0) continue;
1091	–	if( C[iz][j][k] <0){
1092	–	continue;
1093	–	}
1094	–
1095	–	int ndeadflag = ndeadflag_endcap[iz][j][k];
1096	–	if( ndeadflag >=0){
1097	–	C[iz][j][k] /= mean_deadflag[ndeadflag];
1098	–	}
1099	–
1100	–	}
1101	–	}
1102	–	}
1103	–
1104	–
1105	–
1106	–	}
1107	–
1108	–
428
1110	–	void scaleMeanToUnitEndcap(float C[2][101][101]){
1111	–
1112	–
1113	–	double mean = 0;
1114	–	int n = 0;
1115	–
1116	–
1117	–	for(int iz =0; iz < 2; iz++){
1118	–	for(int j=0; j< 101; j++){
1119	–	for(int k=0; k<101; k++){
1120	–	if( validRecHitEndCap[iz][j][k] ==0) continue;
1121	–	if( C[iz][j][k] <0){
1122	–	continue;
1123	–	}
1124	–	mean += C[iz][j][k];
1125	–	n ++;
1126	–	}
1127	–	}
1128	–	}
1129	–
1130	–	mean /= n;
1131	–	cout<<"mean " << mean<<" "<<n<<endl;
1132	–
1133	–
1134	–	for(int iz =0; iz < 2; iz++){
1135	–	for(int j=0; j< 101; j++){
1136	–	for(int k=0; k<101; k++){
1137	–	if( validRecHitEndCap[iz][j][k] ==0) continue;
1138	–	if( C[iz][j][k] <0){
1139	–	continue;
1140	–	}
1141	–
1142	–	C[iz][j][k] /= mean;
1143	–	}
1144	–	}
1145	–	}
1146	–
1147	–	}

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines