22 #include "TGraphErrors.h"
35 int getDVInRun(
int nInRun,
int *runMacro,
unsigned int *runTime,
double *runDVE){
37 const double maxRMS = 0.001;
38 const int minFilesInRun=5;
42 int nToUseInRun = nInRun-1;
43 int iFirstToUseInRun=0;
44 int *runDVEOrderI =
new int[nToUseInRun];
45 double *runDVEOrderV =
new double[nToUseInRun];
47 for(
int iInRun=0; iInRun<nToUseInRun; iInRun++){
48 runDVEOrderV[iInRun]=runDVE[iInRun];
49 runDVEOrderI[iInRun]=iInRun;
56 for(
int iInRun=0; iInRun<nToUseInRun-1; iInRun++){
57 if(runDVEOrderV[iInRun]>runDVEOrderV[iInRun+1]){
58 int tempMinI = runDVEOrderI[iInRun+1];
59 double tempMinV = runDVEOrderV[iInRun+1];
60 runDVEOrderI[iInRun+1]=runDVEOrderI[iInRun];
61 runDVEOrderI[iInRun]=tempMinI;
62 runDVEOrderV[iInRun+1]=runDVEOrderV[iInRun];
63 runDVEOrderV[iInRun]=tempMinV;
69 for(
int iInRun=0; iInRun<nInRun; iInRun++){
70 cout <<
"In run: " << iInRun <<
"/" << nInRun
71 <<
" " << runMacro[iInRun]
72 <<
" " << runTime[iInRun]
73 <<
" " << runDVE[iInRun] << endl;
76 for(
int iInRun=0; iInRun<nToUseInRun; iInRun++){
77 cout <<
"Order: " << iInRun <<
"/" << nToUseInRun
78 <<
" " << runDVEOrderV[iInRun]
79 <<
" " << runDVEOrderI[iInRun] << endl;
84 for(
int iInRun=0; iInRun<nToUseInRun; iInRun++){
85 if(runDVEOrderI[iInRun]==iInRun){
88 if(runDVEOrderI[iInRun]==nToUseInRun-iInRun-1){
93 if(nToUseInRun>=minFilesInRun-1){
94 if(nInOrderU==nToUseInRun || nInOrderD==nToUseInRun){
95 double *meanDiff =
new double[nToUseInRun-2];
96 if(nInOrderU==nToUseInRun){
97 cout <<
"Drift Velocities East are ordered (U) in time. Check convergence ..." << endl;
98 for(
int i=0; i<nToUseInRun-2; i++){
100 for(
int j=i; j<nToUseInRun-1; j++){
101 meanDiff[i]+=(runDVEOrderV[j+1]-runDVEOrderV[j]);
103 meanDiff[i]/=(nToUseInRun-i-1);
106 if(nInOrderD==nToUseInRun){
107 cout <<
"Drift Velocities East are ordered (D) in time. Check convergence ..." << endl;
108 for(
int i=0; i<nToUseInRun-2; i++){
110 for(
int j=nToUseInRun-i-1; j>0; j--){
111 meanDiff[i]+=runDVEOrderV[j-1]-runDVEOrderV[j];
113 meanDiff[i]/=(nToUseInRun-i-1);
116 int nDiffInOrder = 0;
117 for(
int i=0; i<nToUseInRun-3; i++){
118 if(meanDiff[i]>meanDiff[i+1]){
122 if(nDiffInOrder>=nToUseInRun-3-1){
123 if(nInOrderU==nToUseInRun) iFirstToUseInRun++;
124 if(nInOrderD==nToUseInRun) nToUseInRun--;
125 cout <<
"Convergence found! Will remove earliest entry and continue." << endl;
128 cout <<
"No convergence found! This run should be checked. No macro will be kept!" << endl;
135 int iNotValidRMS1 = -1;
136 int iNotValidRMS2 = -1;
138 if(nToUseInRun==minFilesInRun-2){
142 while(tryThreeRMS<3){
147 for(
int iInRun=iFirstToUseInRun; iInRun<nToUseInRun; iInRun++){
148 if(iInRun!=iNotValidRMS1 && iInRun!=iNotValidRMS2){
149 runMean+=runDVEOrderV[iInRun];
150 runRMS+=runDVEOrderV[iInRun]*runDVEOrderV[iInRun];
156 runRMS-=runMean*runMean;
157 runRMS=TMath::Sqrt(runRMS);
159 cout <<
"This run RMS = " << runRMS << endl;
162 int iNotValidLow = -1;
167 while(tryTwiceLow<2){
169 double diffToMean = 100.;
170 for(
int iInRun=iFirstToUseInRun; iInRun<nToUseInRun; iInRun++){
171 double lDiffToMean = TMath::Abs(runDVEOrderV[iInRun]-runMean);
172 if(lDiffToMean<diffToMean && iInRun!=iNotValidLow &&
173 iInRun!=iNotValidRMS1 && iInRun!=iNotValidRMS2){
174 diffToMean=lDiffToMean;
178 int lowThird = nInRun/3;
179 if(runDVEOrderI[iMeanOrder]/lowThird!=0 || (tryTwiceLow==1 && tryThreeRMS==2)){
180 int iKeepMacro = runMacro[runDVEOrderI[iMeanOrder]];
181 cout <<
"Keep Macro: " << iKeepMacro <<
" " << runTime[runDVEOrderI[iMeanOrder]]
182 <<
" " << runDVE[runDVEOrderI[iMeanOrder]] << endl;
183 cout << runDVEOrderI[iMeanOrder] <<
"/" << nInRun <<
" " << iMeanOrder << endl;
189 iNotValidLow=iMeanOrder;
190 cout <<
"Drift Velocity East value closer to mean was found in the first third of the run, is safer to try next closest value." << endl;
191 cout <<
"tryTwiceLow= " << tryTwiceLow <<
" iNotValidLow= " << iNotValidLow << endl;
196 double diffToMean = 0.;
197 for(
int iInRun=iFirstToUseInRun; iInRun<nToUseInRun; iInRun++){
198 double lDiffToMean = TMath::Abs(runDVEOrderV[iInRun]-runMean);
199 if(lDiffToMean>diffToMean && iInRun!=iNotValidRMS1 && iInRun!=iNotValidRMS2){
200 diffToMean=lDiffToMean;
206 iNotValidRMS1=iFarOrder;
208 cout <<
"Will remove Further most value and try again." << endl;
210 cout <<
"RMS of Drift Velocity East is too large! Will remove Further most value." << endl;
211 cout <<
"RMS " << runRMS <<
" > " << maxRMS << endl;
212 cout <<
"tryThreeRMS= " << tryThreeRMS <<
" iNotValidRMS1= " << iNotValidRMS1 << endl;
215 iNotValidRMS2=iFarOrder;
217 cout <<
"Will remove Further most value and try again." << endl;
219 cout <<
"RMS of Drift Velocity East is too large! Will remove Further most value." << endl;
220 cout <<
"RMS " << runRMS <<
" > " << maxRMS << endl;
221 cout <<
"tryThreeRMS= " << tryThreeRMS <<
" iNotValidRMS2= " << iNotValidRMS2 << endl;
225 cout <<
"Couldn't find a good macro to be kept for this run! Run should be checked!" << endl;
230 int getDriftVelocityDB(
unsigned int funixTime,
unsigned int& uitimedb,
double& ldvedb ){
232 char *fDbName =
"Calibrations_tpc";
233 char *fTableName =
"tpcDriftVelocity";
234 char *fFlavorName =
"laserDV";
235 char *fTimestamp =
"2004-10-31 00:00:00";
243 dbtable = configNode->addDbTable(fTableName);
247 cout <<
" No Table : " << fTableName << endl;
251 if (fFlavorName != 0){
252 dbtable -> setFlavor(fFlavorName);
253 cout <<
"Flavor is set as " << fFlavorName <<
" by StDbTable::setFlavor." << endl;
255 cout <<
"Flavor is NOT assigned. Default value is set as 'ofl'. " << endl;
256 dbtable -> setFlavor(
"ofl");
260 mgr->setRequestTime(funixTime);
262 mgr->setRequestTime(fTimestamp);
265 mgr->fetchDbTable(dbtable);
266 void* cstruct = dbtable->GetTableCpy();
267 Int_t nrows = dbtable->GetNRows();
269 St_tpcDriftVelocity * tpcDV =
270 (St_tpcDriftVelocity*)
TTable::New(fTableName,fTableName,cstruct,nrows);
272 if(!(tpcDV && tpcDV->HasData())) {
273 cout <<
"No TTable returned!" << endl;
276 tpcDriftVelocity_st *ltpcDVs = (tpcDriftVelocity_st*)tpcDV->
At(0);
278 uitimedb = dbtable->getBeginTime();
279 ldvedb = ltpcDVs->laserDriftVelocityEast;
290 void LoadLaserDriftVelocityToDb(
const char* dirName,
const char* listOfMacros,
const char* baseName,
int mode=0){
302 const int nMacrosMax = 800;
305 TString fullListName(dirName);
307 fullListName +=listOfMacros;
309 const char* tmpLines =
"/tmp/nLines.txt";
310 TString myCount(
"/usr/bin/wc -l ");
311 myCount += fullListName;
319 cout <<
"Executing: " << myCount.Data() << endl;
321 gSystem->Exec(myCount.Data());
322 ifstream dvc(tmpLines);
324 TString myRm(
"/bin/rm ");
326 gSystem->Exec(myRm.Data());
328 cout <<
"Number of macros in list: " << nMacros << endl;
329 nMacros = TMath::Min(nMacros,nMacrosMax);
330 cout <<
"Number of macros to read: " << nMacros << endl;
332 const double maxdDVdt = 3*0.000002;
338 gSystem->Load(
"St_base");
339 gSystem->Load(
"StUtilities");
340 gSystem->Load(
"StarClassLibrary");
341 gSystem->Load(
"St_Tables");
342 gSystem->Load(
"StDbLib");
345 double *xLaserDVE =
new double[nMacros];
346 double *yLaserDVE =
new double[nMacros];
348 double *xLaserDVW =
new double[nMacros];
349 double *yLaserDVW =
new double[nMacros];
351 int *macroStatus =
new int[nMacros];
354 TString fullNameMacro;
355 unsigned int uidatetime;
357 St_tpcDriftVelocity *tpcDV;
358 tpcDriftVelocity_st *tpcDVs;
361 theList = fopen(fullListName.Data(),
"r");
367 const int maxDeltaT=120;
368 const int maxFilesInRun=30;
369 const int minFilesInRun=5;
374 unsigned int lastTime=0;
375 int *macroToKeep =
new int[nMacros];
376 int *runMacro =
new int[maxFilesInRun];
377 unsigned int *runTime =
new unsigned int[maxFilesInRun];
378 double *runDVE =
new double[maxFilesInRun];
381 for(
int iMacro=0; iMacro<nMacros; iMacro++){
383 cout <<
"\n" <<
"iMacro = " << iMacro << endl;
384 macroStatus[iMacro]=6;
385 fscanf(theList,
"%s",&nameMacro);
386 if(strstr(nameMacro,baseName)){
387 fullNameMacro = TString(dirName);
389 fullNameMacro+=nameMacro;
390 printf(
"Reading macro: %s\n",fullNameMacro.Data());
392 char* sdatetime=strstr(nameMacro,
"."); sdatetime++;
393 TString tsdatetime(sdatetime);
394 tsdatetime = tsdatetime.ReplaceAll(
".C",
"");
395 tsdatetime = tsdatetime.ReplaceAll(
".",
" ");
396 tsdatetime = tsdatetime.Insert(4,
"-");
397 tsdatetime = tsdatetime.Insert(7,
"-");
398 tsdatetime = tsdatetime.Insert(13,
":");
399 tsdatetime = tsdatetime.Insert(16,
":");
400 sdatetime = tsdatetime.Data();
402 mgr->setRequestTime(sdatetime);
403 uidatetime = mgr->getUnixRequestTime();
404 cout <<
" " << sdatetime <<
" " << uidatetime << endl;
405 if (uidatetime == 0){
406 cout <<
" Macro error " << endl;
411 gROOT->LoadMacro(fullNameMacro.Data());
412 tpcDV = (St_tpcDriftVelocity*)CreateTable();
413 if(tpcDV && tpcDV->HasData()) {
414 cout <<
" Macro succesfully loaded" << endl;
416 tpcDVs = (tpcDriftVelocity_st*)tpcDV->At(0);
418 xLaserDVE[iMacro]=uidatetime;
419 yLaserDVE[iMacro]=tpcDVs->laserDriftVelocityEast;
421 xLaserDVW[iMacro]=uidatetime;
422 yLaserDVW[iMacro]=tpcDVs->laserDriftVelocityWest;
424 cout <<
" Laser Drift Velocity East: " <<tpcDVs->laserDriftVelocityEast << endl;
432 macroToKeep[nMacroToKeep++] = iMacro;
433 macroStatus[iMacro] = iMacro;
439 if (nInRun==0){ lastTime = uidatetime;}
441 if ((uidatetime-lastTime)<maxDeltaT && nInRun<maxFilesInRun ){
443 cout <<
" Delta time passed & not enough count "
444 << nInRun <<
" <=> " << maxFilesInRun <<
" yet" << endl;
445 runMacro[nInRun]=iMacro;
446 runTime[nInRun]=uidatetime;
447 runDVE[nInRun]=tpcDVs->laserDriftVelocityEast;
450 lastTime = uidatetime;
452 }
else if (nInRun>=minFilesInRun-1){
455 int keepIt = getDVInRun(nInRun,runMacro,runTime,runDVE);
457 macroToKeep[nMacroToKeep++] = keepIt;
458 for (
int iInRun=0; iInRun<nInRun; iInRun++){
459 macroStatus[runMacro[iInRun]]=keepIt;
463 for (
int iInRun=0; iInRun<nInRun; iInRun++){
464 macroStatus[runMacro[iInRun]]=2;
470 runMacro[nInRun]=iMacro;
471 runTime[nInRun]=uidatetime;
472 runDVE[nInRun]=tpcDVs->laserDriftVelocityEast;
474 lastTime = uidatetime;
478 cout <<
" Not Enough macros!" << nInRun <<
" " << minFilesInRun << endl;
479 for (
int iInRun=0; iInRun<nInRun; iInRun++){
480 macroStatus[runMacro[iInRun]]=1;
488 runMacro[nInRun]=iMacro;
489 runTime[nInRun]=uidatetime;
490 runDVE[nInRun]=tpcDVs->laserDriftVelocityEast;
492 lastTime = uidatetime;
500 cout <<
"\n" <<
"Processing information" << endl;
502 if (nInRun>=minFilesInRun-1 || mode==1){
504 int keepIt = mode || getDVInRun(nInRun,runMacro,runTime,runDVE);
506 macroToKeep[nMacroToKeep++] = keepIt;
507 for (
int iInRun=0; iInRun<nInRun; iInRun++){
508 macroStatus[runMacro[iInRun]]=keepIt;
512 for (
int iInRun=0; iInRun<nInRun; iInRun++){
513 macroStatus[runMacro[iInRun]]=2;
518 cout <<
"We concluded we do not have Enough macros at this stage" << endl;
519 for (
int iInRun=0; iInRun<nInRun; iInRun++){
520 macroStatus[runMacro[iInRun]]=1;
524 cout <<
"Selected " << nMacroToKeep <<
" out of " << nTryToKeep
525 <<
" runs with enough macros." << endl;
529 TString LoadDone(dirName); LoadDone +=
"/Load/Done/";
530 TString LoadFailed(dirName); LoadFailed +=
"/Load/Failed/";
531 TString LoadOthers(dirName); LoadOthers +=
"/Load/Others/";
532 TString CheckBadRun(dirName); CheckBadRun+=
"/Check/BadRun/";
533 TString CheckVarSel(dirName); CheckVarSel+=
"/Check/VarSel/";
534 TString CheckVarOth(dirName); CheckVarOth+=
"/Check/VarOth/";
541 unsigned int uidatetimedb =0;
542 double laseredvdb =0.;
547 int iMacroToKeep = 0;
550 int *macroKept =
new int[nMacros];
557 for(
int iMacro=0; iMacro<nMacros; iMacro++){
558 cout <<
"iMacro = " << iMacro << endl;
560 fscanf(theList,
"%s",&nameMacro);
561 cout <<
"macroToKeep[" << iMacroToKeep <<
"]=" << macroToKeep[iMacroToKeep] << endl;
562 if (iMacro==macroToKeep[iMacroToKeep]){
563 if(strstr(nameMacro,baseName)){
564 fullNameMacro = TString(dirName);
566 fullNameMacro+=nameMacro;
567 printf(
"Reading macro: %s\n",fullNameMacro.Data());
569 char* sdatetime=strstr(nameMacro,
"."); sdatetime++;
570 TString tsdatetime(sdatetime);
571 tsdatetime = tsdatetime.ReplaceAll(
".C",
"");
572 tsdatetime = tsdatetime.ReplaceAll(
".",
" ");
573 tsdatetime = tsdatetime.Insert(4,
"-");
574 tsdatetime = tsdatetime.Insert(7,
"-");
575 tsdatetime = tsdatetime.Insert(13,
":");
576 tsdatetime = tsdatetime.Insert(16,
":");
577 sdatetime = tsdatetime.Data();
578 mgr->setRequestTime(sdatetime);
579 uidatetime = mgr->getUnixRequestTime();
580 cout <<
"sdatetime=" << sdatetime <<
" uidatetime=" << uidatetime << endl;
584 getDriftVelocityDB(uidatetime,uidatetimedb,laseredvdb);
586 dt = -(double) uidatetimedb;
591 gROOT->LoadMacro(fullNameMacro.Data());
592 tpcDV = (St_tpcDriftVelocity*)CreateTable();
593 if(!(tpcDV && tpcDV->HasData())) {
594 cout <<
"No TTable returned!" << endl;
597 tpcDVs = (tpcDriftVelocity_st*)tpcDV->At(0);
600 cout <<
"dt=" << dt <<
" uidatetime=" << uidatetime <<
" -> " << (double) uidatetime << endl;
601 dt += (double) uidatetime;
602 cout <<
"Now dt=" << dt << endl;
605 dDV += tpcDVs->laserDriftVelocityEast;
609 if(dDV==0) dDVdt = 0;
610 else dDVdt = maxdDVdt+1;
612 dDVdt = TMath::Abs(dDV/dt);
615 cout <<
"dDV/dt=" << dDV <<
"/" << dt <<
"=" << dDVdt << endl;
618 macroKept[nMacroKept++] = iMacro;
620 while(macroStatus[++jMacro]==iMacro){
621 macroStatus[jMacro]=3;
624 while(macroStatus[--jMacro]==iMacro){
625 macroStatus[jMacro]=3;
627 macroStatus[iMacro]=0;
629 dt=-(double)uidatetime;
630 dDV=-tpcDVs->laserDriftVelocityEast;
634 while(macroStatus[++jMacro]==iMacro){
635 macroStatus[jMacro]=4;
638 while(macroStatus[--jMacro]==iMacro){
639 macroStatus[jMacro]=4;
641 macroStatus[iMacro]=5;
643 dt-=(double)uidatetime;
644 dDV-=tpcDVs->laserDriftVelocityEast;
647 gInterpreter->ResetGlobals();
667 StDbModifier *dm =
new StDbModifier();
668 dm->SetDbName(
"Calibrations_tpc");
669 dm->SetFlavor(
"laserDV");
670 dm->SetTableName(
"tpcDriftVelocity");
675 for(
int iMacro=0; iMacro<nMacros; iMacro++){
676 cout <<
"iMacro = " << iMacro << endl;
678 fscanf(theList,
"%s",&nameMacro);
679 cout <<
"macroStatus[" << iMacro <<
"]=" << macroStatus[iMacro] << endl;
680 if(strstr(nameMacro,baseName)){
681 fullNameMacro = TString(dirName);
683 fullNameMacro+=nameMacro;
685 switch (macroStatus[iMacro]) {
687 dm->SetInputFileName(fullNameMacro.Data());
688 char* sdatetime=strstr(nameMacro,
"."); sdatetime++;
689 TString tsdatetime(sdatetime);
690 tsdatetime = tsdatetime.ReplaceAll(
".C",
"");
691 tsdatetime = tsdatetime.ReplaceAll(
".",
" ");
692 tsdatetime = tsdatetime.Insert(4,
"-");
693 tsdatetime = tsdatetime.Insert(7,
"-");
694 tsdatetime = tsdatetime.Insert(13,
":");
695 tsdatetime = tsdatetime.Insert(16,
":");
696 sdatetime = tsdatetime.Data();
697 dm->SetDateTime(sdatetime);
698 cout <<
"Loading " << fullNameMacro.Data() <<
" ... ";
699 if(dm->WriteDataToDB()==1){
700 cout <<
" Succeeded!" << endl;
705 cout <<
" Failed!" << endl;
712 newFile = TString(dirName);
717 newFile+=CheckBadRun;
725 newFile+=CheckVarOth;
729 newFile+=CheckVarSel;
733 newFile = TString(dirName);
739 cout << newFile.Data() << endl;
740 if(rename(fullNameMacro.Data(),newFile.Data()))
741 cout <<
"Error moving file to " << newFile.Data() << endl;
743 if( (FI = fopen(fullNameMacro.Data(),
"r")) != NULL){
744 cout <<
"File " << fullNameMacro.Data() <<
" remained at its original place " << endl;
752 double *xLaserDVEk =
new double[nMacroKept];
753 double *yLaserDVEk =
new double[nMacroKept];
755 double *xLaserDVWk =
new double[nMacroKept];
756 double *yLaserDVWk =
new double[nMacroKept];
758 for(
int iMacro=0; iMacro<nMacroKept; iMacro++){
759 xLaserDVEk[iMacro] = xLaserDVE[macroKept[iMacro]];
760 yLaserDVEk[iMacro] = yLaserDVE[macroKept[iMacro]];
761 xLaserDVWk[iMacro] = xLaserDVW[macroKept[iMacro]];
762 yLaserDVWk[iMacro] = yLaserDVW[macroKept[iMacro]];
765 TGraph *tpcLaserDriftVelocityEast =
new TGraph(nMacros,xLaserDVE,yLaserDVE);
766 TGraph *tpcLaserDriftVelocityWest =
new TGraph(nMacros,xLaserDVW,yLaserDVW);
768 TGraph *tpcLaserDriftVelocityEastK =
new TGraph(nMacroKept,xLaserDVEk,yLaserDVEk);
769 TGraph *tpcLaserDriftVelocityWestK =
new TGraph(nMacroKept,xLaserDVWk,yLaserDVWk);
771 TString Load(dirName); Load+=
"/Load/";
772 Load+=
"tpcLaserDV_Auto.root";
773 TFile *hFile =
new TFile(Load.Data(),
"RECREATE");
775 tpcLaserDriftVelocityEast->Write(
"tpcLaserDriftVelocityEast");
776 tpcLaserDriftVelocityWest->Write(
"tpcLaserDriftVelocityWest");
777 tpcLaserDriftVelocityEastK->Write(
"tpcLaserDriftVelocityEastK");
778 tpcLaserDriftVelocityWestK->Write(
"tpcLaserDriftVelocityWestK");
785 TH1F *oneDHisto =
new TH1F(
"oneDHisto",
"oneDHisto",100,xLaserDVE[0]-100,xLaserDVE[nMacros-1]+100);
787 oneDHisto->SetMaximum(5.6);
788 oneDHisto->SetMinimum(5.4);
789 oneDHisto->SetTitle(
"TPC Laser Drift Velocity");
790 oneDHisto->SetXTitle(
"Time (unix)");
791 oneDHisto->SetYTitle(
"dv (cm/#mus)");
794 tpcLaserDriftVelocityWest->SetMarkerStyle(23);
795 tpcLaserDriftVelocityWest->SetMarkerColor(2);
798 tpcLaserDriftVelocityEast->SetMarkerStyle(22);
799 tpcLaserDriftVelocityEast->SetMarkerColor(4);
800 tpcLaserDriftVelocityEast->Draw(
"PL");
802 tpcLaserDriftVelocityEastK->SetMarkerStyle(22);
803 tpcLaserDriftVelocityEastK->Draw(
"PL");
807 Load.ReplaceAll(
".root",
".eps");
808 cvn1.SaveAs(Load.Data());
810 cout << endl << endl <<
" ******************************************** " << endl;
811 cout <<
" Number of Macros attempted = " << nMacros << endl;
812 cout << nSuccess <<
" loaded to DB and " << nFailed <<
" failed to be load!" << endl;
static TTable * New(const Char_t *name, const Char_t *type, void *array, UInt_t size)
This static method creates a new TTable object if provided.
static StDbManager * Instance()
strdup(..) is not ANSI
const void * At(Int_t i) const
Returns a pointer to the i-th row of the table.