tofil.py 5.27 KB
Newer Older
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
1
2
3
4
5
6
7
#!/usr/bin/env python2
from __future__ import print_function

import numpy as np
import filterbank
import sigproc
from astropy.time import Time
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
8
from astropy import units as u
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
9
10

import os.path
11
import os
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
12
13
14
15
16
import sys
import re

import argparse

17
18
19
20
21
22
23
telescope_found = False
try:
    from telescope import telescope
    telescope_found = True
except:
    print("Could not connect to console, no RA/DEC recorded")

24
def create_filterbank(infile=None, outfile=None, date=None, source=None):
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
25
26
27
28
29
30
31
32
33
34
35
36
37
    """
    Create a filterbank file from the CAMRAS-format, which is generated
    by pulsar_record.

    infile:   name of the file generated by pulsar_record
    outfile: name of the filterbank file to be created; if not given, a name
              will be created from the infile
    date:     the date to be used; default is the file creation date of infile
              if given, the date must be in a format the astropy.time understands
    source:   name of the source to be used; default is guessed from filename
    """

    if not outfile:
38
39
        if not infile:
            raise RuntimeError("Specify at least one of inputfile or outputfile")
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
40
41
42
        outfile = re.split('[+-\.]', infile)[0] + ".fil"

    if not source:
43
44
        if not infile:
            raise RuntimeError("If no input file is given, you must specify source")
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
45
46
47
48
49
50
51
        match = re.match("[BJ][0-9]*[+-][0-9]*", infile)
        if match:
            source = "PSR " + match.group(0)
        else:
            raise ValueError("Could not guess source name from filename")

    if not date:
52
53
54
55
        if infile:
            date = Time(os.path.getctime(infile), format='unix')
        else:
            date = Time.now()
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
56
57
58
    else:
        date = Time(date, format='isot')

59
60
61
62
63
    az_start = 0.
    za_start = 0.
    src_raj = 0.
    src_dej = 0.
    if telescope_found:
64
65
66
67
        dt = telescope(consoleHost='console')
        az_start = (180*u.deg + dt.getAzEl(waitForUpdate=True)[0]).to(u.deg).value
        print(az_start)
        za_start = (90*u.deg - dt.getAzEl()[1]).to(u.deg).value
68
69
70
        src_raj = dt.getRaDec().ra.to_string(unit=u.hour, sep='')
        src_dej = dt.getRaDec().dec.to_string(unit=u.degree, sep='')

71
72
    print("Input file: {}".format(infile))
    print("Output file: {}".format(outfile))
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
73
    print("Source name: {}".format(source))
74
    print("Observation time: {}".format(date.isot))
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
75
76

    date = date.mjd
77
78
79
80
81

    if infile:
        rawdatafile = infile
    else:
        rawdatafile = ""
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
82
83
84
85
86
87
 
    # Generate filterbank header
    fil_header = {}
    fil_header["telescope_id"]  = sigproc.telescope_ids["Effelsberg"]
    fil_header["machine_id"]    = sigproc.machine_ids["FAKE"]
    fil_header["data_type"]     = 1
88
    fil_header["rawdatafile"]   = rawdatafile
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
89
90
91
    fil_header["source_name"]   = source
    fil_header["barycentric"]   = 0
    fil_header["pulsarcentric"] = 0
92
93
94
95
    fil_header["az_start"]      = az_start
    fil_header["za_start"]      = za_start
    fil_header["src_raj"]       = src_raj
    fil_header["src_dej"]       = src_dej
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
96
97
    fil_header["tstart"]        = date
    fil_header["tsamp"]         = 0.00046811428571414528 # 1.0/35e6*256 chans*64 decimation
98
    fil_header["nbits"]         = 16
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
99
100
101
102
103
104
    fil_header["fch1"]          = 441.4
    fil_header["foff"]          = -35.0/256.0
    fil_header["nchans"]        = 256
    fil_header["nifs"]          = 1

    # Write header
105
    out = filterbank.create_filterbank_file(outfile, fil_header, nbits=16)
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
106

107
108
109
110
    if not infile:
        out.close()
        return

Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
111
    # Read file
112
    f = open(infile)
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
113

114
115
    filesize = os.stat(infile).st_size
    offset = 0
116
    CHUNKSIZE = 30000
117
    while offset < filesize:
118
        # Loop through file in chunks of 256*10000 numbers
119
        # i.e. 256*CHUNKSIZE*4 bytes (10Mb) to avoid memory problems
120
        f.seek(offset, os.SEEK_SET)
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
121

122
        data = np.fliplr(np.fromfile(f, dtype='>u4', count=256*CHUNKSIZE).reshape(-1,256))
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
123

124
        # Check that no packets were dropped
125
126
127
        if np.sum(data[1:,-1] - data[:-1,-1] - 1) > 0.1:
            print("Packets were dropped:", np.count_nonzero(data[:-1,-1] - data[1:,-1] - 1),
                                           np.sum(data[1:,-1] - data[:-1,-1] - 1))
128
129

        # Write data
130
        out.append_spectra(data.astype('u2'))
131

132
        offset += 256*CHUNKSIZE*4
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
133
134
135
136
137
138

    # Close file
    out.close()

if __name__ == "__main__":
    parser = argparse.ArgumentParser(description="Convert a camras pulsar file (output of pulsar_record) to filterbank format")
139
    parser.add_argument("-i", "--infile", help="Input CAMRAS pulsar file. If empty, create just the header.")
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
140
    parser.add_argument("-o", "--outfile", help="Name of output file, e.g. 'B0329.fil' (default: guessed from input filename)", default=None)
141
    parser.add_argument("-d", "--date", help="Start date/time of observation, in 'isot' format (defaults to creation date of file, or now if no infile given)", default=None)
Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
142
143
144
145
    parser.add_argument("-s", "--source", help="Name of the source, e.g. 'PSR B0329+54' (default: guessed from input filename)", default=None)

    args = parser.parse_args()

Tammo Jan Dijkema's avatar
Tammo Jan Dijkema committed
146
147
148
149
150
    if args.infile is None and args.outfile is None:
        parser.error("One of --infile or --outfile is required")
    if args.infile is None and args.source is None:
        parser.error("Without infile, --source must be specified")
 
151
    create_filterbank(infile=args.infile, outfile=args.outfile, date=args.date, source=args.source)