POTSpillRate.py
Go to the documentation of this file.
1 #!/usr/bin/env python
2 """
3  Compute average delivered POT during a given time range
4  Justin Vasel (jvasel@fnal.gov)
5 
6  """
7 
8 #.......................................................................
9 # INCLUDES
10 #
11 #.......................................................................
12 # Basics
13 import sys
14 import math
15 import urllib2 # Retrieves contents of file at URL
16 import csv # CSV file i/o
17 import getopt
18 
19 # Date and time stuff
20 from datetime import date, datetime, timedelta
21 
22 # Functions from `POTGetData.py'
23 from POTGetData import database_connect
24 from POTGetData import get_runs
25 from POTGetData import get_active_chans
26 #from POTGetData import get_spills
27 from POTGetData import unix_timestamp
28 
29 # PyRoot
30 from ROOT import TFile, TTree
31 from ROOT import gROOT, AddressOf
32 
33 #.......................................................................
34 # DEFINITIONS
35 #
36 #.......................................................................
37 
38 # Get system arguments
39 view = sys.argv[1]
40 beam = str("numi")
41 if len(sys.argv) > 2:
42  beam = sys.argv[2]
43 
44 # Initialize stuff
45 ndays = -1
46 interval = -1
47 
48 BASE_DIR = "/nusoft/app/web/htdoc/nova/datacheck/nearline/data/"
49 
50 # PyRoot stuff
51 gROOT.ProcessLine("struct BeamStruct {Double_t Time; Double_t Beam;};")
52 from ROOT import BeamStruct
53 beamStruct = BeamStruct()
54 
55 # Get arguments and set options
56 if beam == str("numi"):
57  if view == "Day" or view == "day":
58  f = TFile(BASE_DIR + 'POTSpillRateDay.root','RECREATE')
59  interval = 10 # in minutes
60  ndays = 1 # in days
61  elif view == "Week" or view == "week":
62  f = TFile(BASE_DIR + 'POTSpillRateWeek.root','RECREATE')
63  interval = 60 # in minutes
64  ndays = 7 # in days
65  elif view == "Month" or view == "month":
66  f = TFile(BASE_DIR + 'POTSpillRateMonth.root','RECREATE')
67  interval = 3 * 60 # in minutes
68  ndays = 31 # in days
69  elif view == "6Months" or view == "6months":
70  f = TFile(BASE_DIR + 'POTSpillRate6Months.root','RECREATE')
71  interval = 24 * 60 # in minutes
72  ndays = 182 # in days
73  elif view == "Year" or view == "year":
74  f = TFile(BASE_DIR + 'POTSpillRateYear.root','RECREATE')
75  interval = 24 * 60 # in minutes
76  ndays = 365 # in days
77  else:
78  print "ERROR: You did not specify a valid time range. Possibilities are 'Day', 'Week', 'Month', '6Months', or 'Year'."
79  sys.exit()
80 
81 elif beam == str("bnb"):
82  if view == "Day" or view == "day":
83  f = TFile(BASE_DIR + 'BNBPOTSpillRateDay.root','RECREATE')
84  interval = 10 # in minutes
85  ndays = 1 # in days
86  elif view == "Week" or view == "week":
87  f = TFile(BASE_DIR + 'BNBPOTSpillRateWeek.root','RECREATE')
88  interval = 60 # in minutes
89  ndays = 7 # in days
90  elif view == "Month" or view == "month":
91  f = TFile(BASE_DIR + 'BNBPOTSpillRateMonth.root','RECREATE')
92  interval = 3 * 60 # in minutes
93  ndays = 31 # in days
94  elif view == "6Months" or view == "6months":
95  f = TFile(BASE_DIR + 'BNBPOTSpillRate6Months.root','RECREATE')
96  interval = 24 * 60 # in minutes
97  ndays = 182 # in days
98  elif view == "Year" or view == "year":
99  f = TFile(BASE_DIR + 'BNBPOTSpillRateYear.root','RECREATE')
100  interval = 24 * 60 # in minutes
101  ndays = 365 # in days
102  else:
103  print "ERROR: You did not specify a valid time range. Possibilities are 'Day', 'Week', 'Month', '6Months', or 'Year'."
104  sys.exit()
105 
106 # More PyRoot stuff
107 tree = TTree('BeamMetrics','POTtree')
108 tree.Branch('Time', AddressOf(beamStruct, 'Time'), 'Time/D')
109 tree.Branch('SpillRate', AddressOf(beamStruct, 'Beam'), 'Beam/D')
110 
111 
112 #.......................................................................
113 def get_spills(day, deltaDay, beam):
114  """ Get info about spills during a particular day
115 
116  Args:
117  day : A datetime object
118  deltaDay : Number of days to consider
119 
120  Returns:
121  spills : A list of spills (time in sec, pot number in 1E+12)
122  """
123 
124  # Initialize list of spills
125  spills = []
126 
127  # Convert start and end of day to UNIX timestamp. We must explicitly specify
128  # the UTC timezone (+00:00) because the ifb data server uses CST (-06:00) as
129  # the default
130  t0 = str(day.isoformat('T')) + "+00:00"
131  t1 = str((day + timedelta(days = deltaDay)).isoformat('T')) + "+00:00"
132 
133  print "Getting spills for the period between {} and {}".format(t0,t1)
134 
135  # URL where the spill data live
136  if beam == str("numi"):
137  url = 'http://ifb-data.fnal.gov:8099/ifbeam/data/data?'\
138  'v=E:TRTGTD&e=e,a9&'\
139  't0={}&t1={}&&f=csv'.format(t0, t1)
140  elif beam == str("bnb"):
141  url = 'http://ifb-data.fnal.gov:8099/ifbeam/data/data?'\
142  'v=E:TOR860&e=e,1d&'\
143  't0={}&t1={}&&f=csv'.format(t0, t1)
144 
145  print "URL: {}".format(url)
146 
147  # Store the HTTP GET response here
148  response = None
149 
150  print "Trying for a response...",
151 
152  # Try to retrieve the URL contents. If this fails three times, then
153  # give up and return an empty array.
154  try:
155  response = urllib2.urlopen(url, None, 200)
156  except:
157  print "No response. Trying again..."
158  try:
159  response = urllib2.urlopen(url, None, 200)
160  except:
161  print "No response. Trying again..."
162  try:
163  response = urllib2.urlopen(url, None, 200)
164  except:
165  print "No response. Will not try again."
166  sys.exit()
167 
168  if (response == None):
169  return spills
170 
171  # Use csv module to parse the response as CSV
172  csvdata = csv.reader(response)
173  print "SUCCESS"
174 
175  # The response will look like this:
176  # "e,a9", E:TRTGTD, 1406851201047, E12, 22.0598346211
177  # We are interested in the third column (unix timestamp in ms)
178  # and the fifth column (value). The list will be populated with these
179  # values.
180  print "Importing Spills..."
181  if beam == str("numi"):
182  for row in csvdata:
183  if (row.count('E:TRTGTD') > 0):
184  if (row[4] != 'null'):
185  spills.append( (float(row[2]) / 1000., float(row[4])) )
186 
187  elif beam == str("bnb"):
188  for row in csvdata:
189  if (row.count('E:TOR860') > 0):
190  if (row[4] != 'null'):
191  spills.append( (float(row[2]) / 1000., float(row[4])) )
192  print "DONE"
193 
194  # Return spill info
195  return spills
196 
197 #.......................................................................
198 def num_spills_between(t0, t1, spills):
199  """ Determines the number of POT in a given time interval
200 
201  Args:
202  t0 : A datetime object
203  t1 : A datetime object
204  spills : A list of spills
205 
206  Returns:
207  potTotal : int
208  """
209 
210  # Initiate connection to psql database
211  database_connect(beam)
212 
213  print "Couting spills between {} and {}: ".format(str(t0.isoformat('T')),str(t1.isoformat('T'))),
214 
215  # Convert to unix timestamps
216  t0 = unix_timestamp(t0)
217  t1 = unix_timestamp(t1)
218 
219  # initialize counter
220  counter = 0
221 
222  # Loop over spills, but only count the ones that fall into the
223  # appropriate time window.
224  for (time, value) in spills:
225  if ((time >= t0) & (time <= t1)):
226  counter = counter + 1
227 
228  return counter
229 
230 #.......................................................................
231 def main():
232 
233  # Set startTime and currentTime to now
234  # (startTime is a constant; currenTime is variable).
235  startTime = datetime.utcnow()
236  currentTime = startTime
237 
238  # Loop over days
239  for i in range(1, int(ndays) + 1):
240  print "DAY {} of {}".format(i, ndays)
241  print "====================="
242 
243  # Get two days worth of spills (we'll only use a subset of this)
244  spills = get_spills(startTime - timedelta(days = i), 2, beam)
245 
246  # Get number of time steps
247  nsteps = int(math.floor(ndays * 24 * 60 / interval))
248 
249  print "Begin looping over {} intervals of {} minutes each".format(nsteps, interval)
250 
251  # Loop over steps in day
252  for j in range(0, int(nsteps / ndays)):
253  # Find total number of spills in this timestep
254  spillTotal = num_spills_between(currentTime - timedelta(minutes = interval), currentTime, spills)
255 
256  # Calculate spill rate (spills per second)
257  spillRate = spillTotal * 1.0 / (interval * 60)
258 
259  # Pring spill rate to screen
260  print spillRate
261 
262  # Fill ROOT branches/tree
263  beamStruct.Time = int((currentTime - datetime(1970,1,1)).total_seconds())
264  beamStruct.Beam = float(spillRate)
265  tree.Fill()
266 
267  # Advance to next timestep
268  currentTime = currentTime - timedelta(minutes = interval)
269 
270  print "\n"
271 
272  f.Write()
273  f.Close()
274 
275 
276 
277 # Prevent execution on import
278 if __name__ == "__main__":
279  main()
280 ########################################################################
def get_spills(day, deltaDay, beam)
std::string format(const int32_t &value, const int &ndigits=8)
Definition: HexUtils.cpp:14
def unix_timestamp(input)
Definition: POTGetData.py:46
def database_connect()
def num_spills_between(t0, t1, spills)