MARK 5B

The Mark 5B format is the output format of the Mark 5B disk-based VLBI data system. It is described in its design specifications.

File Structure

Each data frame consists of a header consisting of four 32-bit words (16 bytes) followed by a payload of 2500 32-bit words (10000 bytes). The header contains a sync word, frame number, and timestamp (accurate to 1 ms), as well as user-specified data; see Sec. 1 of the design specifications for details. The payload supports \(2^n\) bit streams, for \(0 \leq n \leq 5\), and the first sample of each stream corresponds precisely to the header time. elementary samples may be 1 or 2 bits in size, with the latter being stored in two successive bit streams. The number of channels is equal to the number of bit-streams divided by the number of bits per elementary sample (Baseband currently only supports files where all bit-streams are active). Files begin at a header (unlike for Mark 4), and an integer number of frames fit within 1 second.

The Mark 5B system also outputs files with the active bit-stream mask, number of frames per second, and observational metadata (Sec. 1.3 of the design specifications). Baseband does not yet use these files, and instead requires the user specify, for example, the sample rate.

Usage

This section covers reading and writing Mark 5B files with Baseband; general usage can be found under the Using Baseband section. For situations in which one is unsure of a file’s format, Baseband features the general baseband.open and baseband.file_info functions, which are also discussed in Using Baseband. The examples below use the small sample file baseband/data/sample.m5b, and the numpy, astropy.units, astropy.time.Time, and baseband.mark5b modules:

>>> import numpy as np
>>> import astropy.units as u
>>> from astropy.time import Time
>>> from baseband import mark5b
>>> from baseband.data import SAMPLE_MARK5B

Opening a Mark 5B file with open in binary mode provides a normal file reader extended with methods to read a Mark5BFrame. The number of channels, kiloday (thousands of MJD) and number of bits per sample must all be passed when using read_frame:

>>> fb = mark5b.open(SAMPLE_MARK5B, 'rb', kday=56000, nchan=8)
>>> frame = fb.read_frame()
>>> frame.shape
(5000, 8)
>>> fb.close()

Our sample file has 2-bit component samples, which is also the default for read_frame, so it does not need to be passed. Also, we may pass a reference Time object within 500 days of the observation start time to ref_time, rather than kday.

Opening as a stream wraps the low-level routines such that reading and writing is in units of samples. It also provides access to header information. Here, we also must provide nchan, sample_rate, and ref_time or kday:

>>> fh = mark5b.open(SAMPLE_MARK5B, 'rs', sample_rate=32*u.MHz, nchan=8,
...                  ref_time=Time('2014-06-13 12:00:00'))
>>> fh
<Mark5BStreamReader name=... offset=0
    sample_rate=32.0 MHz, samples_per_frame=5000,
    sample_shape=SampleShape(nchan=8), bps=2,
    start_time=2014-06-13T05:30:01.000000000>
>>> header0 = fh.header0    # To be used for writing, below.
>>> d = fh.read(10000)
>>> d.shape
(10000, 8)
>>> d[0, :3]    
array([-3.316505, -1.      ,  1.      ], dtype=float32)
>>> fh.close()

When writing to file, we again need to pass in sample_rate and nchan, though time can either be passed explicitly or inferred from the header:

>>> fw = mark5b.open('test.m5b', 'ws', header0=header0,
...                  sample_rate=32*u.MHz, nchan=8)
>>> fw.write(d)
>>> fw.close()
>>> fh = mark5b.open('test.m5b', 'rs', sample_rate=32*u.MHz,
...                  kday=57000, nchan=8)
>>> np.all(fh.read() == d)
True
>>> fh.close()

Reference/API

baseband.mark5b Package

Mark5B VLBI data reader.

Code inspired by Walter Brisken’s mark5access. See https://github.com/difx/mark5access.

Also, for the Mark5B design, see https://www.haystack.mit.edu/tech/vlbi/mark5/mark5_memos/019.pdf

Functions

info(name, **kwargs)

Collect Mark5B file information.

open(name[, mode])

Open Mark5B file(s) for reading or writing.

Classes

Mark5BFrame(header, payload[, valid, verify])

Representation of a Mark 5B frame, consisting of a header and payload.

Mark5BHeader(words[, kday, ref_time, verify])

Decoder/encoder of a Mark5B Frame Header.

Mark5BPayload(words, *[, header, ...])

Container for decoding and encoding VDIF payloads.

Class Inheritance Diagram

Inheritance diagram of baseband.mark5b.frame.Mark5BFrame, baseband.mark5b.header.Mark5BHeader, baseband.mark5b.payload.Mark5BPayload

baseband.mark5b.header Module

Definitions for VLBI Mark5B Headers.

Implements a Mark5BHeader class used to store header words, and decode/encode the information therein.

For the specification, see https://www.haystack.mit.edu/tech/vlbi/mark5/docs/Mark%205B%20users%20manual.pdf

Classes

Mark5BHeader(words[, kday, ref_time, verify])

Decoder/encoder of a Mark5B Frame Header.

Variables

CRC16

CRC polynomial used for Mark 5B Headers, as a check on the time code.

crc16

Cyclic Redundancy Check.

Class Inheritance Diagram

Inheritance diagram of baseband.mark5b.header.Mark5BHeader

baseband.mark5b.payload Module

Definitions for VLBI Mark 5B payloads.

Implements a Mark5BPayload class used to store payload words, and decode to or encode from a data array.

For the specification, see https://www.haystack.mit.edu/tech/vlbi/mark5/docs/Mark%205B%20users%20manual.pdf

Functions

init_luts()

Set up the look-up tables for levels as a function of input byte.

decode_1bit(words)

decode_2bit(words)

encode_1bit(values)

Encodes values using 1 bit per sample, packing the result into bytes.

encode_2bit(values)

Generic encoder for data stored using two bits.

Classes

Mark5BPayload(words, *[, header, ...])

Container for decoding and encoding VDIF payloads.

Class Inheritance Diagram

Inheritance diagram of baseband.mark5b.payload.Mark5BPayload

baseband.mark5b.frame Module

Definitions for VLBI Mark 5B frames.

Implements a Mark5BFrame class that can be used to hold a header and a payload, providing access to the values encoded in both.

For the specification, see https://www.haystack.mit.edu/tech/vlbi/mark5/docs/Mark%205B%20users%20manual.pdf

Classes

Mark5BFrame(header, payload[, valid, verify])

Representation of a Mark 5B frame, consisting of a header and payload.

Class Inheritance Diagram

Inheritance diagram of baseband.mark5b.frame.Mark5BFrame

baseband.mark5b.file_info Module

The Mark5BFileReaderInfo property.

Includes information about what is needed to calcuate times.

Classes

Mark5BFileReaderInfo([parent])

Class Inheritance Diagram

Inheritance diagram of baseband.mark5b.file_info.Mark5BFileReaderInfo

baseband.mark5b.base Module

Functions

open(name[, mode])

Open Mark5B file(s) for reading or writing.

info(name, **kwargs)

Collect Mark5B file information.

Classes

Mark5BFileReader(fh_raw[, kday, ref_time, ...])

Simple reader for Mark 5B files.

Mark5BFileWriter(fh_raw)

Simple writer for Mark 5B files.

Mark5BStreamBase()

Provides sample shape maker and fast time and index getting/setting.

Mark5BStreamReader(fh_raw[, sample_rate, ...])

VLBI Mark 5B format reader.

Mark5BStreamWriter(fh_raw[, header0, ...])

VLBI Mark 5B format writer.

Class Inheritance Diagram

Inheritance diagram of baseband.mark5b.base.Mark5BFileReader, baseband.mark5b.base.Mark5BFileWriter, baseband.mark5b.base.Mark5BStreamBase, baseband.mark5b.base.Mark5BStreamReader, baseband.mark5b.base.Mark5BStreamWriter