root / gr-radio-astronomy / src / python / usrp_ra_receiver.py @ 38fe0d40
History | View | Annotate | Download (40.7 kB)
| 1 | #!/usr/bin/env python
|
|---|---|
| 2 | #
|
| 3 | # Copyright 2004,2005,2007 Free Software Foundation, Inc.
|
| 4 | #
|
| 5 | # This file is part of GNU Radio
|
| 6 | #
|
| 7 | # GNU Radio is free software; you can redistribute it and/or modify
|
| 8 | # it under the terms of the GNU General Public License as published by
|
| 9 | # the Free Software Foundation; either version 3, or (at your option)
|
| 10 | # any later version.
|
| 11 | #
|
| 12 | # GNU Radio is distributed in the hope that it will be useful,
|
| 13 | # but WITHOUT ANY WARRANTY; without even the implied warranty of
|
| 14 | # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
| 15 | # GNU General Public License for more details.
|
| 16 | #
|
| 17 | # You should have received a copy of the GNU General Public License
|
| 18 | # along with GNU Radio; see the file COPYING. If not, write to
|
| 19 | # the Free Software Foundation, Inc., 51 Franklin Street,
|
| 20 | # Boston, MA 02110-1301, USA.
|
| 21 | #
|
| 22 | |
| 23 | from gnuradio import gr, gru |
| 24 | from gnuradio import usrp |
| 25 | from usrpm import usrp_dbid |
| 26 | from gnuradio import eng_notation |
| 27 | from gnuradio.eng_option import eng_option |
| 28 | from gnuradio.wxgui import stdgui2, ra_fftsink, ra_stripchartsink, ra_waterfallsink, form, slider |
| 29 | from optparse import OptionParser |
| 30 | import wx |
| 31 | import sys |
| 32 | import Numeric |
| 33 | import time |
| 34 | import numpy.fft |
| 35 | import ephem |
| 36 | |
| 37 | class app_flow_graph(stdgui2.std_top_block): |
| 38 | def __init__(self, frame, panel, vbox, argv): |
| 39 | stdgui2.std_top_block.__init__(self, frame, panel, vbox, argv)
|
| 40 | |
| 41 | self.frame = frame
|
| 42 | self.panel = panel
|
| 43 | |
| 44 | parser = OptionParser(option_class=eng_option) |
| 45 | parser.add_option("-R", "--rx-subdev-spec", type="subdev", default=(0, 0), |
| 46 | help="select USRP Rx side A or B (default=A)")
|
| 47 | parser.add_option("-d", "--decim", type="int", default=16, |
| 48 | help="set fgpa decimation rate to DECIM [default=%default]")
|
| 49 | parser.add_option("-f", "--freq", type="eng_float", default=None, |
| 50 | help="set frequency to FREQ", metavar="FREQ") |
| 51 | parser.add_option("-a", "--avg", type="eng_float", default=1.0, |
| 52 | help="set spectral averaging alpha")
|
| 53 | parser.add_option("-i", "--integ", type="eng_float", default=1.0, |
| 54 | help="set integration time")
|
| 55 | parser.add_option("-g", "--gain", type="eng_float", default=None, |
| 56 | help="set gain in dB (default is midpoint)")
|
| 57 | parser.add_option("-l", "--reflevel", type="eng_float", default=30.0, |
| 58 | help="Set Total power reference level")
|
| 59 | parser.add_option("-y", "--division", type="eng_float", default=0.5, |
| 60 | help="Set Total power Y division size")
|
| 61 | parser.add_option("-e", "--longitude", type="eng_float", default=-76.02,help="Set Observer Longitude") |
| 62 | parser.add_option("-c", "--latitude", type="eng_float", default=44.85,help="Set Observer Latitude") |
| 63 | parser.add_option("-o", "--observing", type="eng_float", default=0.0, |
| 64 | help="Set observing frequency")
|
| 65 | parser.add_option("-x", "--ylabel", default="dB", help="Y axis label") |
| 66 | parser.add_option("-z", "--divbase", type="eng_float", default=0.025, help="Y Division increment base") |
| 67 | parser.add_option("-v", "--stripsize", type="eng_float", default=2400, help="Size of stripchart, in 2Hz samples") |
| 68 | parser.add_option("-F", "--fft_size", type="eng_float", default=1024, help="Size of FFT") |
| 69 | parser.add_option("-N", "--decln", type="eng_float", default=999.99, help="Observing declination") |
| 70 | parser.add_option("-X", "--prefix", default="./") |
| 71 | parser.add_option("-M", "--fft_rate", type="eng_float", default=8.0, help="FFT Rate") |
| 72 | parser.add_option("-A", "--calib_coeff", type="eng_float", default=1.0, help="Calibration coefficient") |
| 73 | parser.add_option("-B", "--calib_offset", type="eng_float", default=0.0, help="Calibration coefficient") |
| 74 | parser.add_option("-W", "--waterfall", action="store_true", default=False, help="Use Waterfall FFT display") |
| 75 | parser.add_option("-S", "--setimode", action="store_true", default=False, help="Enable SETI processing of spectral data") |
| 76 | parser.add_option("-K", "--setik", type="eng_float", default=1.5, help="K value for SETI analysis") |
| 77 | parser.add_option("-T", "--setibandwidth", type="eng_float", default=12500, help="Instantaneous SETI observing bandwidth--must be divisor of 250Khz") |
| 78 | parser.add_option("-Q", "--seti_range", type="eng_float", default=1.0e6, help="Total scan width, in Hz for SETI scans") |
| 79 | parser.add_option("-Z", "--dual_mode", action="store_true", |
| 80 | default=False, help="Dual-polarization mode") |
| 81 | parser.add_option("-I", "--interferometer", action="store_true", default=False, help="Interferometer mode") |
| 82 | parser.add_option("-D", "--switch_mode", action="store_true", default=False, help="Dicke Switching mode") |
| 83 | parser.add_option("-P", "--reference_divisor", type="eng_float", default=1.0, help="Reference Divisor") |
| 84 | parser.add_option("-U", "--ref_fifo", default="@@@@") |
| 85 | parser.add_option("-n", "--notches", action="store_true", |
| 86 | default=False, help="Notch frequencies after all other args") |
| 87 | (options, args) = parser.parse_args() |
| 88 | |
| 89 | self.setimode = options.setimode
|
| 90 | self.dual_mode = options.dual_mode
|
| 91 | self.interferometer = options.interferometer
|
| 92 | self.normal_mode = False |
| 93 | self.switch_mode = options.switch_mode
|
| 94 | self.switch_state = 0 |
| 95 | self.reference_divisor = options.reference_divisor
|
| 96 | self.ref_fifo = options.ref_fifo
|
| 97 | |
| 98 | self.NOTCH_TAPS = 16 |
| 99 | self.notches = Numeric.zeros(self.NOTCH_TAPS,Numeric.Float64) |
| 100 | # Get notch locations
|
| 101 | j = 0
|
| 102 | for i in args: |
| 103 | self.notches[j] = float(i) |
| 104 | j = j + 1
|
| 105 | |
| 106 | self.use_notches = options.notches
|
| 107 | |
| 108 | if (self.ref_fifo != "@@@@"): |
| 109 | self.ref_fifo_file = open (self.ref_fifo, "w") |
| 110 | |
| 111 | modecount = 0
|
| 112 | for modes in (self.dual_mode, self.interferometer): |
| 113 | if (modes == True): |
| 114 | modecount = modecount + 1
|
| 115 | |
| 116 | if (modecount > 1): |
| 117 | print "must select only 1 of --dual_mode, or --interferometer" |
| 118 | sys.exit(1)
|
| 119 | |
| 120 | self.chartneeded = True |
| 121 | |
| 122 | if (self.setimode == True): |
| 123 | self.chartneeded = False |
| 124 | |
| 125 | if (self.setimode == True and self.interferometer == True): |
| 126 | print "can't pick both --setimode and --interferometer" |
| 127 | sys.exit(1)
|
| 128 | |
| 129 | if (self.setimode == True and self.switch_mode == True): |
| 130 | print "can't pick both --setimode and --switch_mode" |
| 131 | sys.exit(1)
|
| 132 | |
| 133 | if (self.interferometer == True and self.switch_mode == True): |
| 134 | print "can't pick both --interferometer and --switch_mode" |
| 135 | sys.exit(1)
|
| 136 | |
| 137 | if (modecount == 0): |
| 138 | self.normal_mode = True |
| 139 | |
| 140 | self.show_debug_info = True |
| 141 | |
| 142 | # Pick up waterfall option
|
| 143 | self.waterfall = options.waterfall
|
| 144 | |
| 145 | # SETI mode stuff
|
| 146 | self.setimode = options.setimode
|
| 147 | self.seticounter = 0 |
| 148 | self.setik = options.setik
|
| 149 | self.seti_fft_bandwidth = int(options.setibandwidth) |
| 150 | |
| 151 | # Calculate binwidth
|
| 152 | binwidth = self.seti_fft_bandwidth / options.fft_size
|
| 153 | |
| 154 | # Use binwidth, and knowledge of likely chirp rates to set reasonable
|
| 155 | # values for SETI analysis code. We assume that SETI signals will
|
| 156 | # chirp at somewhere between 0.10Hz/sec and 0.25Hz/sec.
|
| 157 | #
|
| 158 | # upper_limit is the "worst case"--that is, the case for which we have
|
| 159 | # to wait the longest to actually see any drift, due to the quantizing
|
| 160 | # on FFT bins.
|
| 161 | upper_limit = binwidth / 0.10
|
| 162 | self.setitimer = int(upper_limit * 2.00) |
| 163 | self.scanning = True |
| 164 | |
| 165 | # Calculate the CHIRP values based on Hz/sec
|
| 166 | self.CHIRP_LOWER = 0.10 * self.setitimer |
| 167 | self.CHIRP_UPPER = 0.25 * self.setitimer |
| 168 | |
| 169 | # Reset hit counters to 0
|
| 170 | self.hitcounter = 0 |
| 171 | self.s1hitcounter = 0 |
| 172 | self.s2hitcounter = 0 |
| 173 | self.avgdelta = 0 |
| 174 | # We scan through 2Mhz of bandwidth around the chosen center freq
|
| 175 | self.seti_freq_range = options.seti_range
|
| 176 | # Calculate lower edge
|
| 177 | self.setifreq_lower = options.freq - (self.seti_freq_range/2) |
| 178 | self.setifreq_current = options.freq
|
| 179 | # Calculate upper edge
|
| 180 | self.setifreq_upper = options.freq + (self.seti_freq_range/2) |
| 181 | |
| 182 | # Maximum "hits" in a line
|
| 183 | self.nhits = 20 |
| 184 | |
| 185 | # Number of lines for analysis
|
| 186 | self.nhitlines = 4 |
| 187 | |
| 188 | # We change center frequencies based on nhitlines and setitimer
|
| 189 | self.setifreq_timer = self.setitimer * (self.nhitlines * 5) |
| 190 | |
| 191 | # Create actual timer
|
| 192 | self.seti_then = time.time()
|
| 193 | |
| 194 | # The hits recording array
|
| 195 | self.hits_array = Numeric.zeros((self.nhits,self.nhitlines), Numeric.Float64) |
| 196 | self.hit_intensities = Numeric.zeros((self.nhits,self.nhitlines), Numeric.Float64) |
| 197 | # Calibration coefficient and offset
|
| 198 | self.calib_coeff = options.calib_coeff
|
| 199 | self.calib_offset = options.calib_offset
|
| 200 | if self.calib_offset < -750: |
| 201 | self.calib_offset = -750 |
| 202 | if self.calib_offset > 750: |
| 203 | self.calib_offset = 750 |
| 204 | |
| 205 | if self.calib_coeff < 1: |
| 206 | self.calib_coeff = 1 |
| 207 | if self.calib_coeff > 100: |
| 208 | self.calib_coeff = 100 |
| 209 | |
| 210 | self.integ = options.integ
|
| 211 | self.avg_alpha = options.avg
|
| 212 | self.gain = options.gain
|
| 213 | self.decln = options.decln
|
| 214 | |
| 215 | # Set initial values for datalogging timed-output
|
| 216 | self.continuum_then = time.time()
|
| 217 | self.spectral_then = time.time()
|
| 218 | |
| 219 | |
| 220 | # build the graph
|
| 221 | |
| 222 | self.subdev = [(0, 0), (0,0)] |
| 223 | |
| 224 | #
|
| 225 | # If SETI mode, we always run at maximum USRP decimation
|
| 226 | #
|
| 227 | if (self.setimode): |
| 228 | options.decim = 256
|
| 229 | |
| 230 | if (self.dual_mode == False and self.interferometer == False): |
| 231 | self.u = usrp.source_c(decim_rate=options.decim)
|
| 232 | self.u.set_mux(usrp.determine_rx_mux_value(self.u, options.rx_subdev_spec)) |
| 233 | # determine the daughterboard subdevice we're using
|
| 234 | self.subdev[0] = usrp.selected_subdev(self.u, options.rx_subdev_spec) |
| 235 | self.subdev[1] = self.subdev[0] |
| 236 | self.cardtype = self.subdev[0].dbid() |
| 237 | else:
|
| 238 | self.u=usrp.source_c(decim_rate=options.decim, nchan=2) |
| 239 | self.subdev[0] = usrp.selected_subdev(self.u, (0, 0)) |
| 240 | self.subdev[1] = usrp.selected_subdev(self.u, (1, 0)) |
| 241 | self.cardtype = self.subdev[0].dbid() |
| 242 | self.u.set_mux(0x32103210) |
| 243 | c1 = self.subdev[0].name() |
| 244 | c2 = self.subdev[1].name() |
| 245 | if (c1 != c2):
|
| 246 | print "Must have identical cardtypes for --dual_mode or --interferometer" |
| 247 | sys.exit(1)
|
| 248 | #
|
| 249 | # Set 8-bit mode
|
| 250 | #
|
| 251 | width = 8
|
| 252 | shift = 8
|
| 253 | format = self.u.make_format(width, shift)
|
| 254 | r = self.u.set_format(format)
|
| 255 | |
| 256 | # Set initial declination
|
| 257 | self.decln = options.decln
|
| 258 | |
| 259 | input_rate = self.u.adc_freq() / self.u.decim_rate() |
| 260 | self.bw = input_rate
|
| 261 | #
|
| 262 | # Set prefix for data files
|
| 263 | #
|
| 264 | self.prefix = options.prefix
|
| 265 | |
| 266 | #
|
| 267 | # The lower this number, the fewer sample frames are dropped
|
| 268 | # in computing the FFT. A sampled approach is taken to
|
| 269 | # computing the FFT of the incoming data, which reduces
|
| 270 | # sensitivity. Increasing sensitivity inreases CPU loading.
|
| 271 | #
|
| 272 | self.fft_rate = options.fft_rate
|
| 273 | |
| 274 | self.fft_size = int(options.fft_size) |
| 275 | |
| 276 | # This buffer is used to remember the most-recent FFT display
|
| 277 | # values. Used later by self.write_spectral_data() to write
|
| 278 | # spectral data to datalogging files, and by the SETI analysis
|
| 279 | # function.
|
| 280 | #
|
| 281 | self.fft_outbuf = Numeric.zeros(self.fft_size, Numeric.Float64) |
| 282 | |
| 283 | #
|
| 284 | # If SETI mode, only look at seti_fft_bandwidth
|
| 285 | # at a time.
|
| 286 | #
|
| 287 | if (self.setimode): |
| 288 | self.fft_input_rate = self.seti_fft_bandwidth |
| 289 | |
| 290 | #
|
| 291 | # Build a decimating bandpass filter
|
| 292 | #
|
| 293 | self.fft_input_taps = gr.firdes.complex_band_pass (1.0, |
| 294 | input_rate, |
| 295 | -(int(self.fft_input_rate/2)), int(self.fft_input_rate/2), 200, |
| 296 | gr.firdes.WIN_HAMMING, 0)
|
| 297 | |
| 298 | #
|
| 299 | # Compute required decimation factor
|
| 300 | #
|
| 301 | decimation = int(input_rate/self.fft_input_rate) |
| 302 | self.fft_bandpass = gr.fir_filter_ccc (decimation,
|
| 303 | self.fft_input_taps)
|
| 304 | else:
|
| 305 | self.fft_input_rate = input_rate
|
| 306 | |
| 307 | # Set up FFT display
|
| 308 | if self.waterfall == False: |
| 309 | self.scope = ra_fftsink.ra_fft_sink_c (panel,
|
| 310 | fft_size=int(self.fft_size), sample_rate=self.fft_input_rate, |
| 311 | fft_rate=int(self.fft_rate), title="Spectral", |
| 312 | ofunc=self.fft_outfunc, xydfunc=self.xydfunc) |
| 313 | else:
|
| 314 | self.scope = ra_waterfallsink.waterfall_sink_c (panel,
|
| 315 | fft_size=int(self.fft_size), sample_rate=self.fft_input_rate, |
| 316 | fft_rate=int(self.fft_rate), title="Spectral", ofunc=self.fft_outfunc, size=(1100, 600), xydfunc=self.xydfunc, ref_level=0, span=10) |
| 317 | |
| 318 | # Set up ephemeris data
|
| 319 | self.locality = ephem.Observer()
|
| 320 | self.locality.long = str(options.longitude) |
| 321 | self.locality.lat = str(options.latitude) |
| 322 | |
| 323 | # We make notes about Sunset/Sunrise in Continuum log files
|
| 324 | self.sun = ephem.Sun()
|
| 325 | self.sunstate = "??" |
| 326 | |
| 327 | # Set up stripchart display
|
| 328 | tit = "Continuum"
|
| 329 | if (self.dual_mode != False): |
| 330 | tit = "H+V Continuum"
|
| 331 | if (self.interferometer != False): |
| 332 | tit = "East x West Correlation"
|
| 333 | self.stripsize = int(options.stripsize) |
| 334 | if self.chartneeded == True: |
| 335 | self.chart = ra_stripchartsink.stripchart_sink_f (panel,
|
| 336 | stripsize=self.stripsize,
|
| 337 | title=tit, |
| 338 | xlabel="LMST Offset (Seconds)",
|
| 339 | scaling=1.0, ylabel=options.ylabel,
|
| 340 | divbase=options.divbase) |
| 341 | |
| 342 | # Set center frequency
|
| 343 | self.centerfreq = options.freq
|
| 344 | |
| 345 | # Set observing frequency (might be different from actual programmed
|
| 346 | # RF frequency)
|
| 347 | if options.observing == 0.0: |
| 348 | self.observing = options.freq
|
| 349 | else:
|
| 350 | self.observing = options.observing
|
| 351 | |
| 352 | # Remember our input bandwidth
|
| 353 | self.bw = input_rate
|
| 354 | |
| 355 | #
|
| 356 | #
|
| 357 | # The strip chart is fed at a constant 1Hz rate
|
| 358 | #
|
| 359 | |
| 360 | #
|
| 361 | # Call constructors for receive chains
|
| 362 | #
|
| 363 | |
| 364 | if (self.dual_mode == True): |
| 365 | self.setup_dual (self.setimode,self.use_notches) |
| 366 | |
| 367 | if (self.interferometer == True): |
| 368 | self.setup_interferometer(self.setimode,self.use_notches) |
| 369 | |
| 370 | if (self.normal_mode == True): |
| 371 | self.setup_normal(self.setimode,self.use_notches) |
| 372 | |
| 373 | if (self.setimode == True): |
| 374 | self.setup_seti()
|
| 375 | |
| 376 | self._build_gui(vbox)
|
| 377 | |
| 378 | # Make GUI agree with command-line
|
| 379 | self.integ = options.integ
|
| 380 | if self.setimode == False: |
| 381 | self.myform['integration'].set_value(int(options.integ)) |
| 382 | self.myform['offset'].set_value(self.calib_offset) |
| 383 | self.myform['dcgain'].set_value(self.calib_coeff) |
| 384 | self.myform['average'].set_value(int(options.avg)) |
| 385 | |
| 386 | |
| 387 | if self.setimode == False: |
| 388 | # Make integrator agree with command line
|
| 389 | self.set_integration(int(options.integ)) |
| 390 | |
| 391 | self.avg_alpha = options.avg
|
| 392 | |
| 393 | # Make spectral averager agree with command line
|
| 394 | if options.avg != 1.0: |
| 395 | self.scope.set_avg_alpha(float(1.0/options.avg)) |
| 396 | self.scope.set_average(True) |
| 397 | |
| 398 | if self.setimode == False: |
| 399 | # Set division size
|
| 400 | self.chart.set_y_per_div(options.division)
|
| 401 | # Set reference(MAX) level
|
| 402 | self.chart.set_ref_level(options.reflevel)
|
| 403 | |
| 404 | # set initial values
|
| 405 | |
| 406 | if options.gain is None: |
| 407 | # if no gain was specified, use the mid-point in dB
|
| 408 | g = self.subdev[0].gain_range() |
| 409 | options.gain = float(g[0]+g[1])/2 |
| 410 | |
| 411 | if options.freq is None: |
| 412 | # if no freq was specified, use the mid-point
|
| 413 | r = self.subdev[0].freq_range() |
| 414 | options.freq = float(r[0]+r[1])/2 |
| 415 | |
| 416 | # Set the initial gain control
|
| 417 | self.set_gain(options.gain)
|
| 418 | |
| 419 | if not(self.set_freq(options.freq)): |
| 420 | self._set_status_msg("Failed to set initial frequency") |
| 421 | |
| 422 | # Set declination
|
| 423 | self.set_decln (self.decln) |
| 424 | |
| 425 | |
| 426 | # RF hardware information
|
| 427 | self.myform['decim'].set_value(self.u.decim_rate()) |
| 428 | self.myform['USB BW'].set_value(self.u.adc_freq() / self.u.decim_rate()) |
| 429 | if (self.dual_mode == True): |
| 430 | self.myform['dbname'].set_value(self.subdev[0].name()+'/'+self.subdev[1].name()) |
| 431 | else:
|
| 432 | self.myform['dbname'].set_value(self.subdev[0].name()) |
| 433 | |
| 434 | # Set analog baseband filtering, if DBS_RX
|
| 435 | if self.cardtype in (usrp_dbid.DBS_RX, usrp_dbid.DBS_RX_REV_2_1): |
| 436 | lbw = (self.u.adc_freq() / self.u.decim_rate()) / 2 |
| 437 | if lbw < 1.0e6: |
| 438 | lbw = 1.0e6
|
| 439 | #self.subdev[0].set_bw(lbw)
|
| 440 | #self.subdev[1].set_bw(lbw)
|
| 441 | |
| 442 | # Start the timer for the LMST display and datalogging
|
| 443 | self.lmst_timer.Start(1000) |
| 444 | if (self.switch_mode == True): |
| 445 | self.other_timer.Start(330) |
| 446 | |
| 447 | |
| 448 | def _set_status_msg(self, msg): |
| 449 | self.frame.GetStatusBar().SetStatusText(msg, 0) |
| 450 | |
| 451 | def _build_gui(self, vbox): |
| 452 | |
| 453 | def _form_set_freq(kv): |
| 454 | # Adjust current SETI frequency, and limits
|
| 455 | self.setifreq_lower = kv['freq'] - (self.seti_freq_range/2) |
| 456 | self.setifreq_current = kv['freq'] |
| 457 | self.setifreq_upper = kv['freq'] + (self.seti_freq_range/2) |
| 458 | |
| 459 | # Reset SETI analysis timer
|
| 460 | self.seti_then = time.time()
|
| 461 | # Zero-out hits array when changing frequency
|
| 462 | self.hits_array[:,:] = 0.0 |
| 463 | self.hit_intensities[:,:] = -60.0 |
| 464 | |
| 465 | return self.set_freq(kv['freq']) |
| 466 | |
| 467 | def _form_set_decln(kv): |
| 468 | return self.set_decln(kv['decln']) |
| 469 | |
| 470 | # Position the FFT display
|
| 471 | vbox.Add(self.scope.win, 15, wx.EXPAND) |
| 472 | |
| 473 | if self.setimode == False: |
| 474 | # Position the Total-power stripchart
|
| 475 | vbox.Add(self.chart.win, 15, wx.EXPAND) |
| 476 | |
| 477 | # add control area at the bottom
|
| 478 | self.myform = myform = form.form()
|
| 479 | hbox = wx.BoxSizer(wx.HORIZONTAL) |
| 480 | hbox.Add((7,0), 0, wx.EXPAND) |
| 481 | vbox1 = wx.BoxSizer(wx.VERTICAL) |
| 482 | myform['freq'] = form.float_field(
|
| 483 | parent=self.panel, sizer=vbox1, label="Center freq", weight=1, |
| 484 | callback=myform.check_input_and_call(_form_set_freq, self._set_status_msg))
|
| 485 | |
| 486 | vbox1.Add((4,0), 0, 0) |
| 487 | |
| 488 | myform['lmst_high'] = form.static_text_field(
|
| 489 | parent=self.panel, sizer=vbox1, label="Current LMST", weight=1) |
| 490 | vbox1.Add((4,0), 0, 0) |
| 491 | |
| 492 | if self.setimode == False: |
| 493 | myform['spec_data'] = form.static_text_field(
|
| 494 | parent=self.panel, sizer=vbox1, label="Spectral Cursor", weight=1) |
| 495 | vbox1.Add((4,0), 0, 0) |
| 496 | |
| 497 | vbox2 = wx.BoxSizer(wx.VERTICAL) |
| 498 | if self.setimode == False: |
| 499 | vbox3 = wx.BoxSizer(wx.VERTICAL) |
| 500 | g = self.subdev[0].gain_range() |
| 501 | myform['gain'] = form.slider_field(parent=self.panel, sizer=vbox2, label="RF Gain", |
| 502 | weight=1,
|
| 503 | min=int(g[0]), max=int(g[1]), |
| 504 | callback=self.set_gain)
|
| 505 | |
| 506 | vbox2.Add((4,0), 0, 0) |
| 507 | if self.setimode == True: |
| 508 | max_savg = 100
|
| 509 | else:
|
| 510 | max_savg = 3000
|
| 511 | myform['average'] = form.slider_field(parent=self.panel, sizer=vbox2, |
| 512 | label="Spectral Averaging (FFT frames)", weight=1, min=1, max=max_savg, callback=self.set_averaging) |
| 513 | |
| 514 | # Set up scan control button when in SETI mode
|
| 515 | if (self.setimode == True): |
| 516 | # SETI scanning control
|
| 517 | buttonbox = wx.BoxSizer(wx.HORIZONTAL) |
| 518 | self.scan_control = form.button_with_callback(self.panel, |
| 519 | label="Scan: On ",
|
| 520 | callback=self.toggle_scanning)
|
| 521 | |
| 522 | buttonbox.Add(self.scan_control, 0, wx.CENTER) |
| 523 | vbox2.Add(buttonbox, 0, wx.CENTER)
|
| 524 | |
| 525 | vbox2.Add((4,0), 0, 0) |
| 526 | |
| 527 | if self.setimode == False: |
| 528 | myform['integration'] = form.slider_field(parent=self.panel, sizer=vbox2, |
| 529 | label="Continuum Integration Time (sec)", weight=1, min=1, max=180, callback=self.set_integration) |
| 530 | |
| 531 | vbox2.Add((4,0), 0, 0) |
| 532 | |
| 533 | myform['decln'] = form.float_field(
|
| 534 | parent=self.panel, sizer=vbox2, label="Current Declination", weight=1, |
| 535 | callback=myform.check_input_and_call(_form_set_decln)) |
| 536 | vbox2.Add((4,0), 0, 0) |
| 537 | |
| 538 | if self.setimode == False: |
| 539 | myform['offset'] = form.slider_field(parent=self.panel, sizer=vbox3, |
| 540 | label="Post-Detector Offset", weight=1, min=-750, max=750, |
| 541 | callback=self.set_pd_offset)
|
| 542 | vbox3.Add((2,0), 0, 0) |
| 543 | myform['dcgain'] = form.slider_field(parent=self.panel, sizer=vbox3, |
| 544 | label="Post-Detector Gain", weight=1, min=1, max=100, |
| 545 | callback=self.set_pd_gain)
|
| 546 | vbox3.Add((2,0), 0, 0) |
| 547 | hbox.Add(vbox1, 0, 0) |
| 548 | hbox.Add(vbox2, wx.ALIGN_RIGHT, 0)
|
| 549 | |
| 550 | if self.setimode == False: |
| 551 | hbox.Add(vbox3, wx.ALIGN_RIGHT, 0)
|
| 552 | |
| 553 | vbox.Add(hbox, 0, wx.EXPAND)
|
| 554 | |
| 555 | self._build_subpanel(vbox)
|
| 556 | |
| 557 | self.lmst_timer = wx.PyTimer(self.lmst_timeout) |
| 558 | self.other_timer = wx.PyTimer(self.other_timeout) |
| 559 | |
| 560 | |
| 561 | def _build_subpanel(self, vbox_arg): |
| 562 | # build a secondary information panel (sometimes hidden)
|
| 563 | |
| 564 | # FIXME figure out how to have this be a subpanel that is always
|
| 565 | # created, but has its visibility controlled by foo.Show(True/False)
|
| 566 | |
| 567 | if not(self.show_debug_info): |
| 568 | return
|
| 569 | |
| 570 | panel = self.panel
|
| 571 | vbox = vbox_arg |
| 572 | myform = self.myform
|
| 573 | |
| 574 | #panel = wx.Panel(self.panel, -1)
|
| 575 | #vbox = wx.BoxSizer(wx.VERTICAL)
|
| 576 | |
| 577 | hbox = wx.BoxSizer(wx.HORIZONTAL) |
| 578 | hbox.Add((5,0), 0) |
| 579 | myform['decim'] = form.static_float_field(
|
| 580 | parent=panel, sizer=hbox, label="Decim")
|
| 581 | |
| 582 | hbox.Add((5,0), 1) |
| 583 | myform['USB BW'] = form.static_float_field(
|
| 584 | parent=panel, sizer=hbox, label="USB BW")
|
| 585 | |
| 586 | hbox.Add((5,0), 1) |
| 587 | myform['dbname'] = form.static_text_field(
|
| 588 | parent=panel, sizer=hbox) |
| 589 | |
| 590 | hbox.Add((5,0), 1) |
| 591 | myform['baseband'] = form.static_float_field(
|
| 592 | parent=panel, sizer=hbox, label="Analog BB")
|
| 593 | |
| 594 | hbox.Add((5,0), 1) |
| 595 | myform['ddc'] = form.static_float_field(
|
| 596 | parent=panel, sizer=hbox, label="DDC")
|
| 597 | |
| 598 | hbox.Add((5,0), 0) |
| 599 | vbox.Add(hbox, 0, wx.EXPAND)
|
| 600 | |
| 601 | |
| 602 | |
| 603 | def set_freq(self, target_freq): |
| 604 | """
|
| 605 | Set the center frequency we're interested in.
|
| 606 | |
| 607 | @param target_freq: frequency in Hz
|
| 608 | @rypte: bool
|
| 609 | |
| 610 | Tuning is a two step process. First we ask the front-end to
|
| 611 | tune as close to the desired frequency as it can. Then we use
|
| 612 | the result of that operation and our target_frequency to
|
| 613 | determine the value for the digital down converter.
|
| 614 | """ |
| 615 | #
|
| 616 | # Everything except BASIC_RX should support usrp.tune()
|
| 617 | #
|
| 618 | if not (self.cardtype == usrp_dbid.BASIC_RX): |
| 619 | r = usrp.tune(self.u, self.subdev[0].which(), self.subdev[0], target_freq) |
| 620 | r = usrp.tune(self.u, self.subdev[1].which(), self.subdev[1], target_freq) |
| 621 | else:
|
| 622 | r = self.u.set_rx_freq(0, target_freq) |
| 623 | f = self.u.rx_freq(0) |
| 624 | if abs(f-target_freq) > 2.0e3: |
| 625 | r = 0
|
| 626 | if r:
|
| 627 | self.myform['freq'].set_value(target_freq) # update displayed value |
| 628 | #
|
| 629 | # Make sure calibrator knows our target freq
|
| 630 | #
|
| 631 | |
| 632 | # Remember centerfreq---used for doppler calcs
|
| 633 | delta = self.centerfreq - target_freq
|
| 634 | self.centerfreq = target_freq
|
| 635 | self.observing -= delta
|
| 636 | self.scope.set_baseband_freq (self.observing) |
| 637 | |
| 638 | self.myform['baseband'].set_value(r.baseband_freq) |
| 639 | self.myform['ddc'].set_value(r.dxc_freq) |
| 640 | |
| 641 | if (self.use_notches): |
| 642 | self.compute_notch_taps(self.notches) |
| 643 | if self.dual_mode == False and self.interferometer == False: |
| 644 | self.notch_filt.set_taps(self.notch_taps) |
| 645 | else:
|
| 646 | self.notch_filt1.set_taps(self.notch_taps) |
| 647 | self.notch_filt2.set_taps(self.notch_taps) |
| 648 | |
| 649 | return True |
| 650 | |
| 651 | return False |
| 652 | |
| 653 | def set_decln(self, dec): |
| 654 | self.decln = dec
|
| 655 | self.myform['decln'].set_value(dec) # update displayed value |
| 656 | |
| 657 | def set_gain(self, gain): |
| 658 | self.myform['gain'].set_value(gain) # update displayed value |
| 659 | self.subdev[0].set_gain(gain) |
| 660 | self.subdev[1].set_gain(gain) |
| 661 | self.gain = gain
|
| 662 | |
| 663 | def set_averaging(self, avval): |
| 664 | self.myform['average'].set_value(avval) |
| 665 | self.scope.set_avg_alpha(1.0/(avval)) |
| 666 | self.scope.set_average(True) |
| 667 | self.avg_alpha = avval
|
| 668 | |
| 669 | def set_integration(self, integval): |
| 670 | if self.setimode == False: |
| 671 | self.integrator.set_taps(1.0/((integval)*(self.bw/2))) |
| 672 | self.myform['integration'].set_value(integval) |
| 673 | self.integ = integval
|
| 674 | |
| 675 | #
|
| 676 | # Timeout function
|
| 677 | # Used to update LMST display, as well as current
|
| 678 | # continuum value
|
| 679 | #
|
| 680 | # We also write external data-logging files here
|
| 681 | #
|
| 682 | def lmst_timeout(self): |
| 683 | self.locality.date = ephem.now()
|
| 684 | if self.setimode == False: |
| 685 | x = self.probe.level()
|
| 686 | sidtime = self.locality.sidereal_time()
|
| 687 | # LMST
|
| 688 | s = str(ephem.hours(sidtime)) + " " + self.sunstate |
| 689 | # Continuum detector value
|
| 690 | if self.setimode == False: |
| 691 | sx = "%7.4f" % x
|
| 692 | s = s + "\nDet: " + str(sx) |
| 693 | else:
|
| 694 | sx = "%2d" % self.hitcounter |
| 695 | s1 = "%2d" % self.s1hitcounter |
| 696 | s2 = "%2d" % self.s2hitcounter |
| 697 | sa = "%4.2f" % self.avgdelta |
| 698 | sy = "%3.1f-%3.1f" % (self.CHIRP_LOWER, self.CHIRP_UPPER) |
| 699 | s = s + "\nHits: " + str(sx) + "\nS1:" + str(s1) + " S2:" + str(s2) |
| 700 | s = s + "\nAv D: " + str(sa) + "\nCh lim: " + str(sy) |
| 701 | |
| 702 | self.myform['lmst_high'].set_value(s) |
| 703 | |
| 704 | #
|
| 705 | # Write data out to recording files
|
| 706 | #
|
| 707 | if self.setimode == False: |
| 708 | self.write_continuum_data(x,sidtime)
|
| 709 | self.write_spectral_data(self.fft_outbuf,sidtime) |
| 710 | |
| 711 | else:
|
| 712 | self.seti_analysis(self.fft_outbuf,sidtime) |
| 713 | now = time.time() |
| 714 | if ((self.scanning == True) and ((now - self.seti_then) > self.setifreq_timer)): |
| 715 | self.seti_then = now
|
| 716 | self.setifreq_current = self.setifreq_current + self.fft_input_rate |
| 717 | if (self.setifreq_current > self.setifreq_upper): |
| 718 | self.setifreq_current = self.setifreq_lower |
| 719 | self.set_freq(self.setifreq_current) |
| 720 | # Make sure we zero-out the hits array when changing
|
| 721 | # frequency.
|
| 722 | self.hits_array[:,:] = 0.0 |
| 723 | self.hit_intensities[:,:] = 0.0 |
| 724 | |
| 725 | def other_timeout(self): |
| 726 | if (self.switch_state == 0): |
| 727 | self.switch_state = 1 |
| 728 | |
| 729 | elif (self.switch_state == 1): |
| 730 | self.switch_state = 0 |
| 731 | |
| 732 | if (self.switch_state == 0): |
| 733 | self.mute.set_n(1) |
| 734 | self.cmute.set_n(int(1.0e9)) |
| 735 | |
| 736 | elif (self.switch_state == 1): |
| 737 | self.mute.set_n(int(1.0e9)) |
| 738 | self.cmute.set_n(1) |
| 739 | |
| 740 | if (self.ref_fifo != "@@@@"): |
| 741 | self.ref_fifo_file.write(str(self.switch_state)+"\n") |
| 742 | self.ref_fifo_file.flush()
|
| 743 | |
| 744 | self.avg_reference_value = self.cprobe.level() |
| 745 | |
| 746 | #
|
| 747 | # Set reference value
|
| 748 | #
|
| 749 | self.reference_level.set_k(-1.0 * (self.avg_reference_value/self.reference_divisor)) |
| 750 | |
| 751 | def fft_outfunc(self,data,l): |
| 752 | self.fft_outbuf=data
|
| 753 | |
| 754 | def write_continuum_data(self,data,sidtime): |
| 755 | |
| 756 | # Create localtime structure for producing filename
|
| 757 | foo = time.localtime() |
| 758 | pfx = self.prefix
|
| 759 | filenamestr = "%s/%04d%02d%02d%02d" % (pfx, foo.tm_year,
|
| 760 | foo.tm_mon, foo.tm_mday, foo.tm_hour) |
| 761 | |
| 762 | # Open the data file, appending
|
| 763 | continuum_file = open (filenamestr+".tpdat","a") |
| 764 | |
| 765 | flt = "%6.3f" % data
|
| 766 | inter = self.decln
|
| 767 | integ = self.integ
|
| 768 | fc = self.observing
|
| 769 | fc = fc / 1000000
|
| 770 | bw = self.bw
|
| 771 | bw = bw / 1000000
|
| 772 | ga = self.gain
|
| 773 | |
| 774 | now = time.time() |
| 775 | |
| 776 | #
|
| 777 | # If time to write full header info (saves storage this way)
|
| 778 | #
|
| 779 | if (now - self.continuum_then > 20): |
| 780 | self.sun.compute(self.locality) |
| 781 | enow = ephem.now() |
| 782 | sunset = self.locality.next_setting(self.sun) |
| 783 | sunrise = self.locality.next_rising(self.sun) |
| 784 | sun_insky = "Down"
|
| 785 | self.sunstate = "Dn" |
| 786 | if ((sunrise < enow) and (enow < sunset)): |
| 787 | sun_insky = "Up"
|
| 788 | self.sunstate = "Up" |
| 789 | self.continuum_then = now
|
| 790 | |
| 791 | continuum_file.write(str(ephem.hours(sidtime))+" "+flt+" Dn="+str(inter)+",") |
| 792 | continuum_file.write("Ti="+str(integ)+",Fc="+str(fc)+",Bw="+str(bw)) |
| 793 | continuum_file.write(",Ga="+str(ga)+",Sun="+str(sun_insky)+"\n") |
| 794 | else:
|
| 795 | continuum_file.write(str(ephem.hours(sidtime))+" "+flt+"\n") |
| 796 | |
| 797 | continuum_file.close() |
| 798 | return(data)
|
| 799 | |
| 800 | def write_spectral_data(self,data,sidtime): |
| 801 | |
| 802 | now = time.time() |
| 803 | delta = 10
|
| 804 | |
| 805 | # If time to write out spectral data
|
| 806 | # We don't write this out every time, in order to
|
| 807 | # save disk space. Since the spectral data are
|
| 808 | # typically heavily averaged, writing this data
|
| 809 | # "once in a while" is OK.
|
| 810 | #
|
| 811 | if (now - self.spectral_then >= delta): |
| 812 | self.spectral_then = now
|
| 813 | |
| 814 | # Get localtime structure to make filename from
|
| 815 | foo = time.localtime() |
| 816 | |
| 817 | pfx = self.prefix
|
| 818 | filenamestr = "%s/%04d%02d%02d%02d" % (pfx, foo.tm_year,
|
| 819 | foo.tm_mon, foo.tm_mday, foo.tm_hour) |
| 820 | |
| 821 | # Open the file
|
| 822 | spectral_file = open (filenamestr+".sdat","a") |
| 823 | |
| 824 | # Setup data fields to be written
|
| 825 | r = data |
| 826 | inter = self.decln
|
| 827 | fc = self.observing
|
| 828 | fc = fc / 1000000
|
| 829 | bw = self.bw
|
| 830 | bw = bw / 1000000
|
| 831 | av = self.avg_alpha
|
| 832 | |
| 833 | # Write those fields
|
| 834 | spectral_file.write("data:"+str(ephem.hours(sidtime))+" Dn="+str(inter)+",Fc="+str(fc)+",Bw="+str(bw)+",Av="+str(av)) |
| 835 | spectral_file.write (" [ ")
|
| 836 | for r in data: |
| 837 | spectral_file.write(" "+str(r)) |
| 838 | |
| 839 | spectral_file.write(" ]\n")
|
| 840 | spectral_file.close() |
| 841 | return(data)
|
| 842 | |
| 843 | return(data)
|
| 844 | |
| 845 | def seti_analysis(self,fftbuf,sidtime): |
| 846 | l = len(fftbuf)
|
| 847 | x = 0
|
| 848 | hits = [] |
| 849 | hit_intensities = [] |
| 850 | if self.seticounter < self.setitimer: |
| 851 | self.seticounter = self.seticounter + 1 |
| 852 | return
|
| 853 | else:
|
| 854 | self.seticounter = 0 |
| 855 | |
| 856 | # Run through FFT output buffer, computing standard deviation (Sigma)
|
| 857 | avg = 0
|
| 858 | # First compute average
|
| 859 | for i in range(0,l): |
| 860 | avg = avg + fftbuf[i] |
| 861 | avg = avg / l |
| 862 | |
| 863 | sigma = 0.0
|
| 864 | # Then compute standard deviation (Sigma)
|
| 865 | for i in range(0,l): |
| 866 | d = fftbuf[i] - avg |
| 867 | sigma = sigma + (d*d) |
| 868 | |
| 869 | sigma = Numeric.sqrt(sigma/l) |
| 870 | |
| 871 | #
|
| 872 | # Snarfle through the FFT output buffer again, looking for
|
| 873 | # outlying data points
|
| 874 | |
| 875 | start_f = self.observing - (self.fft_input_rate/2) |
| 876 | current_f = start_f |
| 877 | l = len(fftbuf)
|
| 878 | f_incr = self.fft_input_rate / l
|
| 879 | hit = -1
|
| 880 | |
| 881 | # -nyquist to DC
|
| 882 | for i in range(l/2,l): |
| 883 | #
|
| 884 | # If current FFT buffer has an item that exceeds the specified
|
| 885 | # sigma
|
| 886 | #
|
| 887 | if ((fftbuf[i] - avg) > (self.setik * sigma)): |
| 888 | hits.append(current_f) |
| 889 | hit_intensities.append(fftbuf[i]) |
| 890 | current_f = current_f + f_incr |
| 891 | |
| 892 | # DC to nyquist
|
| 893 | for i in range(0,l/2): |
| 894 | #
|
| 895 | # If current FFT buffer has an item that exceeds the specified
|
| 896 | # sigma
|
| 897 | #
|
| 898 | if ((fftbuf[i] - avg) > (self.setik * sigma)): |
| 899 | hits.append(current_f) |
| 900 | hit_intensities.append(fftbuf[i]) |
| 901 | current_f = current_f + f_incr |
| 902 | |
| 903 | # No hits
|
| 904 | if (len(hits) <= 0): |
| 905 | return
|
| 906 | |
| 907 | |
| 908 | #
|
| 909 | # OK, so we have some hits in the FFT buffer
|
| 910 | # They'll have a rather substantial gauntlet to run before
|
| 911 | # being declared a real "hit"
|
| 912 | #
|
| 913 | |
| 914 | # Update stats
|
| 915 | self.s1hitcounter = self.s1hitcounter + len(hits) |
| 916 | |
| 917 | # Weed out buffers with an excessive number of
|
| 918 | # signals above Sigma
|
| 919 | if (len(hits) > self.nhits): |
| 920 | return
|
| 921 | |
| 922 | |
| 923 | # Weed out FFT buffers with apparent multiple narrowband signals
|
| 924 | # separated significantly in frequency. This means that a
|
| 925 | # single signal spanning multiple bins is OK, but a buffer that
|
| 926 | # has multiple, apparently-separate, signals isn't OK.
|
| 927 | #
|
| 928 | last = hits[0]
|
| 929 | ns2 = 1
|
| 930 | for i in range(1,len(hits)): |
| 931 | if ((hits[i] - last) > (f_incr*3.0)): |
| 932 | return
|
| 933 | last = hits[i] |
| 934 | ns2 = ns2 + 1
|
| 935 | |
| 936 | self.s2hitcounter = self.s2hitcounter + ns2 |
| 937 | |
| 938 | #
|
| 939 | # Run through all available hit buffers, computing difference between
|
| 940 | # frequencies found there, if they're all within the chirp limits
|
| 941 | # declare a good hit
|
| 942 | #
|
| 943 | good_hit = False
|
| 944 | f_ds = Numeric.zeros(self.nhitlines, Numeric.Float64)
|
| 945 | avg_delta = 0
|
| 946 | k = 0
|
| 947 | for i in range(0,min(len(hits),len(self.hits_array[:,0]))): |
| 948 | f_ds[0] = abs(self.hits_array[i,0] - hits[i]) |
| 949 | for j in range(1,len(f_ds)): |
| 950 | f_ds[j] = abs(self.hits_array[i,j] - self.hits_array[i,0]) |
| 951 | avg_delta = avg_delta + f_ds[j] |
| 952 | k = k + 1
|
| 953 | |
| 954 | if (self.seti_isahit (f_ds)): |
| 955 | good_hit = True
|
| 956 | self.hitcounter = self.hitcounter + 1 |
| 957 | break
|
| 958 | |
| 959 | if (avg_delta/k < (self.seti_fft_bandwidth/2)): |
| 960 | self.avgdelta = avg_delta / k
|
| 961 | |
| 962 | # Save 'n shuffle hits
|
| 963 | # Old hit buffers percolate through the hit buffers
|
| 964 | # (there are self.nhitlines of these buffers)
|
| 965 | # and then drop off the end
|
| 966 | # A consequence is that while the nhitlines buffers are filling,
|
| 967 | # you can get some absurd values for self.avgdelta, because some
|
| 968 | # of the buffers are full of zeros
|
| 969 | for i in range(self.nhitlines,1): |
| 970 | self.hits_array[:,i] = self.hits_array[:,i-1] |
| 971 | self.hit_intensities[:,i] = self.hit_intensities[:,i-1] |
| 972 | |
| 973 | for i in range(0,len(hits)): |
| 974 | self.hits_array[i,0] = hits[i] |
| 975 | self.hit_intensities[i,0] = hit_intensities[i] |
| 976 | |
| 977 | # Finally, write the hits/intensities buffer
|
| 978 | if (good_hit):
|
| 979 | self.write_hits(sidtime)
|
| 980 | |
| 981 | return
|
| 982 | |
| 983 | def seti_isahit(self,fdiffs): |
| 984 | truecount = 0
|
| 985 | |
| 986 | for i in range(0,len(fdiffs)): |
| 987 | if (fdiffs[i] >= self.CHIRP_LOWER and fdiffs[i] <= self.CHIRP_UPPER): |
| 988 | truecount = truecount + 1
|
| 989 | |
| 990 | if truecount == len(fdiffs): |
| 991 | return (True) |
| 992 | else:
|
| 993 | return (False) |
| 994 | |
| 995 | def write_hits(self,sidtime): |
| 996 | # Create localtime structure for producing filename
|
| 997 | foo = time.localtime() |
| 998 | pfx = self.prefix
|
| 999 | filenamestr = "%s/%04d%02d%02d%02d" % (pfx, foo.tm_year,
|
| 1000 | foo.tm_mon, foo.tm_mday, foo.tm_hour) |
| 1001 | |
| 1002 | # Open the data file, appending
|
| 1003 | hits_file = open (filenamestr+".seti","a") |
| 1004 | |
| 1005 | # Write sidtime first
|
| 1006 | hits_file.write(str(ephem.hours(sidtime))+" "+str(self.decln)+" ") |
| 1007 | |
| 1008 | #
|
| 1009 | # Then write the hits/hit intensities buffers with enough
|
| 1010 | # "syntax" to allow parsing by external (not yet written!)
|
| 1011 | # "stuff".
|
| 1012 | #
|
| 1013 | for i in range(0,self.nhitlines): |
| 1014 | hits_file.write(" ")
|
| 1015 | for j in range(0,self.nhits): |
| 1016 | hits_file.write(str(self.hits_array[j,i])+":") |
| 1017 | hits_file.write(str(self.hit_intensities[j,i])+",") |
| 1018 | hits_file.write("\n")
|
| 1019 | hits_file.close() |
| 1020 | return
|
| 1021 | |
| 1022 | def xydfunc(self,func,xyv): |
| 1023 | if self.setimode == True: |
| 1024 | return
|
| 1025 | magn = int(Numeric.log10(self.observing)) |
| 1026 | if (magn == 6 or magn == 7 or magn == 8): |
| 1027 | magn = 6
|
| 1028 | dfreq = xyv[0] * pow(10.0,magn) |
| 1029 | if func == 0: |
| 1030 | ratio = self.observing / dfreq
|
| 1031 | vs = 1.0 - ratio
|
| 1032 | vs *= 299792.0
|
| 1033 | if magn >= 9: |
| 1034 | xhz = "Ghz"
|
| 1035 | elif magn >= 6: |
| 1036 | xhz = "Mhz"
|
| 1037 | elif magn <= 5: |
| 1038 | xhz = "Khz"
|
| 1039 | s = "%.6f%s\n%.3fdB" % (xyv[0], xhz, xyv[1]) |
| 1040 | s2 = "\n%.3fkm/s" % vs
|
| 1041 | self.myform['spec_data'].set_value(s+s2) |
| 1042 | else:
|
| 1043 | tmpnotches = Numeric.zeros(self.NOTCH_TAPS,Numeric.Float64)
|
| 1044 | delfreq = -1
|
| 1045 | if self.use_notches == True: |
| 1046 | for i in range(0,len(self.notches)): |
| 1047 | if abs(self.notches[i] - dfreq) <= (self.bw/self.NOTCH_TAPS): |
| 1048 | delfreq = i |
| 1049 | break
|
| 1050 | j = 0
|
| 1051 | for i in range(0,len(self.notches)): |
| 1052 | if (i != delfreq):
|
| 1053 | tmpnotches[j] = self.notches[i]
|
| 1054 | j = j + 1
|
| 1055 | if (delfreq == -1): |
| 1056 | for i in range(0,len(tmpnotches)): |
| 1057 | if (int(tmpnotches[i]) == 0): |
| 1058 | tmpnotches[i] = dfreq |
| 1059 | break
|
| 1060 | self.notches = tmpnotches
|
| 1061 | self.compute_notch_taps(self.notches) |
| 1062 | if self.dual_mode == False and self.interferometer == False: |
| 1063 | self.notch_filt.set_taps(self.notch_taps) |
| 1064 | else:
|
| 1065 | self.notch_filt1.set_taps(self.notch_taps) |
| 1066 | self.notch_filt2.set_taps(self.notch_taps) |
| 1067 | |
| 1068 | def xydfunc_waterfall(self,pos): |
| 1069 | lower = self.observing - (self.seti_fft_bandwidth / 2) |
| 1070 | upper = self.observing + (self.seti_fft_bandwidth / 2) |
| 1071 | binwidth = self.seti_fft_bandwidth / 1024 |
| 1072 | s = "%.6fMHz" % ((lower + (pos.x*binwidth)) / 1.0e6) |
| 1073 | self.myform['spec_data'].set_value(s) |
| 1074 | |
| 1075 | def toggle_cal(self): |
| 1076 | if (self.calstate == True): |
| 1077 | self.calstate = False |
| 1078 | self.u.write_io(0,0,(1<<15)) |
| 1079 | self.calibrator.SetLabel("Calibration Source: Off") |
| 1080 | else:
|
| 1081 | self.calstate = True |
| 1082 | self.u.write_io(0,(1<<15),(1<<15)) |
| 1083 | self.calibrator.SetLabel("Calibration Source: On") |
| 1084 | |
| 1085 | def toggle_annotation(self): |
| 1086 | if (self.annotate_state == True): |
| 1087 | self.annotate_state = False |
| 1088 | self.annotation.SetLabel("Annotation: Off") |
| 1089 | else:
|
| 1090 | self.annotate_state = True |
| 1091 | self.annotation.SetLabel("Annotation: On") |
| 1092 | #
|
| 1093 | # Turn scanning on/off
|
| 1094 | # Called-back by "Recording" button
|
| 1095 | #
|
| 1096 | def toggle_scanning(self): |
| 1097 | # Current scanning? Flip state
|
| 1098 | if (self.scanning == True): |
| 1099 | self.scanning = False |
| 1100 | self.scan_control.SetLabel("Scan: Off") |
| 1101 | # Not scanning
|
| 1102 | else:
|
| 1103 | self.scanning = True |
| 1104 | self.scan_control.SetLabel("Scan: On ") |
| 1105 | |
| 1106 | def set_pd_offset(self,offs): |
| 1107 | self.myform['offset'].set_value(offs) |
| 1108 | self.calib_offset=offs
|
| 1109 | x = self.calib_coeff / 100.0 |
| 1110 | self.cal_offs.set_k(offs*(x*8000)) |
| 1111 | |
| 1112 | def set_pd_gain(self,gain): |
| 1113 | self.myform['dcgain'].set_value(gain) |
| 1114 | self.cal_mult.set_k(gain*0.01) |
| 1115 | self.calib_coeff = gain
|
| 1116 | x = gain/100.0
|
| 1117 | self.cal_offs.set_k(self.calib_offset*(x*8000)) |
| 1118 | |
| 1119 | def compute_notch_taps(self,notchlist): |
| 1120 | tmptaps = Numeric.zeros(self.NOTCH_TAPS,Numeric.Complex64)
|
| 1121 | binwidth = self.bw / self.NOTCH_TAPS |
| 1122 | |
| 1123 | for i in range(0,self.NOTCH_TAPS): |
| 1124 | tmptaps[i] = complex(1.0,0.0) |
| 1125 | |
| 1126 | for i in notchlist: |
| 1127 | diff = i - self.observing
|
| 1128 | if int(i) == 0: |
| 1129 | break
|
| 1130 | if (diff > 0): |
| 1131 | idx = diff / binwidth |
| 1132 | idx = round(idx)
|
| 1133 | idx = int(idx)
|
| 1134 | if (idx < 0 or idx > (self.NOTCH_TAPS/2)): |
| 1135 | break
|
| 1136 | tmptaps[idx] = complex(0.0, 0.0) |
| 1137 | |
| 1138 | if (diff < 0): |
| 1139 | idx = -diff / binwidth |
| 1140 | idx = round(idx)
|
| 1141 | idx = (self.NOTCH_TAPS/2) - idx |
| 1142 | idx = int(idx+(self.NOTCH_TAPS/2)) |
| 1143 | if (idx < 0 or idx > (self.NOTCH_TAPS)): |
| 1144 | break
|
| 1145 | tmptaps[idx] = complex(0.0, 0.0) |
| 1146 | |
| 1147 | self.notch_taps = numpy.fft.ifft(tmptaps)
|
| 1148 | |
| 1149 | #
|
| 1150 | # Setup common pieces of radiometer mode
|
| 1151 | #
|
| 1152 | def setup_radiometer_common(self,n): |
| 1153 | # The IIR integration filter for post-detection
|
| 1154 | self.integrator = gr.single_pole_iir_filter_ff(1.0) |
| 1155 | self.integrator.set_taps (1.0/self.bw) |
| 1156 | |
| 1157 | if (self.use_notches == True): |
| 1158 | self.compute_notch_taps(self.notches) |
| 1159 | if (n == 2): |
| 1160 | self.notch_filt1 = gr.fft_filter_ccc(1, self.notch_taps) |
| 1161 | self.notch_filt2 = gr.fft_filter_ccc(1, self.notch_taps) |
| 1162 | else:
|
| 1163 | self.notch_filt = gr.fft_filter_ccc(1, self.notch_taps) |
| 1164 | |
| 1165 | |
| 1166 | # Signal probe
|
| 1167 | self.probe = gr.probe_signal_f()
|
| 1168 | |
| 1169 | #
|
| 1170 | # Continuum calibration stuff
|
| 1171 | #
|
| 1172 | x = self.calib_coeff/100.0 |
| 1173 | self.cal_mult = gr.multiply_const_ff(self.calib_coeff/100.0) |
| 1174 | self.cal_offs = gr.add_const_ff(self.calib_offset*(x*8000)) |
| 1175 | |
| 1176 | #
|
| 1177 | # Mega decimator after IIR filter
|
| 1178 | #
|
| 1179 | if (self.switch_mode == False): |
| 1180 | self.keepn = gr.keep_one_in_n(gr.sizeof_float, self.bw) |
| 1181 | else:
|
| 1182 | self.keepn = gr.keep_one_in_n(gr.sizeof_float, int(self.bw/2)) |
| 1183 | |
| 1184 | #
|
| 1185 | # For the Dicke-switching scheme
|
| 1186 | #
|
| 1187 | #self.switch = gr.multiply_const_ff(1.0)
|
| 1188 | |
| 1189 | #
|
| 1190 | if (self.switch_mode == True): |
| 1191 | self.vector = gr.vector_sink_f()
|
| 1192 | self.swkeep = gr.keep_one_in_n(gr.sizeof_float, int(self.bw/3)) |
| 1193 | self.mute = gr.keep_one_in_n(gr.sizeof_float, 1) |
| 1194 | self.cmute = gr.keep_one_in_n(gr.sizeof_float, int(1.0e9)) |
| 1195 | self.cintegrator = gr.single_pole_iir_filter_ff(1.0/(self.bw/2)) |
| 1196 | self.cprobe = gr.probe_signal_f()
|
| 1197 | else:
|
| 1198 | self.mute = gr.multiply_const_ff(1.0) |
| 1199 | |
| 1200 | |
| 1201 | self.avg_reference_value = 0.0 |
| 1202 | self.reference_level = gr.add_const_ff(0.0) |
| 1203 | |
| 1204 | #
|
| 1205 | # Setup ordinary single-channel radiometer mode
|
| 1206 | #
|
| 1207 | def setup_normal(self, setimode): |
| 1208 | |
| 1209 | self.setup_radiometer_common(1) |
| 1210 | |
| 1211 | self.head = self.u |
| 1212 | if (self.use_notches == True): |
| 1213 | self.shead = self.notch_filt |
| 1214 | else:
|
| 1215 | self.shead = self.u |
| 1216 | |
| 1217 | if setimode == False: |
| 1218 | |
| 1219 | self.detector = gr.complex_to_mag_squared()
|
| 1220 | self.connect(self.shead, self.scope) |
| 1221 | |
| 1222 | if (self.use_notches == False): |
| 1223 | self.connect(self.head, self.detector, self.mute, self.reference_level, |
| 1224 | self.integrator, self.keepn, self.cal_mult, self.cal_offs, self.chart) |
| 1225 | else:
|
| 1226 | self.connect(self.head, self.notch_filt, self.detector, self.mute, self.reference_level, |
| 1227 | self.integrator, self.keepn, self.cal_mult, self.cal_offs, self.chart) |
| 1228 | |
| 1229 | self.connect(self.cal_offs, self.probe) |
| 1230 | |
| 1231 | #
|
| 1232 | # Add a side-chain detector chain, with a different integrator, for sampling
|
| 1233 | # The reference channel data
|
| 1234 | # This is used to derive the offset value for self.reference_level, used above
|
| 1235 | #
|
| 1236 | if (self.switch_mode == True): |
| 1237 | self.connect(self.detector, self.cmute, self.cintegrator, self.swkeep, self.cprobe) |
| 1238 | |
| 1239 | return
|
| 1240 | |
| 1241 | #
|
| 1242 | # Setup dual-channel (two antenna, usual orthogonal polarity probes in the same waveguide)
|
| 1243 | #
|
| 1244 | def setup_dual(self, setimode,notches): |
| 1245 | |
| 1246 | self.setup_radiometer_common(2) |
| 1247 | |
| 1248 | self.di = gr.deinterleave(gr.sizeof_gr_complex)
|
| 1249 | self.addchans = gr.add_cc ()
|
| 1250 | self.detector = gr.add_ff ()
|
| 1251 | self.h_power = gr.complex_to_mag_squared()
|
| 1252 | self.v_power = gr.complex_to_mag_squared()
|
| 1253 | self.connect (self.u, self.di) |
| 1254 | |
| 1255 | if (self.use_notches == True): |
| 1256 | self.connect((self.di, 0), self.notch_filt1, (self.addchans, 0)) |
| 1257 | self.connect((self.di, 1), self.notch_filt2, (self.addchans, 1)) |
| 1258 | else:
|
| 1259 | #
|
| 1260 | # For spectral, adding the two channels works, assuming no gross
|
| 1261 | # phase or amplitude error
|
| 1262 | self.connect ((self.di, 0), (self.addchans, 0)) |
| 1263 | self.connect ((self.di, 1), (self.addchans, 1)) |
| 1264 | |
| 1265 | #
|
| 1266 | # Connect heads of spectral and total-power chains
|
| 1267 | #
|
| 1268 | if (self.use_notches == False): |
| 1269 | self.head = self.di |
| 1270 | else:
|
| 1271 | self.head = (self.notch_filt1, self.notch_filt2) |
| 1272 | |
| 1273 | self.shead = self.addchans |
| 1274 | |
| 1275 | if (setimode == False): |
| 1276 | #
|
| 1277 | # For dual-polarization mode, we compute the sum of the
|
| 1278 | # powers on each channel, after they've been detected
|
| 1279 | #
|
| 1280 | self.detector = gr.add_ff()
|
| 1281 | |
| 1282 | #
|
| 1283 | # In dual-polarization mode, we compute things a little differently
|
| 1284 | # In effect, we have two radiometer chains, terminating in an adder
|
| 1285 | #
|
| 1286 | if self.use_notches == True: |
| 1287 | self.connect(self.notch_filt1, self.h_power) |
| 1288 | self.connect(self.notch_filt2, self.v_power) |
| 1289 | else:
|
| 1290 | self.connect((self.head, 0), self.h_power) |
| 1291 | self.connect((self.head, 1), self.v_power) |
| 1292 | self.connect(self.h_power, (self.detector, 0)) |
| 1293 | self.connect(self.v_power, (self.detector, 1)) |
| 1294 | self.connect(self.detector, self.mute, self.reference_level, |
| 1295 | self.integrator, self.keepn, self.cal_mult, self.cal_offs, self.chart) |
| 1296 | self.connect(self.cal_offs, self.probe) |
| 1297 | self.connect(self.shead, self.scope) |
| 1298 | |
| 1299 | #
|
| 1300 | # Add a side-chain detector chain, with a different integrator, for sampling
|
| 1301 | # The reference channel data
|
| 1302 | # This is used to derive the offset value for self.reference_level, used above
|
| 1303 | #
|
| 1304 | if (self.switch_mode == True): |
| 1305 | self.connect(self.detector, self.cmute, self.cintegrator, self.swkeep, self.cprobe) |
| 1306 | return
|
| 1307 | |
| 1308 | #
|
| 1309 | # Setup correlating interferometer mode
|
| 1310 | #
|
| 1311 | def setup_interferometer(self, setimode): |
| 1312 | self.setup_radiometer_common(2) |
| 1313 | |
| 1314 | self.di = gr.deinterleave(gr.sizeof_gr_complex)
|
| 1315 | self.connect (self.u, self.di) |
| 1316 | self.corr = gr.multiply_cc()
|
| 1317 | self.c2f = gr.complex_to_float()
|
| 1318 | |
| 1319 | self.shead = (self.di, 0) |
| 1320 | |
| 1321 | # Channel 0 to multiply port 0
|
| 1322 | # Channel 1 to multiply port 1
|
| 1323 | if (self.use_notches == False): |
| 1324 | self.connect((self.di, 0), (self.corr, 0)) |
| 1325 | self.connect((self.di, 1), (self.corr, 1)) |
| 1326 | else:
|
| 1327 | self.connect((self.di, 0), self.notch_filt1, (self.corr, 0)) |
| 1328 | self.connect((self.di, 1), self.notch_filt2, (self.corr, 0)) |
| 1329 | |
| 1330 | #
|
| 1331 | # Multiplier (correlator) to complex-to-float, followed by integrator, etc
|
| 1332 | #
|
| 1333 | self.connect(self.corr, self.c2f, self.integrator, self.keepn, self.cal_mult, self.cal_offs, self.chart) |
| 1334 | |
| 1335 | #
|
| 1336 | # FFT scope gets only 1 channel
|
| 1337 | # FIX THIS, by cross-correlating the *outputs* of two different FFTs, then display
|
| 1338 | # Funky!
|
| 1339 | #
|
| 1340 | self.connect(self.shead, self.scope) |
| 1341 | |
| 1342 | #
|
| 1343 | # Output of correlator/integrator chain to probe
|
| 1344 | #
|
| 1345 | self.connect(self.cal_offs, self.probe) |
| 1346 | |
| 1347 | return
|
| 1348 | |
| 1349 | #
|
| 1350 | # Setup SETI mode
|
| 1351 | #
|
| 1352 | def setup_seti(self): |
| 1353 | self.connect (self.shead, self.fft_bandpass, self.scope) |
| 1354 | return
|
| 1355 | |
| 1356 | |
| 1357 | |
| 1358 | def main (): |
| 1359 | app = stdgui2.stdapp(app_flow_graph, "RADIO ASTRONOMY SPECTRAL/CONTINUUM RECEIVER: $Revision$", nstatus=1) |
| 1360 | app.MainLoop() |
| 1361 | |
| 1362 | if __name__ == '__main__': |
| 1363 | main () |