API tools faq
paste
Guest User

Untitled

a guest
Jul 23rd, 2018
161
0
Never
Add comment
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 13.53 KB | None | 0 0
raw download clone embed print report
  1. /*
  2. plotter.c
  3. Giulio Candreva (giulio.candreva@cern.ch)
  4. ------------------------------
  5. the purpose is to plot the data from the summary files, generated by the script WIP.c
  6. the script is divided in various functions:
  7.  
  8. - split(string s, char delim, vector <string> *elements)
  9. Takes a string s, splits it in chunks using the 'delim' char as divider.
  10. The chunks are put in the *elements vector
  11.  
  12. - strToTime(string date)
  13. Takes a string in the format dd/mm/yyyy and returns the UNIX timestamp
  14.  
  15. - plot(vector <string> files, string column, string start, float min, float max)
  16. Reads all the files contained in the vector "file", extract the column "column" and the column "Date"
  17. it plots all the points in a plot that has the time in the X axis
  18. the X axis is limited by start and stop (they are two strings containing two dates in the format dd/mm/yyyy) stop is the last analysis date
  19. the Y axis is limited by min and max
  20.  
  21. - plotter()
  22. The main function of the script. It calls the function plot with various parameteres in orderd to generate
  23. the desired plots
  24.  
  25. - drawLogo()
  26. helpful function that draws a overlaying logo in the canvas
  27.  
  28. !!IF YOU OPENED THIS FILE YOU PROBABLY WANT TO EDIT "plot" AND "plotter" FUNCTIONS SINCE THEY'RE THE CORE OF THE SCRIPT!!
  29.  
  30. */
  31.  
  32. #include <string>
  33. #include <stdlib>
  34. #include <fstream>
  35. #include <iostream>
  36. #include <TGaxis.h>
  37. #include "TCanvas.h"
  38. #include <TMath.h>
  39. #include "TColor.h"
  40.  
  41. /* User defined function that splits a string s using a delimitator delim and puts the chuncks in the elements vector */
  42. /* You probably don't need to know how it works */
  43. void split(string s, char delim, vector <string> * elements) {
  44. string item = "";
  45. int i=0;
  46. while(i<s.size()) {
  47. if (s.at(i) == delim) {
  48. (*elements).push_back(item);
  49. item = "";
  50. } else {
  51. item += s.at(i);
  52. }
  53. i++;
  54. }
  55. (*elements).push_back(item);
  56. }
  57.  
  58.  
  59. /* user defined function that converts a date in the format dd/mm/yyyy to a UNIX timestamp */
  60. int strToTime(string datestring) {
  61.  
  62. int day, year, month;
  63. vector <string> elements;
  64. split(datestring, '/', &elements);
  65.  
  66. day = atoi(elements[0].c_str());
  67. month = atoi(elements[1].c_str());
  68. year = atoi(elements[2].c_str());
  69.  
  70. TDatime * date = new TDatime(year, month, day, 0, 0, 0);
  71. UInt_t ttt;
  72. ttt = date->Convert();
  73. int pTime = (int)ttt;
  74. return pTime;
  75. }
  76.  
  77.  
  78. // User defined function that draws lines for acceptable limits
  79. void drawLimits(string det, string column, float startday, float stopday) {
  80.  
  81. float suplim=0, inflim, values[2], times_l[2];
  82. times_l[0] = startday;
  83. times_l[1] = stopday;
  84.  
  85. if (det == "RICH2") {
  86. if (column == "CF4") {
  87. suplim = 93;
  88. inflim = 91;
  89. } else if (column == "CO2") {
  90. suplim = 9;
  91. inflim = 7;
  92. } else if (column == "Oxygen") {
  93. suplim = 3000;
  94. inflim = -1;
  95. } else if (column == "Nitrogen") {
  96. suplim = 12000;
  97. inflim = -1;
  98. } else if (column == "Oxygen_2018") {
  99. suplim = 3000;
  100. inflim = -1;
  101. } else if (column == "Nitrogen_2018") {
  102. suplim = 12000;
  103. inflim = -1;
  104. }
  105. } else if (det == "MWPC") {
  106. if (column == "CO2") {
  107. suplim = 58;
  108. inflim = 56;
  109. } else if (column == "Ar") {
  110. suplim = 39;
  111. inflim = 37;
  112. } else if (column == "CF4") {
  113. suplim = 5.5;
  114. inflim = 4.5;
  115. }
  116. } else if (det == "RICH1") {
  117. if (column == "C4F10") {
  118. suplim = 100;
  119. inflim = 100;
  120. } else if (column == "CO2") {
  121. suplim = 0.5;
  122. inflim = -1.0;
  123. } else if (column == "Air") {
  124. suplim = 2;
  125. inflim = -1;
  126. } else if (column == "Oxygen") {
  127. suplim = 2000;
  128. inflim = -1;
  129. } else if (column == "Nitrogen") {
  130. suplim = 8000;
  131. inflim = -1;
  132. }
  133. }
  134.  
  135. //MARA 30/04 modified for OT plots
  136. //ANNA change limits for lines on the plot
  137. else if (det == "OT") {
  138. if (column == "CO2") {
  139. suplim = 0.5;
  140. inflim = -1.0;
  141. } else if (column == "Argon") {
  142. suplim = 2;
  143. inflim = -1;
  144. } else if (column == "Oxygen") {
  145. suplim = 2000;
  146. inflim = -1;
  147. } else if (column == "Nitrogen") {
  148. suplim = 8000;
  149. inflim = -1;
  150. }
  151. }
  152. cout << "limits " << det << " " << column << " " << suplim << " " << inflim << endl;
  153. if (suplim != 0) {
  154. // Add line 1 (inferior limit)
  155. values[0] = inflim;
  156. values[1] = inflim;
  157. TGraph *Line = new TGraph(2, times_l, values);
  158. Line->SetMarkerColor(4); //(kMagenta-2);
  159. Line->SetMarkerStyle(8);
  160. Line->SetLineColor(kRed);
  161. Line->SetLineWidth(2);
  162. Line->SetMarkerSize(0.2);
  163. Line->Draw("CP");
  164. // Add line 2 (superior limit)
  165. values[0] = suplim;
  166. values[1] = suplim;
  167. TGraph *LineS = new TGraph(2, times_l, values);
  168. LineS->SetMarkerColor(kAzure+1); //(kMagenta-2);
  169. LineS->SetMarkerStyle(8);
  170. LineS->SetLineColor(kRed);
  171. LineS->SetLineWidth(2);
  172. LineS->SetMarkerSize(0.2);
  173. LineS->Draw("CP");
  174. }
  175.  
  176. }
  177.  
  178.  
  179. // User defined function that draws a logo based on an image
  180. void drawLogo() {
  181. TImage *img = TImage::Open("epdtlogo.jpg");
  182.  
  183. if (!img) {
  184. printf("Could not create an image... exit\n");
  185. return;
  186. }
  187. // TPad *l = new TPad("l", "l",0.1209553,0.7906977,0.220339,0.8901809);
  188. //gPad->cd(0);
  189. TPad *l = new TPad("l", "l",0.1209553,0.7906977,0.3466872,0.8901809);
  190. l->Draw();
  191. l->cd();
  192. img->Draw();
  193. }
  194.  
  195. /* function that actually plots */
  196. void plot(vector <string> names, string column, string start, float min, float max, bool proportional=false){
  197.  
  198. //------- GENERAL STYLE INSTRUCTIONS
  199. gROOT->SetStyle("Plain");
  200. // background is no longer mouse-dropping white
  201. gStyle->SetCanvasColor(kWhite);
  202. // blue to red false color palette. Use 9 for b/w
  203. gStyle->SetPalette(1,0);
  204. // turn off canvas borders
  205. gStyle->SetCanvasBorderMode(0);
  206. gStyle->SetPadBorderMode(0);
  207. // What precision to put numbers if plotted with "TEXT"
  208. gStyle->SetPaintTextFormat("5.2f");
  209.  
  210. // For publishing:
  211. gStyle->SetLineWidth(1.5);
  212. gStyle->SetTextSize(1.1);
  213. gStyle->SetLabelSize(0.03,"xy");
  214. gStyle->SetTitleSize(0.04,"xy");
  215. gStyle->SetTitleOffset(1.1,"x");
  216. gStyle->SetTitleOffset(0.9,"y");
  217. gStyle->SetPadTopMargin(0.1);
  218. gStyle->SetPadRightMargin(0.05);
  219. gStyle->SetPadBottomMargin(0.1);
  220. gStyle->SetPadLeftMargin(0.1);
  221.  
  222. //------ END GENERAL STYLE INSTRUCTIONS
  223.  
  224. // Initialize the Canvas, graphs and legend
  225. TCanvas *c;
  226. char namecanvas[100]; //the TCanvas constructor accepts only char*
  227. sprintf(namecanvas, "%s_%d_%d", column,min,max);
  228. c = new TCanvas(namecanvas,"plot",0,0,1300,800); // It doesn't accept string, so I did the conversion
  229. c->SetFillColor(0);
  230. c->GetFrame()->SetBorderSize(0);
  231. c->SetGridy();
  232. TGraph *Graph[3]; // I will need at most 3 graphs in the same canvas
  233. //leg = new TLegend(0.6687211,0.7622739,0.9306626,0.872093,NULL,"brNDC");
  234. leg = new TLegend(0.801877,0.7769784,0.9442127,0.8893885,NULL,"brNDC");
  235. leg->SetFillColor(0);
  236. leg->SetBorderSize(1);
  237. leg->SetTextSize(0.03);
  238. // end canvas and graphs
  239.  
  240. string name; // this will contain the name of the processed file
  241. int np = 0; // a support variable that counts the plots that I draw
  242. char *title; // this will be the title of the graph
  243. string titolo; // I'll use these two variables later
  244.  
  245. for (np=0;np < names.size();np++) { // I draw many plots as many as files
  246.  
  247. // open file
  248. name = names[np];
  249. ifstream thaFile(name.c_str(),ios::in);
  250.  
  251. int pos=0; // It will contain the number of the column
  252. int i=0; // line counter
  253. string dates[10000]; // this will contain all the dates
  254. float areas[10000]; // this will contain all the concentrations (although it's called 'areas')
  255. float errors_y[10000]; // this will contain the errors
  256. float errors_x[10000]; // this will contain the errors on the x axis (it will be filled with zeros)
  257. float times[100000]; // this will contain all the unix timestamps (converted from the dates)
  258. float areas_sum; // this will contain the sum over all areas (possibly used as calib. value, depending on the "proportional" parameter)
  259.  
  260. string lineString; // this will contain the current line of the file
  261. vector <string> elements; // every line will be split into chunks and put in elements
  262.  
  263. cout << name << endl; // I print the current file that is being processed
  264.  
  265. while(getline(thaFile, lineString)) { // for each line
  266. elements.clear();
  267. split(lineString, '\t', &elements); // I split the line as said before
  268.  
  269. if (elements[0] == "Date") { // <=> If it's the first line
  270.  
  271. for (int j=0;j<elements.size();j++) { // I cycle over all the columns name
  272. if (column == elements[j]) { // till I find my column, specified in the function arguments
  273. (proportional==true) ? pos = j : pos = j+2; // And I save its position (+2 cause that's the column that contains the concentration, +0 if I want the raw area)
  274. }
  275. }
  276. continue;
  277. }
  278.  
  279. // Now I want to collect the sum of all the areas for the given row, because it will be used as calibration value if the parameter "proportional" is set to true
  280. areas_sum = 1;
  281. for (int j=6;j<elements.size();j++) { // I start from 6 (that becomes 9 due to the mod4 ) because the first 8 values are relative to MS column
  282. if ((j-1)%4 == 0) {
  283. areas_sum = areas_sum + atof(elements.at(j).c_str());
  284. }
  285. }
  286.  
  287. times[i] = (float) strToTime(elements[0]); // I convert the date (elements(0)) into a timestamp
  288. areas[i] = atof(elements.at(pos).c_str()); // I convert the area (elements(pos)) into a float number
  289. areas[i] = (areas[i] > 0) ? areas[i] : -1; // If the concentration is zero, then it means that the data is missing: I avoid to plot it, putting it below the x-axis (I assign it a negative value)
  290. errors_y[i] = atof(elements.at(pos+1).c_str()); // the uncertainty is the next column
  291. errors_x[i] = 0; // I don't want any error for my x-axis
  292.  
  293. if (proportional) { // If proportional is set to true, I divide the area by the sum of all the others
  294. areas[i] = 100.0 * areas[i] / areas_sum;
  295. errors_y[i] = 100.0 * errors_y[i] / areas_sum;
  296. }
  297. cout << i << " : " << (int) times[i] << " -> " << areas[i] << " +- " << errors_y[i] << "( areas sum = " << areas_sum << ")" << endl; // just for debug
  298. i++;
  299. } // end reading the file
  300.  
  301. // Now I use the name of the file in order to construct the title of the plot
  302. vector <string> parts;
  303. split(name, '_', &parts);
  304. titolo = parts.at(1) + " " + column;
  305. title = titolo;
  306. char *line;
  307. line = parts.at(2); // the second part of the filename should contain the gas line
  308. string detector = parts.at(1);
  309. // End title management
  310.  
  311.  
  312. // I create the Graph and define a lot of graphic things
  313. Graph[np] = new TGraphErrors(i, times, areas, errors_x, errors_y);
  314. Graph[np]->SetTitle(title);
  315. Graph[np]->GetXaxis()->SetTitle("Time");
  316. Graph[np]->GetYaxis()->SetTitle((column == "Oxygen" || column == "Nitrogen") ? "Concentration (ppm) " : "Concentration (%)");
  317. Graph[np]->GetXaxis()->SetNdivisions(608,kFALSE); //605
  318. Graph[np]->GetXaxis()->SetTimeDisplay(1);
  319. Graph[np]->GetXaxis()->SetTimeFormat("%d\/%m\/%y %F1970-01-01 00:00:00");
  320. Graph[np]->SetMarkerStyle(8);
  321. Graph[np]->SetMarkerSize(1);
  322.  
  323. // Down here I set the boundaries of the axis
  324. int startday = strToTime(start);
  325. //int stopday = strToTime(stop);
  326. int stopday = (int) times[i-1];
  327. Graph[np]->GetXaxis()->SetLimits(startday,stopday);
  328. Graph[np]->SetMinimum(min);
  329. Graph[np]->SetMaximum(max);
  330.  
  331. // I add an entry in the legend, corresponding to this graph (line is the name of the gas line)
  332. leg->AddEntry(Graph[np],line,"p");
  333.  
  334. // Since I'm drawing more graphs on the top of each other, I change the colours of the points
  335. TColor *colors[3] = {gROOT->GetColor(kAzure), gROOT->GetColor(kGreen+2), gROOT->GetColor(kRed-2)};
  336. if (np == 0) {
  337. Graph[np]->SetMarkerColor(colors[0]->GetNumber());
  338. Graph[np]->Draw("AP");
  339. } else {
  340. Graph[np]->SetMarkerColor(colors[np]->GetNumber());
  341. Graph[np]->Draw("P");
  342. }
  343.  
  344.  
  345. } // End for each file processed
  346. drawLimits(detector, column, startday, stopday);
  347. leg->Draw();
  348. drawLogo();
  349. // let's save the canvas in a png file
  350. char outfile[100];
  351. sprintf(outfile, "%s_%d.png", titolo, proportional);
  352. c->Print(outfile);
  353. } // end plot() function
  354.  
  355. void plotter() {
  356. vector <string> files;
  357.  
  358. // MWPC
  359. files.push_back("Summaries\\150_MWPC_Mixer_Sum.DIF");
  360. files.push_back("Summaries\\151_MWPC_ExhToDis_Sum.DIF");
  361. files.push_back("Summaries\\148_MWPC_BefPurif_Sum.DIF");
  362. plot(files, "CO2", "1/4/2017", 45, 65);
  363. plot(files, "Ar", "1/4/2017", 30, 44);
  364. plot(files, "CF4", "1/4/2017", 4, 7);
  365. plot(files,"Oxygen", "1/4/2017", 0, 200);
  366. plot(files,"Nitrogen", "1/4/2017", 0, 3000);
  367. files.clear();
  368.  
  369. // OT
  370. files.push_back("Summaries\\155_OT_Exh_Sum.DIF");
  371. //files.push_back("Summaries\\154_OT_Mixer_Sum.DIF");
  372. //plot(nameoffilewithdata, "nameofvariable", "startday", minY, maxY)
  373. plot(files, "CO2", "1/4/2018", 25, 35);
  374. plot(files, "Ar", "1/4/2018", 60, 80);
  375. plot(files,"Oxygen", "1/4/2018", 0, 20000);
  376. plot(files,"Nitrogen", "1/4/2018", 0, 500);
  377. files.clear(); //frees the memory of program to go on with other files
  378.  
  379.  
  380. // RICH2
  381. files.push_back("Summaries\\159_RICH2_BefPurif_Sum.DIF");
  382. files.push_back("Summaries\\160_RICH2_AftPurif_Sum.DIF");
  383. plot(files, "CO2", "1/1/2016", 7, 11);
  384. plot(files, "CF4", "1/1/2016", 85, 95);
  385. // I want to plot Ratios for CF4 e CO2
  386. plot(files, "CF4", "1/1/2016", 90, 100, true);
  387. plot(files, "CO2", "1/1/2016", 7, 11, true);
  388.  
  389. //uncomment to get the complete plot from 2016 to 2018
  390. //plot(files, "Oxygen", "1/1/2016", 0, 6000);
  391. //plot(files, "Nitrogen", "1/1/2016", 0, 20000);
  392. plot(files, "Oxygen", "1/3/2018", 4700, 5200); //Anna
  393. plot(files, "Nitrogen", "1/3/2018", 16600, 18300); //Anna
  394. files.clear();
  395.  
  396. // RICH1
  397. files.push_back("Summaries\\146_RICH1_AftPump_Sum.DIF");
  398. plot(files, "C4F10", "1/1/2016", 90, 100, true);
  399. plot(files, "CO2", "1/1/2016", 0, 10, true);
  400. plot(files, "Oxygen", "1/1/2016", 0, 12000);
  401. plot(files, "Nitrogen", "1/1/2016", 0, 25000);
  402. //plot(files, "Oxygen", "1/3/2018", 0, 12000);
  403. //plot(files, "Nitrogen", "1/3/2018", 0, 25000);
  404. plot(files, "Air", "1/1/2016", 0, 1.5, true);
  405. files.clear();
  406.  
  407. }
Add Comment
Please, Sign In to add comment
Public Pastes
We use cookies for various purposes including analytics. By continuing to use Pastebin, you agree to our use of cookies as described in the Cookies Policy.  OK, I Understand
Not a member of Pastebin yet?
Sign Up, it unlocks many cool features!