glue.lal

1 # Copyright (C) 2006-2013 Kipp Cannon 2 # 3 # This program is free software; you can redistribute it and/or modify it 4 # under the terms of the GNU General Public License as published by the 5 # Free Software Foundation; either version 3 of the License, or (at your 6 # option) any later version. 7 # 8 # This program is distributed in the hope that it will be useful, but 9 # WITHOUT ANY WARRANTY; without even the implied warranty of 10 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU General 11 # Public License for more details. 12 # 13 # You should have received a copy of the GNU General Public License along 14 # with this program; if not, write to the Free Software Foundation, Inc., 15 # 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA. 16 17 # 18 # ============================================================================= 19 # 20 # Preamble 21 # 22 # ============================================================================= 23 # 24 25 26 """ 27 This module contains bits and pieces of use when interacting with LAL and 28 LAL-derived code (eg. LALApps programs) 29 """ 30 31 32 import fnmatch 33 import math 34 import os 35 import re 36 import sys 37 from six.moves import urllib 38 import warnings 39 import six 40 41 try: # python < 3 42 long 43 except NameError: # python >= 3 44 long = int 45 46 from glue import git_version 47 from glue import segments 48 49 50 __author__ = "Kipp Cannon <kipp.cannon@ligo.org>" 51 __version__ = "git id %s" % git_version.id 52 __date__ = git_version.date 53 54 55 # 56 # ============================================================================= 57 # 58 # High precision time object 59 # 60 # ============================================================================= 61 # 62 63 64 # 65 # Python version in case LAL isn't available 66 # 67 68

69 -class LIGOTimeGPS(object):

70 """ 71 An object for storing times with nanosecond resolution. LAL defines an 72 equivalent object which is used through-out the search algorithms to 73 represent times. Many LALApps routines input and output times in a 74 manner that meshes well with this object. 75 76 Internally the time is represented as a signed integer "seconds" part 77 and an unsigned integer "nanoseconds" part. The actual time is always 78 constructed by adding the nanoseconds to the seconds. So -0.5 s is 79 represented by setting seconds = -1, and nanoseconds to 500000000. 80 That's the way LAL does it. 81 """ 82 83 # basic class methods 84

85 - def __init__(self, seconds, nanoseconds = 0):

86 """ 87 Create a LIGOTimeGPS instance. The first parameter is the 88 count of seconds, and the second (optional) parameter is the 89 count of nanoseconds. If the nanoseconds parameter is not 90 supplied, it is assumed to be 0. Either parameter can be 91 a numeric type or an ASCII string. 92 93 Example: 94 95 >>> LIGOTimeGPS(100.5) 96 LIGOTimeGPS(100, 500000000) 97 >>> LIGOTimeGPS("100.5") 98 LIGOTimeGPS(100, 500000000) 99 >>> LIGOTimeGPS(100, 500000000) 100 LIGOTimeGPS(100, 500000000) 101 >>> LIGOTimeGPS(0, 100500000000) 102 LIGOTimeGPS(100, 500000000) 103 >>> LIGOTimeGPS(100.2, 300000000) 104 LIGOTimeGPS(100, 500000000) 105 >>> LIGOTimeGPS("0.000000001") 106 LIGOTimeGPS(0, 1) 107 >>> LIGOTimeGPS("0.0000000012") 108 LIGOTimeGPS(0, 1) 109 >>> LIGOTimeGPS("0.0000000018") 110 LIGOTimeGPS(0, 2) 111 >>> LIGOTimeGPS("-0.8") 112 LIGOTimeGPS(-1, 200000000) 113 >>> LIGOTimeGPS("-1.2") 114 LIGOTimeGPS(-2, 800000000) 115 """ 116 if not isinstance(nanoseconds, (float, int, long)): 117 try: 118 nanoseconds = float(nanoseconds) 119 except: 120 raise TypeError(nanoseconds) 121 if isinstance(seconds, float): 122 ns, seconds = math.modf(seconds) 123 seconds = int(seconds) 124 nanoseconds += ns * 1e9 125 elif not isinstance(seconds, six.integer_types): 126 if isinstance(seconds, (six.binary_type, six.text_type)): 127 sign = -1 if seconds.lstrip().startswith("-") else +1 128 try: 129 if "." in seconds: 130 seconds, ns = seconds.split(".") 131 ns = round(sign * float("." + ns) * 1e9) 132 else: 133 ns = 0 134 seconds = int(seconds) 135 except: 136 raise TypeError("invalid literal for LIGOTimeGPS(): %s" % seconds) 137 nanoseconds += ns 138 elif hasattr(seconds, "gpsSeconds") and hasattr(seconds, "gpsNanoSeconds"): 139 # handle LIGOTimeGPS(x) where x is an 140 # object with gpsSeconds and gpsNanoSeconds 141 # fields. 142 nanoseconds += seconds.gpsNanoSeconds 143 seconds = seconds.gpsSeconds 144 elif hasattr(seconds, "seconds") and hasattr(seconds, "nanoseconds"): 145 # handle LIGOTimeGPS(x) where x is an 146 # object with seconds and nanoseconds 147 # fields. 148 nanoseconds += seconds.nanoseconds 149 seconds = seconds.seconds 150 else: 151 raise TypeError(seconds) 152 self.__seconds = seconds + int(nanoseconds // 1000000000) 153 self.__nanoseconds = int(nanoseconds % 1000000000)

154 155 seconds = gpsSeconds = property(lambda self: self.__seconds) 156 nanoseconds = gpsNanoSeconds = property(lambda self: self.__nanoseconds) 157

158 - def __repr__(self):

159 return "LIGOTimeGPS(%d, %u)" % (self.__seconds, self.__nanoseconds)

160

161 - def __str__(self):

162 """ 163 Return an ASCII string representation of a LIGOTimeGPS. 164 """ 165 if (self.__seconds >= 0) or (self.__nanoseconds == 0): 166 s = "%d.%09u" % (self.__seconds, self.__nanoseconds) 167 elif self.__seconds < -1: 168 s = "%d.%09u" % (self.__seconds + 1, 1000000000 - self.__nanoseconds) 169 else: 170 s = "-0.%09u" % (1000000000 - self.__nanoseconds) 171 return s.rstrip("0").rstrip(".")

172 173 # type conversion 174

175 - def __float__(self):

176 """ 177 Convert a LIGOTimeGPS to seconds as a float. 178 179 Example: 180 181 >>> float(LIGOTimeGPS(100.5)) 182 100.5 183 """ 184 return self.__seconds + self.__nanoseconds * 1e-9

185

186 - def __int__(self):

187 """ 188 Return the integer part (seconds) of a LIGOTimeGPS as an int. 189 190 Example: 191 192 >>> int(LIGOTimeGPS(100.5)) 193 100 194 """ 195 return self.__seconds

196

197 - def __long__(self):

198 """ 199 Return the integer part (seconds) of a LIGOTimeGPS as a long. 200 201 Example: 202 203 >>> long(LIGOTimeGPS(100.5)) 204 100L 205 """ 206 return long(self.__seconds)

207

208 - def ns(self):

209 """ 210 Convert a LIGOTimeGPS to a count of nanoseconds as a long. 211 212 Example: 213 214 >>> LIGOTimeGPS(100.5).ns() 215 100500000000 216 """ 217 return self.__seconds * 1000000000 + self.__nanoseconds

218 219 # comparison 220

221 - def __hash__(self):

222 return self.__seconds ^ self.__nanoseconds

223

224 - def __cmp__(self, other):

225 """ 226 Compare a value to a LIGOTimeGPS. If the value being compared 227 to the LIGOTimeGPS is not also a LIGOTimeGPS, then an attempt 228 is made to convert it to a LIGOTimeGPS. 229 230 Example: 231 232 >>> LIGOTimeGPS(100.5) < LIGOTimeGPS(200) 233 True 234 >>> LIGOTimeGPS(100.5) < 200 235 True 236 >>> LIGOTimeGPS(100.5) < "200" 237 True 238 """ 239 if not isinstance(other, LIGOTimeGPS): 240 try: 241 other = LIGOTimeGPS(other) 242 except TypeError: 243 return NotImplemented 244 return cmp(self.__seconds, other.seconds) or cmp(self.__nanoseconds, other.nanoseconds)

245

246 - def __nonzero__(self):

247 """ 248 Return True if the LIGOTimeGPS is nonzero. 249 250 Example: 251 252 >>> bool(LIGOTimeGPS(100.5)) 253 True 254 """ 255 return self.__seconds or self.__nanoseconds

256 257 # arithmetic 258

259 - def __add__(self, other):

260 """ 261 Add a value to a LIGOTimeGPS. If the value being added to the 262 LIGOTimeGPS is not also a LIGOTimeGPS, then an attempt is made 263 to convert it to a LIGOTimeGPS. 264 265 Example: 266 267 >>> LIGOTimeGPS(100.5) + LIGOTimeGPS(3) 268 LIGOTimeGPS(103, 500000000) 269 >>> LIGOTimeGPS(100.5) + 3 270 LIGOTimeGPS(103, 500000000) 271 >>> LIGOTimeGPS(100.5) + "3" 272 LIGOTimeGPS(103, 500000000) 273 """ 274 if not isinstance(other, LIGOTimeGPS): 275 other = LIGOTimeGPS(other) 276 return LIGOTimeGPS(self.__seconds + other.seconds, self.__nanoseconds + other.nanoseconds)

277 278 # addition is commutative. 279 __radd__ = __add__ 280

281 - def __sub__(self, other):

282 """ 283 Subtract a value from a LIGOTimeGPS. If the value being 284 subtracted from the LIGOTimeGPS is not also a LIGOTimeGPS, then 285 an attempt is made to convert it to a LIGOTimeGPS. 286 287 Example: 288 289 >>> LIGOTimeGPS(100.5) - LIGOTimeGPS(3) 290 LIGOTimeGPS(97, 500000000) 291 >>> LIGOTimeGPS(100.5) - 3 292 LIGOTimeGPS(97, 500000000) 293 >>> LIGOTimeGPS(100.5) - "3" 294 LIGOTimeGPS(97, 500000000) 295 """ 296 if not isinstance(other, LIGOTimeGPS): 297 other = LIGOTimeGPS(other) 298 return LIGOTimeGPS(self.__seconds - other.seconds, self.__nanoseconds - other.nanoseconds)

299

300 - def __rsub__(self, other):

301 """ 302 Subtract a LIGOTimeGPS from a value. 303 """ 304 if not isinstance(other, LIGOTimeGPS): 305 other = LIGOTimeGPS(other) 306 return LIGOTimeGPS(other.seconds - self.__seconds, other.nanoseconds - self.__nanoseconds)

307

308 - def __mul__(self, other):

309 """ 310 Multiply a LIGOTimeGPS by a number. 311 312 Example: 313 314 >>> LIGOTimeGPS(100.5) * 2 315 LIGOTimeGPS(201, 0) 316 """ 317 seconds = self.__seconds 318 nanoseconds = self.__nanoseconds 319 320 if seconds < 0 and nanoseconds > 0: 321 seconds += 1 322 nanoseconds -= 1000000000 323 elif seconds > 0 and nanoseconds < 0: 324 seconds -=1 325 nanoseconds += 1000000000 326 327 slo = seconds % 131072 328 shi = seconds - slo 329 olo = other % 2**(int(math.log(other, 2)) - 26) if other else 0 330 ohi = other - olo 331 332 nanoseconds *= float(other) 333 seconds = 0. 334 for addend in (slo * olo, shi * olo, slo * ohi, shi * ohi): 335 n, s = math.modf(addend) 336 seconds += s 337 nanoseconds += n * 1e9 338 339 return LIGOTimeGPS(seconds, round(nanoseconds))

340 341 # multiplication is commutative 342 __rmul__ = __mul__ 343

344 - def __div__(self, other):

345 """ 346 Divide a LIGOTimeGPS by a number. 347 348 Example: 349 350 >>> LIGOTimeGPS(100.5) / 2 351 LIGOTimeGPS(50, 250000000) 352 """ 353 quotient = LIGOTimeGPS(float(self) / float(other)) 354 for n in range(100): 355 residual = float(self - quotient * other) / float(other) 356 quotient += residual 357 if abs(residual) <= 0.5e-9: 358 break 359 return quotient

360

361 - def __mod__(self, other):

362 """ 363 Compute the remainder when a LIGOTimeGPS is divided by a number. 364 365 Example: 366 367 >>> LIGOTimeGPS(100.5) % 3 368 LIGOTimeGPS(1, 500000000) 369 """ 370 quotient = int(self / other) 371 return self - quotient * other

372 373 # unary arithmetic 374

375 - def __pos__(self):

376 return self

377

378 - def __neg__(self):

379 return LIGOTimeGPS(0, -self.ns())

380

381 - def __abs__(self):

382 if self.__seconds >= 0: 383 return self 384 return -self

385 386 387 # 388 # ============================================================================= 389 # 390 # LAL Cache File Manipulation 391 # 392 # ============================================================================= 393 # 394 395 396 import imp 397 lal = imp.load_module("lal", *imp.find_module("lal", sys.path[1:])) 398 lal.utils = imp.load_module("lal.utils", *imp.find_module("utils", lal.__path__))

399 -class CacheEntry(lal.utils.CacheEntry):

400 - def __init__(self, *args, **kwargs):

401 warnings.warn("glue.lal.CacheEntry is deprecated, use lal.utils.CacheEntry instead", DeprecationWarning) 402 return super(CacheEntry, self).__init__(*args, **kwargs)

403 404 405 # 406 # An object representing a LAL cache file 407 # 408 409

410 -class Cache(list):

411 """ 412 An object representing a LAL cache file. Currently it is possible to 413 add anything to a Cache. This method should check that the thing you 414 are adding is a CacheEntry and throw and error if it is not. 415 """ 416 entry_class = CacheEntry 417 418 # methods to create new Cache objects

419 - def fromfile(cls, fileobj, coltype=LIGOTimeGPS):

420 """ 421 Return a Cache object whose entries are read from an open file. 422 """ 423 c = [cls.entry_class(line, coltype=coltype) for line in fileobj] 424 return cls(c)

425 fromfile = classmethod(fromfile) 426

427 - def fromfilenames(cls, filenames, coltype=LIGOTimeGPS):

428 """ 429 Read Cache objects from the files named and concatenate the results into a 430 single Cache. 431 """ 432 cache = cls() 433 for filename in filenames: 434 cache.extend(cls.fromfile(open(filename), coltype=coltype)) 435 return cache

436 fromfilenames = classmethod(fromfilenames) 437

438 - def from_urls(cls, urllist, coltype=LIGOTimeGPS):

439 """ 440 Return a Cache whose entries are inferred from the URLs 441 in urllist, if possible. PFN lists will also work; for PFNs, the path 442 will be absolutized and "file://" and "localhost" will be assumed 443 for the schemes and hosts. 444 445 The filenames must be in the format set forth by DASWG in T050017-00. 446 """ 447 def pfn_to_url(url): 448 scheme, host, path, dummy, dummy = urllib.parse.urlsplit(url) 449 if scheme == "": path = os.path.abspath(path) 450 return urllib.parse.urlunsplit((scheme or "file", host or "localhost", 451 path, "", ""))

452 return cls([cls.entry_class.from_T050017(pfn_to_url(f), coltype=coltype) \ 453 for f in urllist])

454 from_urls = classmethod(from_urls) 455 456 # some set arithmetic to make life better

457 - def __isub__(self, other):

458 """ 459 Remove elements from self that are in other. 460 """ 461 end = len(self) - 1 462 for i, elem in enumerate(self[::-1]): 463 if elem in other: 464 del self[end - i] 465 return self

466

467 - def __sub__(self, other):

468 """ 469 Return a Cache containing the entries of self that are not in other. 470 """ 471 return self.__class__([elem for elem in self if elem not in other])

472

473 - def __ior__(self, other):

474 """ 475 Append entries from other onto self without introducing (new) duplicates. 476 """ 477 self.extend(other - self) 478 return self

479

480 - def __or__(self, other):

481 """ 482 Return a Cache containing all entries of self and other. 483 """ 484 return self.__class__(self[:]).__ior__(other)

485

486 - def __iand__(self, other):

487 """ 488 Remove elements in self that are not in other. 489 """ 490 end = len(self) - 1 491 for i, elem in enumerate(self[::-1]): 492 if elem not in other: 493 del self[end - i] 494 return self

495

496 - def __and__(self, other):

497 """ 498 Return a Cache containing the entries of self that are also in other. 499 """ 500 return self.__class__([elem for elem in self if elem in other])

501

502 - def unique(self):

503 """ 504 Return a Cache which has every element of self, but without 505 duplication. Preserve order. Does not hash, so a bit slow. 506 """ 507 seen = set() 508 return self.__class__([x for x in self if x not in seen and not seen.add(x)])

509 #new = self.__class__([]) 510 #for elem in self: 511 # if elem not in new: 512 # new.append(elem) 513 #return new 514 515 # other useful manipulations

516 - def tofile(self, fileobj):

517 """ 518 write a cache object to the fileobj as a lal cache file 519 """ 520 for entry in self: 521 fileobj.write('%s\n' % str(entry)) 522 fileobj.close()

523

524 - def topfnfile(self, fileobj):

525 """ 526 write a cache object to filename as a plain text pfn file 527 """ 528 for entry in self: 529 fileobj.write('%s\n' % entry.path) 530 fileobj.close()

531

532 - def to_segmentlistdict(self):

533 """ 534 Return a segmentlistdict object describing the instruments 535 and times spanned by the entries in this Cache. The return 536 value is coalesced. 537 """ 538 d = segments.segmentlistdict() 539 for entry in self: 540 d |= entry.segmentlistdict 541 return d

542

543 - def sieve(self, ifos=None, description=None, segment=None, 544 segmentlist=None, exact_match=False):

545 """ 546 Return a Cache object with those CacheEntries that 547 contain the given patterns (or overlap, in the case of 548 segment or segmentlist). If exact_match is True, then 549 non-None ifos, description, and segment patterns must 550 match exactly, and a non-None segmentlist must contain 551 a segment which matches exactly). 552 553 It makes little sense to specify both segment and 554 segmentlist arguments, but it is not prohibited. 555 556 Bash-style wildcards (*?) are allowed for ifos and description. 557 """ 558 if exact_match: 559 segment_func = lambda e: e.segment == segment 560 segmentlist_func = lambda e: e.segment in segmentlist 561 else: 562 if ifos is not None: ifos = "*" + ifos + "*" 563 if description is not None: description = "*" + description + "*" 564 segment_func = lambda e: segment.intersects(e.segment) 565 segmentlist_func = lambda e: segmentlist.intersects_segment(e.segment) 566 567 c = self 568 569 if ifos is not None: 570 ifos_regexp = re.compile(fnmatch.translate(ifos)) 571 c = [entry for entry in c if ifos_regexp.match(entry.observatory) is not None] 572 573 if description is not None: 574 descr_regexp = re.compile(fnmatch.translate(description)) 575 c = [entry for entry in c if descr_regexp.match(entry.description) is not None] 576 577 if segment is not None: 578 c = [entry for entry in c if segment_func(entry)] 579 580 if segmentlist is not None: 581 # must coalesce for intersects_segment() to work 582 segmentlist.coalesce() 583 c = [entry for entry in c if segmentlist_func(entry)] 584 585 return self.__class__(c)

586

587 - def pfnlist(self):

588 """ 589 Return a list of physical file names. 590 """ 591 return [entry.path for entry in self]

592

593 - def checkfilesexist(self, on_missing="warn"):

594 ''' 595 Runs through the entries of the Cache() object and checks each entry 596 if the file which it points to exists or not. If the file does exist then 597 it adds the entry to the Cache() object containing found files, otherwise it 598 adds the entry to the Cache() object containing all entries that are missing. 599 It returns both in the follwing order: Cache_Found, Cache_Missed. 600 601 Pass on_missing to control how missing files are handled: 602 "warn": print a warning message saying how many files 603 are missing out of the total checked. 604 "error": raise an exception if any are missing 605 "ignore": do nothing 606 ''' 607 if on_missing not in ("warn", "error", "ignore"): 608 raise ValueError("on_missing must be \"warn\", \"error\", or \"ignore\".") 609 610 c_found = [] 611 c_missed = [] 612 for entry in self: 613 if os.path.isfile(entry.path): 614 c_found.append(entry) 615 else: 616 c_missed.append(entry) 617 618 if len(c_missed) > 0: 619 msg = "%d of %d files in the cache were not found "\ 620 "on disk" % (len(c_missed), len(self)) 621 if on_missing == "warn": 622 sys.stderr.write("warning: %s\n" % msg) 623 elif on_missing == "error": 624 raise ValueError(msg) 625 elif on_missing == "ignore": 626 pass 627 else: 628 raise ValueError("Why am I here? "\ 629 "Please file a bug report!") 630 return self.__class__(c_found), self.__class__(c_missed)

631

632 - def __getslice__(self, i, j):

633 return self.__class__(super(Cache, self).__getslice__(i, j))

634

Source Code for Module glue.lal