Comparison with Cotter
- Dev Null
Owned by Dev Null
The visibility values produced by Birli are thoroughly tested against the equivalent values from Cotter on each commit. Any change which causes these tests to fail is rejected.
Because of a number of small differences between Cotter and Birli, these tests vary slightly in error tolerance, depending which corrections are applied in the test data.
- Corrections in Birli are done using double precision float values, while Cotter only uses single precision.
- Although this results in a slight speed decrease for these particular operations, Birli more than makes up for this in other areas and is significantly faster overall.
- Birli uses a more accurate and higher resolution PFB gain curve which is provided at the highest frequency resolution that the MWAX correlator can use.
- You can use the PFB gains from Cotter with the --pfb-gains cotter command line flag
- Birli uses a more accurate (double precision) value for the array position.
- You can use the --emulate-cotter flag to use the same array position value as Cotter.
- Birli pays close attention to the timestamp values provided in the gpubox files and validates that these values match what is specified in the metafits.
- Cotter only reads the first timestamp value of the first gpubox file, which can result in it processing corrupted visibility files without warning the user.
- When Birli encounters a conflict between the raw files and the metafits, the user may be required to update the values in the metafits or raw files to resolve this conflict.
- On the flip side, there are some cases where the legacy correlator could produce a file that can be preprocessed with Cotter but not Birli.
- Some examples of legacy correlator issues are documented here
- Some Measurement Set and UVFits files created by Birli, Hyperdrive or Marlu can use the DUT1 timescale
- The --ignore-dut1 flag can be used to emulate Cotter's behaviour and write timesteps the the UTC timescale.
- it's best to read the timescale (DUT1 or UTC) specified in the file metadata.
Testing your own observation
One of the most efficient ways to validate that Cotter and Birli produce the same visibilities is using taql, which can query a measurement significantly faster than Python.
This bash script preprocesses the same observation with Cotter and Birli, giving the best direct comparison assumes that the raw files for the observation are placed in ${outdir}/${obsid}/raw
export outdir="" # the root directory containing your observation data
export obsid="" # the specific observation to process
# preprocessor settings
export timeres_s=2
export freqres_khz=80
export edgewidth_khz=280
export metafits="${outdir}/${obsid}/raw/${obsid}.metafits"
export birli_ms="${outdir}/${obsid}/prep/birli_${obsid}_${timeres_s}s_${freqres_khz}kHz.ms"
birli \
--avg-freq-res ${freqres_khz} --avg-time-res ${timeres_s} \
--flag-edge-width ${edgewidth_khz} \
--no-rfi \
--emulate-cotter \
--ignore-dut1 \
--pfb-gains cotter \
-M "${birli_ms}" \
-m "${metafits}" \
${outdir}/${obsid}/raw/${obsid}_2*.fits
export cotter_ms="${outdir}/${obsid}/prep/cotter_${obsid}_${timeres_s}s_${freqres_khz}kHz.ms"
cotter \
-freqres ${freqres_khz} -timeres ${timeres_s} \
-edgewidth ${edgewidth_khz} \
-o "${cotter_ms}" \
-m "${metafits}" \
-allowmissing \
-noantennapruning \
-noflagautos \
-nostats \
-norfi \
${outdir}/${obsid}/raw/${obsid}_2*.fits
Birli will print some helpful information about the start times provided in the raw correlator files, along with the processing time.
For a detailed explanation of the MWALib concepts of the difference between "Provided", "Common" and "Good" times, check out this wki.
observation name: Sun
Array position: { longitude: 116.6708°, latitude: -26.7033°, height: 377.827m }
Phase centre: (164.9390°, 6.2214°) => (10h59m45.3566s, 6d13m17.1999s)
Pointing centre: (161.4930°, 5.2610°) => (10h45m58.3127s, 5d15m39.4799s)
Scheduled start: 2014-09-07 02:00:00.000 UTC, unix=1410055200.000, gps=1094090416.000, mjd=4916772000.000, lmst=132.7462°, lmst2k=132.5860°, lat2k=-26.6457°
Scheduled end: 2014-09-07 02:04:00.000 UTC, unix=1410055440.000, gps=1094090656.000, mjd=4916772240.000, lmst=133.7489°, lmst2k=133.5882°, lat2k=-26.6447°
Scheduled duration: 240.000s = 480 * 0.500s
Quack duration: 2.500s = 5 * 0.500s
Output duration: 237.000s = 474 * 0.500s
Scheduled Bandwidth: 30.720MHz = 24 * 32 * 40.000kHz
Output Bandwidth: 30.720MHz = 384 * 80.000kHz (2x)
Timestep details (all=480, provided=474, common=472, good=471, select=474, flag=13):
2014-09-07 UTC + unix [s] gps [s] p c g s f
ts0: 02:00:00.000 1410055200.000 1094090416.000 f
ts1: 02:00:00.500 1410055200.500 1094090416.500 f
ts2: 02:00:01.000 1410055201.000 1094090417.000 f
ts3: 02:00:01.500 1410055201.500 1094090417.500 f
ts4: 02:00:02.000 1410055202.000 1094090418.000 p c s f
ts5: 02:00:02.500 1410055202.500 1094090418.500 p c g s f
ts6: 02:00:03.000 1410055203.000 1094090419.000 p c g s f
ts7: 02:00:03.500 1410055203.500 1094090419.500 p c g s f
ts8: 02:00:04.000 1410055204.000 1094090420.000 p c g s f
ts9: 02:00:04.500 1410055204.500 1094090420.500 p c g s
ts10: 02:00:05.000 1410055205.000 1094090421.000 p c g s
...
ts474: 02:03:57.000 1410055437.000 1094090653.000 p c g s
ts475: 02:03:57.500 1410055437.500 1094090653.500 p c g s
ts476: 02:03:58.000 1410055438.000 1094090654.000 p s f
ts477: 02:03:58.500 1410055438.500 1094090654.500 p s f
ts478: 02:03:59.000 1410055439.000 1094090655.000 f
ts479: 02:03:59.500 1410055439.500 1094090655.500 f
Coarse channel details (metafits=24, provided=24, common=24, good=24, select=24, flag=0):
gpu corr rec cen [MHz] p c g s f
cc0: 1 0 57 72.9600 p c g s
cc1: 2 1 58 74.2400 p c g s
cc2: 3 2 59 75.5200 p c g s
cc3: 4 3 60 76.8000 p c g s
cc4: 5 4 61 78.0800 p c g s
cc5: 6 5 62 79.3600 p c g s
cc6: 7 6 63 80.6400 p c g s
cc7: 8 7 64 81.9200 p c g s
cc8: 9 8 65 83.2000 p c g s
cc9: 10 9 66 84.4800 p c g s
cc10: 11 10 67 85.7600 p c g s
cc11: 12 11 68 87.0400 p c g s
cc12: 13 12 69 88.3200 p c g s
cc13: 14 13 70 89.6000 p c g s
cc14: 15 14 71 90.8800 p c g s
cc15: 16 15 72 92.1600 p c g s
cc16: 17 16 73 93.4400 p c g s
cc17: 18 17 74 94.7200 p c g s
cc18: 19 18 75 96.0000 p c g s
cc19: 20 19 76 97.2800 p c g s
cc20: 21 20 77 98.5600 p c g s
cc21: 22 21 78 99.8400 p c g s
cc22: 23 22 79 101.1200 p c g s
cc23: 24 23 80 102.4000 p c g s
Antenna details (all=128, flag=8)
Baseline Details (all=8256, auto=128, select=8256, flag=996):
Estimated memory usage per timestep = 768ch * 8256bl * (32<Jones<f32>> + 4<f32> + 1<bool>) = 0.22 GiB
Estimated memory selected = 474ts * 768ch * 8256bl * (32<Jones<f32>> + 4<f32> + 1<bool>) = 103.56 GiB
Estimated output size = 474ts * 384ch * 8256bl * 4pol * (8<c32> + 4<f32> + 1<bool>) = 72.78 GiB
Preprocessing Context:
Will correct cable lengths.
Will correct digital gains.
Will correct coarse pfb passband gains.
Will not flag with aoflagger
Will correct geometry.
...
[2023-02-03T06:32:01Z INFO birli] init duration: 81.309324237s
[2023-02-03T06:32:01Z INFO birli] correct_passband duration: 20.408972661s
[2023-02-03T06:32:01Z INFO birli] read duration: 14.37518451s
[2023-02-03T06:32:01Z INFO birli] correct_geom duration: 3.533045269s
[2023-02-03T06:32:01Z INFO birli] write duration: 526.83874591s
[2023-02-03T06:32:01Z INFO birli] correct_cable duration: 13.050491358s
[2023-02-03T06:32:01Z INFO birli] correct_digital duration: 9.208166293s
[2023-02-03T06:32:01Z INFO birli] total duration: 668.723930238s
Cotter's output for the same task shows Birli taking a significant lead in processing time.
Running Cotter MWA preprocessing pipeline, version 4.6 (2022-01-20). Flagging is performed by AOFlagger 3.2.1 (2022-05-19). Input filenames succesfully parsed: using 96 files covering 4 timeranges from 24 GPU boxes. Detected 503.7 GB of system memory. Observation covers 72.32-103 MHz. Output resolution: 0.5 s / 80 kHz (time avg: 1x, freq avg: 2x). The first 8 samples (4 s), last 0 samples (0 s) and 7 edge channels will be flagged. Using a-priori subband passband with 32 channels. Using per-input subband gains. Average gains: 3.21436,2.90552,2.58398,2.29688,2.0874,1.97717,1.94177,1.93787,1.92065,1.88062,1.80994,1.74097,1.69104,1.67554,1.6897,1.71851,1.74182,1.7478,1.72314,1.68665,1.64905,1.62866,1.6272,1.64795 Subband order: 0,1,2,3,4,5,6,7,8,9,10,11,12,13,14,15,16,17,18,19,20,21,22,23 Observation's bandwidth is contiguous. All 480 scans fit in memory; no partitioning necessary. WARNING: start time according to raw files is 2014-09-07 02:00:02, but meta files say 2014-09-07 02:00:00 ! Will use start time from raw file, which should be most accurate. Warning: header specifies 480 scans, but there are only 472 in the data. Last 8 scan(s) will be flagged. Wall-clock time in reading: 00:04:42.520529 processing: 00:00:04.020476 writing: 00:25:38.408673
A python script can be used to display information about which timesteps are flagged by each preprocessor, which can vary slightly.
Some Measurement Set files created by Birli, Hyperdrive or Marlu can use the DUT1 timescale, so it's best to read the timescale specified in the column metadata. In this instance, the --ignore-dut1 flag was used to emulate Cotter's behaviour
for path in ['birli_1094090416_0.5s_80kHz.ms', 'cotter_1094090416_0.5s_80kHz.ms']:
print(path)
tb_main = table(os.path.join('/data/dev/1094090416/prep', path))
timescale = tb_main.getcolkeyword('TIME', 'MEASINFO')['Ref'].lower() or 'utc'
all_times = set()
unflagged_times = set()
for row_idx in range(tb_main.nrows()):
time = Time(tb_main.getcell('TIME', row_idx) / 86400.0, format='mjd', scale=timescale)
all_times.add(time.gps)
if np.any(np.logical_not(tb_main.getcell('FLAG', row_idx))):
unflagged_times.add(time.gps)
print(f"number of times: {len(all_times):12d}")
print(f"number unflagged: {len(unflagged_times):12d}")
print(f"first time: {min(all_times):12.4f}")
print(f"first unflagged: {min(unflagged_times):12.4f}")
print(f"last unflagged: {max(unflagged_times):12.4f}")
print(f"last time: {max(all_times):12.4f}")
The results of running this command on both files show that Birli and Cotter flags start and end at different times (these are centroid times, where the Birli logs above are start times).
birli_1094090416_0.5s_80kHz.ms Successful readonly open of default-locked table /data/dev/1094090416/prep/birli_1094090416_0.5s_80kHz.ms: 23 columns, 3913344 rows number of times: 474 number unflagged: 467 first time: 1094090418.2500 first unflagged: 1094090420.7500 last unflagged: 1094090653.7500 last time: 1094090654.7500 cotter_1094090416_0.5s_80kHz.ms Successful readonly open of default-locked table /data/dev/1094090416/prep/cotter_1094090416_0.5s_80kHz.ms: 23 columns, 3962880 rows number of times: 480 number unflagged: 464 first time: 1094090418.2500 first unflagged: 1094090422.2500 last unflagged: 1094090653.7500 last time: 1094090657.7500
Birli flags all times that are within quack time (in this case, 2.5s) after the first provided timestamp, while Cotter flags more timesteps by default.
Taql can also be used to enumerate timestamps in both measurement set files
for proc in birli cotter; do export table="${proc}_${obsid}_${timeres_s}s_${freqres_khz}kHz.ms"; taql -o "${proc}_times.txt" -p "select unique TIME from $table"; done;
From these timestamps, we can construct a taql file to print a table of visibility magnitude values from a group of timesteps in the middle of the observation, looking at the band of frequencies bookended by the edge channels of the first coarse channel, showing only visibilities from the first (XX) polarization.
comparison.taql
SELECT
TIME,
mscal.ant1name() as ant1name,
mscal.ant2name() as ant2name,
STR(FLATTEN(abs(DATA[3:13, 0])),"%9.3f") as ch057XX
FROM $table
WHERE
TIME>=MJD("07-Sep-2014/02:00:10.000")
AND TIME<MJD("07-Sep-2014/02:00:50.000")
AND mscal.ant1name() in ["Tile011","Tile021","Tile031"]
AND mscal.ant2name() in ["Tile012","Tile022","Tile032"]
this taql file can be run on both measurement sets, and the results can be compared with your preferred diff viewer
for proc in birli cotter; do export table="${proc}_${obsid}_${timeres_s}s_${freqres_khz}kHz.ms"; taql -f comparison.taql -o ${proc}_ch057XX.txt -p; done;
diff -u birli_ch057XX.txt cotter_ch057XX.txt | diff-highlight
Although the visibilities are not identical because of the reasons mentioned above, there is negligible difference between the two.
