Collection of scripts and small programs used by the EVN Support Scientists at JIVE during the regular data processing of EVN observations.
You can not select more than 25 topics Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
 
 

75 lines
3.0 KiB

  1. #!/usr/bin/env python3
  2. ### --- scale1bit.py --- ####
  3. ### JMB Oct 2018 ###
  4. ### TACL by HV Oct 2018 ###
  5. ### Scales 1 bit data to account for quantization loss factors
  6. import pyrap.tables
  7. import math
  8. import argparse
  9. import sys #for exit
  10. import os.path
  11. #old glish comments for safekeeping
  12. # value of factor updated 20nov2015: \sqrt{(pi/2)/1.1329552}
  13. ### 1.1329552 = 2bit-2bit correction factor already present in SFXC output
  14. ### pi/2 = what would be the case for 1bit-1bit data
  15. ### \sqrt = to get to 2bit-1bit factor, which forms the basis for the
  16. ### existing logic
  17. ## factor := 1.17631;
  18. ## factor := 1.177480083;
  19. factor1b1b = math.pi /2.0 / 1.1329552
  20. factor1b2b = math.sqrt(factor1b1b)
  21. # debug function for feedback
  22. def debug(*msg, **kwargs):
  23. if args.verbose:
  24. print (*msg, file=sys.stderr, **kwargs)
  25. #handle command line arguments: ms, antenna
  26. parser = argparse.ArgumentParser(description='Scale 1 bit data to correct quantization losses')
  27. parser.add_argument('ms', help="The measurement set to be corrected")
  28. parser.add_argument('ant',nargs='+', help="The antenna(s) to correct (space delimited)")
  29. parser.add_argument('-v', '--verbose', help="Print debug information", action="store_true")
  30. parser.add_argument('-u', '--undo', help="Undo a previous run, take care to specify the same stations!", action="store_true")
  31. parser.add_argument('-w', '--scale-weights', help="Also scale the weights", dest='to_scale', default=['DATA'], action='append_const', const='WEIGHT')
  32. args = parser.parse_args()
  33. if not os.path.exists(args.ms):
  34. print("Error, no such ms found!")
  35. sys.exit(1)
  36. try:
  37. # query antenna table to translate 1 or more antenna names into antenna id's
  38. aList = list(pyrap.tables.taql("""SELECT ROWID() AS ANTENNAS FROM {table}::ANTENNA
  39. WHERE UPCASE(NAME) IN {0}""".format(
  40. list(map(str.upper, args.ant)),table=args.ms)
  41. ).getcol('ANTENNAS'))
  42. except Exception as e:
  43. print(e, "AAAAAAAAAARGH we have no such dishes, abort fail die :'(")
  44. sys.exit(1)
  45. debug("Antenna list:", aList)
  46. # depending on length of the result, use "ANTENNAx == <id>" or "ANTENNAx IN [<id>, <id>,...]"
  47. # because "==" is faster than "IN [<id>]" when only one <id> is present
  48. antcond = ("== {0[0]}".format if len(aList)==1 else "IN {0}".format)(aList)
  49. # and let taql do the real work ...
  50. factor1, factor2 = (factor1b1b, factor1b2b) if not args.undo else (1/factor1b1b, 1/factor1b2b)
  51. scale_it = "{{0}} = {{0}} * IIF(apply.ant1bit AND apply.ant2bit, {factor1b1b}, {factor1b2b})".format(factor1b1b=factor1, factor1b2b=factor2).format
  52. #to_scale = ['DATA'] + (['WEIGHT'] if args.scale_weights else [])
  53. t = pyrap.tables.taql("""
  54. UPDATE {table}
  55. SET {todo}
  56. FROM [SELECT ANTENNA1 {condition} AS ant1bit, ANTENNA2 {condition} AS ant2bit
  57. FROM {table} GIVING AS memory] apply
  58. WHERE ((apply.ant1bit OR apply.ant2bit) and ANTENNA1 != ANTENNA2)
  59. """.format(condition=antcond, table=args.ms, todo=",".join(map(scale_it, args.to_scale))))