I would like to dedicate this thesis to my dadu (maternal grandfather), who
supported me throughout my childhood and always wanted me to pursue higher education. I wish you were here with us now; I shall do my best to make you proud.
i Abstract
Mm and sub-mm wavelengths are crucial for studying the process of star formation and probing the hot matter around supermassive black holes, as these environments tend to be optically thick at visible and near-IR wavelengths. The atmosphere plays an important role in mm and sub-mm astronomy as it can absorb, and scatter radia- tion. Also, varying atmospheric emission provides additional noise. High-resolution, sensitive instruments are required to probe distant objects with small spatial scales (like star-forming cores). In this work, we first test whether sub-mm polarized dust emission can be used to trace magnetic fields in cores by analyzing the polarized ther- mal dust emissions from starless, prestellar, and protostellar cores. Using stacking analysis, we find that cores that have a bright 70 µm source show a higher polarized intensity signal when stacked, compared to ones that do not. This could indicate that cores that have luminous internal objects may produce polarized emissions that can be used to study their magnetic fields. However, we would need a larger sample of cores at different evolutionary stages to statistically determine the correlation be- tween core properties and polarization efficiency. This is one of the primary science goals for upcoming telescopes, such as the CCAT. We use atmospheric models to investigate the sensitivity of the CCAT’s Prime-Cam instrument to changing atmo-
ii iii
spheric loads. We find that a tone-tracking algorithm will be required to prevent
the Microwave Kinetic Inductance Detector (MKID) detectors from saturating. We also find a lower detector responsivity value would be needed for the 850 GHz array compared to 250 and 350 GHz arrays. Finally, we push to an even higher reso- lution by designing the Balloon-borne VLBI Experiment (BVEX). We lay out a plan for a complete backend readout and data storage system using the latest radio frequency system on chip (RFSoC) boards and show a successfully deployed fully working spectrometer. BVEX will be a pathfinder mission showing the possibility of a balloon-borne VLBI station that would be part of a future high-frequency VLBI network. Acknowledgements
I want to take this opportunity to thank my supervisor Prof Laura Fissel for not only supporting me through excellent science advice but also for letting me into the world of astronomy and instrumentation, something which I had always wanted to do! I will be forever grateful to you for this. You are a true inspiration. All the projects have been an excellent learning experience for me and have taught me skills to carry forward into my Ph.D. and later into my professional career. I also thank Prof Sarah Sadavoy for being incredibly supportive and kind. I appreciate all the scientific inputs and feedback you have provided me. It has been beneficial. A big thanks to Adrian Sinclair for your continuous guidance and valuable input on the CCAT project. I want to thank the CASPER collaboration for their support in implementing the spectrometer design to the RFSoC. I would like to thank everyone from the Qstar team especially Felix Thiel for being an excellent BVEX collaborator. Finally, a special thanks to my family, mom, dad, my sister Mahima, and my best friend, Sayaree. Without you all, I would not have been able to do what I am doing now. Thank you for being there for me always.
BLAST Balloon-borne Large Aperture Submillimeter Telescope
BRAM Block Random Acess Memory
CCAT Cerro Chajnantor Atacama Telescope
CMB Cosmic Microwave Background
CFI Canada Foundation for Innovation
CNES Centre National d’Études Spatiales
CSA Canadian Space Agency
List of Figures xv
CASPER Collaboration for Astronomy Signal Processing and Electronics Research
DAC Digital-to-Analog Converter
DSP Digital Signal Processor
EHT Event Horizon Telescope
ESA European Space Agency
ENOB Effective Number of Bits
FPGA Field-Programmable Gate Array
FITS Flexible Image Transport System
FWHM Full Width at Half Maximum
FYST Fred Young Submillimeter Telescope
FFT Fast Fourier Transform
FIR Finite Impulse Response
GSPS Giga Samples Per Second
GMC Giant Molecular Cloud
GPS Global Positioning System
GBT Green Bank Telescope
HFI High Frequency Instrument
List of Figures xvi
ISM Interstellar Medium
IF Intermediate Frequency
JCMT James Clerk Maxwell Telescope
LMT Large Millimeter Telescope
LMC Large Magellanic Cloud
LFI Low-Frequency Instrument
LWP Liquid Water Profile
LNA Low-Noise Amplifier
LSB Least Significant Bit
MKID Microwave Kinetic Inductance Detector
MSB Most Significant Bit
NSD Noise Spectral Density
ngEHT next generation Event Horizon Telescope
NCO Numerically Controlled Oscillator
OCXO Oven-Controlled Crystal Oscillator
PWV Precipitable Water Vapor
PAH Polycyclic Aromatic Hydrocarbon
List of Figures xvii
PACS Photodetector Array Camera and Spectrometer
PIC Programmable Integrated Circuit
PPS Pulse-Per-Second
PFB Polyphase Filterbank
PID Proportional-Integral-Derivative
RFSoC Radio Frequency System-on-Chip
RH Relative Humidity
RAT Radiative Aligmnent Torque
RF Radio Frequency
RMS Root Mean Square
RFI Radio Frequency Interference
RFDC Radio Frequency Digital Downconverter
SMBH Supermassive Black Hole
SPIRE Spectral and Photometric Imaging Receiver
SMC Small Magellanic Cloud
SEFD System Equivalent Flux Density
SFDR Spurious-Free Dynamic Range
List of Figures xviii
SNR Signal-to-Noise Ratio
TMC Taurus Molecular Cloud
TES Transition Edge Sensor
TICC Timestamping/Time Interval Counter
VLBI Very Long Baseline Interferometry
VCXO Voltage-Controlled Crystal Oscillator
WCS World Coordinate System
YSO Young Stellar Object
Chapter 1
Introduction
1.1 Motivation
With advances in theory, to better understand and characterize our universe, the need for new instrumentation is steadily rising. Astronomical instrumentation has come a long way, from verifying general relativity near black holes by imaging their event horizon to studying the early universe by identifying some of the oldest stars and galaxies. To probe star-forming regions or study the photon rings around su- permassive black holes, we need to observe the dust emission at mm and sub-mm wavelengths, as these environments have significant extinction at optical wavelengths. Magnetic fields are thought to play a key role in the process of star formation as they provide support against the gravitational collapse of dense clumps of gas. However, the importance of the magnetic field in the early stages of star formation is still unclear. Cores are dense gas precursors to individual stellar systems, and to resolve them, we need spatial resolution better than 0.1 pc. In general, polarization levels tend to be low, perhaps indicating that the dust grains are not well aligned in cores
1 1.1. Motivation 2
and that most of the polarized radiation observed towards these cores is from the less dense gas envelopes in the surrounding cloud. To trace the magnetic field mor- phology, we first need to determine whether cores emit polarized radiation. Next, we would need a large sample of core polarization maps to statistically determine the correlation between core properties and polarization efficiency. Supermassive black holes are thought to be at the center of almost every galaxy. These black holes have masses above 106 M⊙ and play a crucial role in shaping the evolution of galaxies. The EHT collaboration, using Very Long Baseline Interfer- ometry (VLBI), analyses the radiation from the trapped photons at the inner orbits of these black holes as they carry vital information like size, mass, and spin. The distances to Supermassive Black Hole (SMBH), like 8.5 kpc (Genzel et al., 2010) for Sagittarius A∗ and 16.8 Mpc (Herrnstein et al., 1999) for M87∗ would require us to push to higher angular resolution to image their event horizons. Hence, there is a clear need to improve the angular resolution of our observations. There are two ways to do that; according to the formula for the diffraction limited resolution of a telescope, λ θ≈ (1.1) D
(where λ is the wavelength at which we are observing and D is the diameter of the dish):
1. Increase the diameter/baseline of the dishes (currently limited by the size of
the Earth),
2. Push to shorter wavelengths (sub-mm).
1.1. Motivation 3
Figure 1.1: First image of a black hole shadow taken using the EHT observations of the center of the galaxy M87 (Collaboration et al., 2019).
EHT already achieves an incredible angular resolution using VLBI or very long baseline interferometry. This technique uses intercontinental baselines (∼ 1000 km), combining some of the best mm/sub-mm radio telescopes to achieve this. However, at the current 350 GHz range, we can only push to a maximum angular resolution of θ = c/vD ≈ 17 µarcsec before getting limited by the diameter of the Earth. We need to push to shorter wavelengths to resolve more black holes and study compact extra-galactic star-forming cores. Currently, on one hand, we have space telescopes like P lanck, which provide us with 5 arc-min resolution with no interferometry provision; on the other, we have in- credibly powerful interferometers like the Atacama Large Millimeter Array (ALMA), which go down to a milli-arcsec scale. With increased angular resolution, we tend to lose sensitivity on larger spatial scales, and this could be a problem if we try to resolve nearby star-forming dense cores. There is a gap that upcoming mm/sub-mm telescopes like the Cerro Chajnantor Atacama Telescope (CCAT) shown in Figure 1.2 and the TolTEC camera on the Large Millimeter Telescope (LMT) are trying 1.1. Motivation 4
to fill. These instruments will be able to resolve up to 5′′ , allowing us to probe the dense cores within our local molecular clouds.
Figure 1.2: An artist’s illustration of the upcoming CCAT telescope on the Cerro Chajnantor plateau, CPI Vertex Antennentechnik GmbH (used with permission).
Millimeter and sub-mm measurements require excellent weather, accurate anten-
nas, and sensitive receivers. Hence it is common to see flagship mm/sub-mm radio telescopes at high altitude sites with dry atmospheric conditions. Eventually, we will be limited by the atmosphere as we observe at higher frequencies. We would need to observe from even higher altitudes to surpass most of the atmosphere. This will also allow VLBI to observe at higher frequencies, thereby improving the resolu- tion even further, allowing the EHT to image more black hole shadows and resolve distant compact sources. Balloon-borne astronomy can be a great solution as bal- loons can ascend high enough to escape 99.5% of the atmosphere, allowing us to capture the high-frequency radiation that could not be observed from most ground- based observatories at just a fraction of the cost of space telescopes. With sufficient 1.1. Motivation 5
advancements in new-generation detectors, stable timing reference systems, high-
speed readout systems, and ultra-precise position tracking systems, balloon-borne telescopes capable of doing VLBI in the near future can push mm and sub-mm astronomy to new heights. Submillimeter astronomy observes light from a wavelength of 200 µm to 1 mm, thereby bridging the gap between infrared and radio bands. However, the atmosphere can be a huge hindrance at these wavelengths as water vapor and other molecules absorb most of the incoming mm photons. Therefore, it is common to see mm and sub-mm sites on high-altitude, dry sites. Even at the best ground-based observing sites, it becomes important to characterize the atmosphere above the telescope as the atmospheric transmission varies significantly depending on the weather. Atmo- spheric modeling predicts how much of the incoming radiation will be absorbed by the atmosphere before it reaches the telescope. It also gives us the radiation power from the atmosphere which will dominate the loading of the instrument. This is important because it helps us determine which wavelengths are best suited for ob- servations and also helps us design important aspects of detectors to maximize their performance. With the advent of new-generation superconducting photon detectors like MKID, astronomical observations have got a noticeable boost, especially in far infrared and millimeter. Due to their sensitivity, and fabricating detector focal plane arrays with thousands of detectors, MKIDs have proven to be an essential technology for building high sensitivity mm and sub-mm cameras. With frequency domain multiplexing a single radio frequency coaxial cable can read out thousands of MKID detectors. Each 1.2. Research Problem 6
detector has a resonance frequency and only has a linear response to radiation near this resonance frequency, but the resonant frequency can change as the incident power from the atmosphere or electronics changes. Current readout techniques significantly restrict their potential dynamic range. Hence, to unlock the true potential of MKIDs, we need to actively track the resonant frequency of the detectors within the array by implementing a firmware known as tone-tracking. With new detector technologies and VLBI pushing to even higher frequencies, there will be an exponential rise in data output and bandwidth. For example, the next generation Event Horizon Telescope (ngEHT) will double the EHT’s recording rate of 64 gigabits per second. With such massive data output, we need adequate backend readout systems to efficiently sample signals, process them, and store them. The latest generation Radio Frequency System-on-Chip (RFSoC)s from AMD’s Xil- inx is being heavily pushed in the market. These devices eliminate bulky Analog- to-Digital Converter (ADC) and Digital-to-Analog Converter (DAC) and condense both onto a single chip which can be programmed via a Field-Programmable Gate Array (FPGA) board. These FPGA boards are slowly overtaking existing readout systems all across the world.
1.2 Research Problem
My research attempts to address the following questions:
1. Can we detect polarized dust emission within molecular cores, or do they not produce significant levels of polarized emission, perhaps because of the lack of grain alignment? If the answer is yes, then we can use polarization to study 1.3. Objective 7
magnetic fields within cores.
2. How will next-generation dust polarimeters and cameras, which will survey thousands of cores using sub-mm MKID detectors, operate in ground-based observing locations where the power loading due to the atmosphere changes with telescope elevation?
3. For what observing frequencies will a tone-tracking system be required to make
several degree-sized maps?
4. Can we develop techniques to probe at an even higher resolution in the sub-mm
with balloon-borne VLBI? What could be the technical challenges, and how would we solve them? Can new commercially available radio frequency FPGA boards be used to simplify the design of these telescopes?
1.3 Objective
To answer the questions identified above, we will need to :
1. Analyse polarized dust emissions at the core scale (∼0.05 to 0.1 pc) of nearby molecular clouds as it can be an important tool for understanding the mag- netic field morphology at those scales. We will use catalogs of the different evolutionary stages of cores and extract their stokes parameters and polariza- tion intensities using the data from P lanck and James Clerk Maxwell Tele- scope (JCMT) POL-2 850 µm maps at the location of cores of NGC1333. To determine whether the cores are polarized we will normalize each map and stack 1.3. Objective 8
them together. Further, we will investigate whether the polarization intensity
is higher for the cores that show protostellar activity (contain a bright 70 µm source) when compared to the ones that do not. Stacking will help us improve the signal-to-noise ratio, allowing us to examine the statistical properties of polarized emission of cores.
2. For Prime-Cam to establish the need for tone tracking which will be used to actively track the shift in the resonant frequency of the detectors, we will, first of all, need to characterize the atmosphere above the Cerro Chajnantor site (5600 m above sea level) where the 6m telescope is located using the Atmospheric Model simulator by Paine (2012). Determining the atmospheric transmission, opacity, brightness temperatures, etc., as a function of frequency will be crucial in planning observation runs. Since all of these parameters vary over the year with different Precipitable Water Vapor (PWV) values, it will significantly impact the quality of data and, finally, the performance of MKID detectors. We will also analyze the various frequency bands as this will show how far in elevation the telescope can move without causing the MKID detectors to depart from their regime of linear response and therefore requiring a full retuning of all of the detectors.
3. Lastly, to address the issues of pushing VLBI observations to higher frequen-
cies and resolution, we will need to demonstrate that VLBI is possible from balloon-borne telescopes. BVEX is a pathfinder experiment that will demon- strate exactly that but at K-band (21-23 GHz) frequency. We will need to design a low-powered sufficient backend system, readout electronics, stable 1.4. Organization of thesis 9
timing references, and an ultra-precise position tracking system. To provide
stable operating conditions to highly sensitive Oven-Controlled Crystal Os- cillator (OCXO), we will need to design a vibrationally stable pressure and temperature-controlled vessel. Finally, to handle the massive amounts of raw data 4 Giga Samples Per Second (GSPS), compute the spectrum of the source, and store them efficiently, we will be designing an RFSoC FPGA-based readout system.
1.4 Organization of thesis
This thesis is divided into five chapters. In Chapter 2, we provide a background
on star formation and focus on the core stage. We perform a stacking analysis of the polarization properties of 57 cores, which helps us determine that active cores (bright in 70µm) show a higher polarized intensity than their inactive counterparts. In Chapter 3, we discuss CCAT and how it is crucial for studying star-forming cores at higher resolution. We will also discuss the MKID detectors that CCAT will use and how they can be efficiently read out. Finally, we discuss our work on characterizing the atmosphere above the CCAT site, which will help plan observing strategies and establish a need for a tone-tracking routine for the MKID detectors. In Chapter 4, we introduce BVEX, a pathfinder experiment that will demonstrate simultaneous VLBI observations from a balloon-borne telescope. We present the designs for a timing reference system and RFSoC-based backend readout and spectrometer. Finally, in Chapter 5, we conclude this work and discuss some of the next steps for each of these projects. Chapter 2
Polarized emission of star-forming cores
This chapter will present some of the theory and background information on star formation, especially towards the very early core stage. By looking at the polarized thermal radiation from dust grains, we can map out the magnetic field morphology projected on the plane-of-sky. It is unclear whether or not dust grains within dense cores will be aligned with the magnetic field (since the short wavelength light that is thought to align grains may not be attenuated by dust in these cores). Therefore it isn’t clear whether or not polarization observations can be used to tell us about mag- netic fields within cores. Hence in this chapter, a statistical study will be carried out to test for evidence of polarized emission from cores by stacking linear polarization maps at core locations within molecular clouds.
10 2.1. Background 11
2.1 Background
2.1.1 The hierarchical structure of star formation
Star formation is a multi-scale process. Stars form due to the gravitational collapse of gas inside Giant Molecular Cloud (GMC), which are large-scale structures in the Interstellar Medium (ISM). The ISM is the matter and radiation that exists be- tween star systems in a galaxy. These GMCs are composed of smaller and denser filamentary structures which contain dense cores and envelopes with young proto- stellar objects embedded in them as in Figure 2.1. Gravity, turbulence, and magnetic fields are the main physical mechanisms associated with this multi-scalar spatial dis- tribution of gas. However, a big area of uncertainty in star formation research is how these different physical processes interact with each other.
Figure 2.1: A cartoon showing the hierarchical multi-scale process of star formation, where clouds have the largest spatial scales followed by clumps, filaments, cores, envelopes, and finally disks and stars. The density similarly increases from clouds to disks, Pokhrel et al. (2018)(used with permission from author). 2.1. Background 12
The principal ingredients of the cold, dense clouds in the ISM are dust and gas. By number density, molecular clouds are mostly composed of hydrogen molecules. The largest of these molecular clouds, known as GMC, have temperatures of roughly 10-30 K, average densities of up to 102 -103 particles/cm3 , masses of a few to more than a million solar masses, and sizes of 20 to 200 parsecs (Solomon et al., 1987; Dame et al., 2001). The exact formation process of these GMCs is not straightforward; however, a study by Inutsuka et al. (2015) suggests that repeated interaction of expanding neutral hydrogen (HI) bubbles might be responsible for their formation.
Figure 2.2: Cartoon showing the sequential formation of molecular clouds by the interaction of expanding shells Inutsuka et al. (2015)(used with permission from the author).
In this model, molecular clouds form when multiple compressions of thick HI shells overlap, propelled by expanding bubbles. As we can see in Figure 2.2, they only form in a few specific places: those with excessive compression or where the compressional direction is almost parallel to the local mean magnetic field lines. The interplay of 2.1. Background 13
gravity, magnetism, and turbulence within these molecular clouds is responsible for forming denser clumps and filaments. Just like GMCs, the dominant mechanism of filament formation is still uncertain. In the HI phase of the ISM, filamentary structures are generally believed to form in shock-compressed layers (Inoue et al., 2009) induced by thermal instability. The turbulent shear flow along the magnetic field (Soler and Hennebelle, 2017) stretches these shocks. On the other hand, in the denser regions of the ISM (GMC), filaments may form due to direct compression of interstellar turbulence or gravitational fragmentation (Pineda et al., 2022).
2.1.2 Interstellar dust grains
Dust grains are crucial to the study of the ISM and influence the chemistry and thermodynamics of molecular gas. According to Bernstein et al. (2002), dust in the universe reradiates at least 30% of the energy that stars release as light. Using observations of extinction, scattering, and emission from dust we can infer the abun- dance, composition, and size of dust grains. The primary constituents of dust are silicates, which contribute a significant fraction of dust mass (Li and Draine, 2002), carbonaceous materials including diamond, graphite, Polycyclic Aromatic Hydrocar- bon (PAH)s, aliphatic hydrocarbons, and, carbonates including calcite and dolomites (Kemper et al., 2002). Depending on their composition, interstellar dust grains have a wide distribution of sizes. Most of the silicates and metal-based dust particles typically have diameters between 0.1 µm and 0.5 µm (Draine, 2003). Carbonaceous material-based dust grains on the other hand are thought to have a diameter ≲ 100 Å (Draine, 2003). 2.1. Background 14
2.1.3 Polarized thermal dust emission
GMCs forming stars have ∼ 1% by mass of interstellar dust mixed with the gas. Therefore, dust emission is an excellent probe of star-forming regions. These dust grains absorb some of the incoming unpolarized optical and U V light thereby heating the grains up (Draine and Li, 2001). They then re-emit the light in infrared and radio frequencies ν through a modified Blackbody radiation,
2hν 3 1 I(ν, T ) = 2 hν/kT [1 − e−τ (ν) ], (2.1) c e −1
where c is the speed of light, h is Planck’s constant, ν is the observing frequency, and τ (ν) is the optical depth at a given frequency ν. If τ (ν) ≪ 1 (optically thin regime) 3 then e−τ (ν) ≈ 1−τ (ν) and I(ν, T ) ≈ τ (ν) 2hν c2 1 ehν/kT −1 . On the other hand if τ (ν) ≫ 1 2hν 3 (optically thick regime) then e−τ (ν) ≈ 0 and I(ν, T ) ≈ B(ν, T ) = 1 c2 ehν/kT −1 , where B(ν, T ) is the specific intensity of a blackbody. The optical depth τ is usually modeled as having a power-law dependence in the mm/sub-mm,
τ ∝ νβ (2.2)
where β is the dust emissivity spectral index.
Dust grains are preferentially aligned with their short axis parallel to the direc- tion of the local magnetic field as shown in Figure 2.3 and discussed in Section 2.1.4. The wavelength of thermal emission is typically much larger than the size of a dust grain (∼0.1 µm), so grains are not that efficient at emitting thermal radiation. The dust grains are usually longer in one direction, which effectively makes a dust grain 2.1. Background 15
a better antenna in the direction parallel to its long axis compared to its short axis. Therefore it radiates with the E-field (and therefore polarization direction) prefer- entially aligned parallel to the long axis of the grain (Hildebrand, 1988). Hence, the polarization of the thermal dust emission is preferentially aligned perpendicular to the magnetic field. Therefore, studying polarized thermal dust emission can be an excellent way to map the local magnetic field orientation within molecular clouds. This means that the polarized emission from dust grains will depend on its temper- ature, emissivity, density, how well aligned the grains are with the magnetic field, and efficiency at producing polarization (for example, the degree of anisotropy of the dust grain shapes).
Figure 2.3: A cartoon illustration showing the dust orientation in molecular clouds. In our analysis, we are analyzing the polarization of thermal radiation as a result of dust grains re-emitting the absorbed light.
All of this polarization information is encoded within the Stokes Parameters I,
Q, U , and V . For a plane electromagnetic wave propagating in the z-direction, the electric field vectors can be expressed as, 2.1. Background 16
Ex = ax (t) cos[ω0 (t) − θx (t)] (2.3)
and, Ey = ay (t) cos[ω0 (t) − θy (t)], (2.4)
where ax and ay are the amplitude and θx and θy are the phase angles. In the case of polarized light, the equations mentioned above will have some correlation between the two components. Stokes Parameters can be defined as the different combination time average of the amplitudes,
I ≡ ⟨a2x ⟩ + ⟨a2y ⟩ (2.5)
Q ≡ ⟨a2x ⟩ − ⟨a2y ⟩ (2.6)
U ≡ ⟨2ax ay cos(θx − θy )⟩ (2.7)
V ≡ ⟨2ax ay sin(θx − θy )⟩ (2.8)
Stokes I corresponds to the total intensity, Q corresponds to horizontal or vertical
linearly polarized light, U corresponds to linearly polarized light at an angle of ± 45◦ , and V corresponds to circularly polarized light. Therefore, the linearly polarized light intensity (P ) can be expressed as, 2.1. Background 17
P = (Q2 + U 2 )1/2 . (2.9)
The linear polarization fraction(p) can be calculated from polarized light intensity by,
(Q2 + U 2 )1/2 P p= = (2.10) I I
Finally the polarization angle θ can be expressed as ,
1 θpol = tan−1 (U/Q) (2.11) 2
The polarization angle θpol traces the plane-of-sky orientation of the long axis of the dust grain along the line of sight. Where, θpol = 0 implies that a vector is pointing towards celestial North. As θpol increases the vector rotates counter-clockwise (that is East of North). To get the plane-of-sky magnetic field, we need to add by 90◦ ,
1 B̂pos = θpol + 90◦ = tan−1 (U/Q) + 90◦ (2.12) 2
The plane-of-sky magnetic field B̂pos is not actually describing vectors with mag- nitude and direction. Rather, B̂pos is a ‘pseudovector’ which only gives us the relative orientation of the plane-of-sky magnetic field. 2.1. Background 18
2.1.4 Theory of radiative alignment
The Radiative Aligmnent Torque (RAT) model (Draine and Weingartner, 1996) is the leading candidate to explain grain alignment. We established that the alignment of the dust grains within the molecular clouds gives rise to polarization. A dust grain in interstellar space will initially rotate around a random axis (and wobble). Radiative grain alignment originates in the differential extinction cross sections for the left- and right-hand circular polarization components of the incident light on the irregular grains exposed to an anisotropic radiation field with λ < 2a where 2a is the grain diameter. The asymmetry in the grains results in torques on them, causing them to spin up (increased angular momentum). If the grain material has unpaired spins (paramagnetic, e.g. silicates) it can minimize energy while keeping angular momentum fixed by trading rotation for spin- flips thus causing magnetization. This is known as the Barnett effect. In the presence of an external magnetic field the magnetized grain Larmor precesses around the field direction. The continued action of radiative torques over the precession period causes the angular momentum of the grains to eventually preferentially align with the magnetic field. Figure 2.4 depicts an ellipsoidal-shaped grain aligned with its spin axis J⃗ parallel to its maximum moment of inertia a1 . The grain precesses about the magnetic field ⃗ with an angle ξ. The incident radiation field at an angle ψ with respect to B B ⃗ ⃗ || to J⃗ and F⃗ ⊥ H. causes a net torque on the grain with components H ⃗ Since the
alignment torque F⃗ is perpendicular to J,
⃗ the fast precession of the grain causes
the average value of F⃗ to approach zero as J⃗ approaches B
of ψ angles, for which the grains align with their long axes parallel to the magnetic field, occurs in a state of low angular momentum for the grains (Lazarian and Hoang, 2007). The assumption of J⃗ perfectly parallel to the grain short axis(⃗a) does not hold, due to thermal fluctuations of the grain energy causing the grain spin axis to wobble slightly. However the radiative torques tend to increase J⃗ so that their rotation is superthermal, therefore fluctuations in the grain energy do not cause significant changes in the orientation of J⃗ with respect to ⃗a. Overall the resulting range of wobbling angles at a low-J⃗ point is much wider than the range of stable angles for which the alignment occurs with J⃗ ⊥ B ⃗ (Lazarian and Hoang, 2008). Therefore,
in the presence of internal excitations (wobbling), grains under radiative alignment
torque preferentially such that J⃗ || B. ⃗ 2.1. Background 20
Polarized thermal emission of dust grains has a strong dependence on grain
alignment in the local environment. The radiative alignment torques as we dis- cussed depend on the degree of anisotropy of the local radiation field, the size of the grains, ambient gas temperature, and density distribution. Single-dish observa- tions of molecular clouds (Poidevin et al., 2013; Fissel et al., 2016), starless cores (Alves et al., 2014), and high-resolution interferometric observation of protostellar cores (Hull et al., 2014) have found a significant reduction in polarization fraction (p) with an increase in column density. This drop in polarization fraction has been interpreted as either a) depolarization or loss of polarization due to disorganized magnetic field lines within the volume of gas probed by the telescope beam or b) loss of alignment efficiency of the dust particles due to lack of irradiation and changes in grain characteristics towards high column density regions. Also particularly for cores it is not clear whether thermal dust polarization is a good tracer of magnetic fields because these cores are found within dense filamentary regions and short wavelength photons needed to align the dust grains may not be able to penetrate deeply into the core envelopes (Jones et al., 2014).
2.1.5 Types of star-forming cores
Cores are dense gas precursors to individual stellar systems. The various evolutionary stages of cores can be summarised as follows:
• Starless cores, are systems that may or may not collapse to form a star. They are a transient concentration of molecular gas and dust without any embedded Young Stellar Object (YSO), André et al. 2009. They provide us with a snap- 2.2. Data 21
shot of the initial conditions of star formation, and even of the first stages of the process if they have already started to collapse (Kirk et al., 2007).
• Prestellar cores, are gravitationally bound starless systems that represent a
somewhat denser and more centrally-concentrated population of cores (André et al., 2009). They are typically observed in sub−mm dust continuum emission and dense gas tracers such as N H3 (Caselli et al., 2002). Prestellar cores fre- quently exhibit evidence of infall motion and are characterized by large density contrasts over the local background medium.
• Dense cores or protostellar cores are believed to result from the gravitational collapse of prestellar cores. They are the highest-density structures in molecu- lar clouds and largely control the origin of the stellar initial mass function (IMF; André et al. 2009). These cores may have protostars or compact luminous ob- jects. Protostellar cores also show evidence of excess far-IR or mid-IR emission and may drive molecular outflows (which indicates an accreting protostar).
2.2 Data
For our analysis, we chose nearby molecular clouds that would allow us to resolve the maximum number of cores within them ideally to <0.1 pc resolution. Since we are interested in studying the polarized thermal dust emission from these cores we used the 353 GHz P lanck Legacy Archive PR4 data (Akrami et al., 2020) of the Taurus molecular cloud referred in Figure 2.5 and the 353 GHz (850 µm) JCMT POL-2 (Doi et al., 2020) data of NGC1333, a star-forming region in the Perseus 2.2. Data 22
molecular cloud. To extract the location of cores from these maps we used the “derived cores catalog” for Perseus (Sadavoy et al., 2012) and Taurus (Marsh et al., 2016) from the Herschel Gould Belt Survey archival data and the SCUBA-2 Gould Belt survey (private communication with Kate Pattle). To identify the cores that show protostellar activity (brighter in 70 µm) we also used the 70 µm dust intensity map for the Perseus molecular cloud (Sadavoy et al., 2014; Pezzuto et al., 2021).
Figure 2.5: P lanck 353 GHz map of the Taurus Molecular Cloud. The first figure shows the Stokes I map followed by Stokes Q and U map in units of a thermodynamic temperature (kcmb) (Akrami et al., 2020). 2.2. Data 23
2.2.1 P lanck data
P lanck was a space observatory operated by the European Space Agency (ESA) from 2009 to 2013 at the Sun-Earth Lagrange point (L2). It was an all-sky Cosmic Mi- crowave Background (CMB) mapping experiment at microwave and infrared frequen- cies (27 GHz – 857 GHz) with high sensitivity and angular resolution. P lanck had a 1.9m x 1.5m Gregorian telescope. It had two primary instruments, a Low-Frequency Instrument (LFI), which observes 30, 44, and 70 GHz bands at a resolution of 33’, 24’, and 14’, respectively, with a bandwidth of 20% of the central frequency, and a High Frequency Instrument (HFI) which observes at 100, 143, 217, 353, 545, and 857 GHz bands at a resolution of 10, 7.1, 5.5, 5.0, 5.0 and 5.0 respectively with a 33% bandwidth of the central frequency. The HFI had a cryogenically cooled bolometer array. We used the all-sky 353 GHz P lanck Legacy Archive PR4 data as it covered every star-forming region. We extracted the Stokes I, Q, and U data centered at a Galactic longitude of 174◦ and a latitude of 13.45◦ for the nearby star-forming Taurus Molecular Cloud (TMC) region for our initial analysis. At 140 pc the TMC is one of the closest molecular cloud to us. With the 5’ resolution of the P lanck data, we can resolve structures up to ∼0.2 pc which is a little bit larger than the expected size of cores (0.05 pc to 0.1 pc). We further smoothed the map with a 4-pixel Full Width at Half Maximum (FWHM) Gaussian kernel to reduce the noise, where each pixel corresponds to an angular size of 1.5′ . 2.2. Data 24
2.2.2 JCMT POL-2 data
Figure 2.6: JCMT POL-2 850 µm map of NGC1333 from the BISTRO Survey. The first figure shows the Stokes I map followed by Stokes Q, U map in units of pW (Doi et al., 2020). 2.2. Data 25
The JCMT is a 15m submillimeter telescope at the summit of Mauna Kea operated by the East Asian Observatory. It is primarily used to study our solar system, interstellar and circumstellar gas and dust, and distant galaxies. The JCMT uses heterodyne receivers, continuum bolometer array instruments, and an ancillary linear polarization-sensitive continuum instrument (Holland et al., 2013). We use data from the POL-2 instrument, a linear polarimeter sensitive to Stokes I, Q, and U vectors. The POL-2 works with the SCUBA-2 detectors observing at 850 and 450 microns. The BISTRO survey (Pattle et al., 2017) uses the POL-2 to map the polarization in the densest part of the Gould Belt star-forming regions to trace the magnetic field direction and strength at the centers of all these regions. We use the 850 µm Stokes I, Q, and U data of NGC 1333 (Doi et al., 2020), which covers the entire ∼1.5 pc × 2 pc region. NGC 1333 was observed at two positions, namely, field 1 (α = 03h 29m 03s .450, δ = 31◦ 14′ 34”.70) and field 2 (α = 03h 29m 06s .350, δ = 31◦ 20′ 22”.50) between August 2017 to January 2018 with a total exposure time of 28.3 hours. The spatial distribution of the Stokes I, Q, and U parameters from their time-series measurements was reduced using the STARLINK software package (Parsons et al., 2018). At a distance of 300 pc with 14.4′′ resolution, we can resolve structures down to ∼0.02 pc. Although it is not enough to spatially resolve individual protostellar envelopes (∼ 0.005 pc or ∼ 1000 au), it allows us to resolve the polarized thermal emission from their surrounding core envelopes. 2.2. Data 26
2.2.3 Herschel Gould Belt survey
The Herschel Space Observatory was a 3.5 m infrared telescope built and operated by ESA. It was sensitive to far-infrared and submillimeter wavelengths ranging from 55 to 672 µm. The two imaging instruments, Spectral and Photometric Imaging Receiver (SPIRE) and Photodetector Array Camera and Spectrometer (PACS), are ideal for imaging prestellar cores and young protostars. The Gould Belt Survey (André et al., 2010) was designed to extensively image the densest portion of the Gould Belt (nearby star-forming regions) with SPIRE at 250-500 microns and PACS at 70-170 µm with a 5σ column density sensitivity (NH ∼1021 cm−2 ) to understand the physical mechanisms for the formation of prestellar cores out of the diffuse medium. In this work, we used the derived cores catalog for the Taurus (Marsh et al., 2016) molecular cloud to extract the location of the cores at different evolutionary stages (as discussed in Section 2.1.5) namely: starless, prestellar and protostellar from the P lanck dust polarization maps. The derived cores catalog used the getsources identification algorithm (Men’Shchikov et al., 2012) to identify different star-forming cores within the molecular cloud by analyzing the spatial decompositions of the original images across a wide range of scales and waveband. The derived cores catalog also provides other useful parameters such as the major and minor FWHM of the core sizes (pc), estimated core mass (M⊙ ), dust temperature (K) from spectral energy distribution or SED fitting, and average column amongst others. 2.3. Method 27
2.2.4 JCMT Gould Belt Survey
For NGC1333, we used the “Perseus west” core catalog (Kate Pattle, private com- munication) to extract 57 cores from the POL-2 maps (Doi et al., 2020). The “Perseus west” catalog was made from the JCMT Gould Belt Legacy Survey which used three instruments, SCUBA-2, HARP, and POL-2 to survey nearby (< 500 pc) star-forming regions. The survey covered Gould’s Belt clouds like Orion, Tau- rus, Perseus, Serpens, Ophiuchus, and less well-studied local star-forming regions, to complement similar Gould Belt surveys with the Spitzer and Herschel space telescopes. The SCUBA-2 targeted half-degree regions based on Herschel 500 µm emission to provide a deep 450/850 µm survey. Overall JCMT has a better resolution than Herschel which is important in our case since NGC 1333 is more distant and crowded than Taurus.
2.3 Method
As we can see in Figure 2.7 the polarization intensity (P ) maps centered on each individual cores of NGC1333 are not significant enough to study the differences in polarization levels at the cores. This can also be seen in the cores of Taurus. Therefore by stacking the polarization maps for various core types as illustrated in Figure 2.8, we investigate whether the polarization intensities are different for the various evolutionary stages of the cores for both NGC1333 and Taurus. Further for NGC1333 we also investigate if stacking cores with signs of protostellar activity (brighter in 70 µm) shows us a significantly higher polarization when compared to ones that do not. Overall stacking helps to improve the S/N, thereby allowing us to 2.3. Method 28
examine the statistical properties of polarized emission of cores.
2.3.1 P lanck data for Taurus Molecular Cloud
Figure 2.7: Polarized intensity (P ) maps of individual cores of NGC1333 extracted
using the JCMT Gould Belt Survey catalog and POL-2 data. Each of the individual “postage stamp” maps of the cores are 12 × 12 pixel (∼ 0.24 × 0.24 pc) wide and the core centers have been marked by white cross-hairs.
Firstly we extracted the position of the individual core types from the Herschel Gould Belt Survey derived core catalog in right ascension (RA) and declination (DEC) and converted them to skycoord objects. 2.3. Method 29
Figure 2.8: Stacking analysis performed on protostellar cores in the Taurus molecular cloud where we take the average of all the individual cores for Stokes I, Q, and U maps.
Next, the standard Flexible Image Transport System (FITS) based maps of the respective molecular clouds were imported using the astropy FITS function along with their World Coordinate System (WCS) object.WCS describes the geometric transformations between one set of coordinates and another. We smoothed the P lanck maps of the Taurus molecular cloud to decrease the noise from the data by applying a 2-D Gaussian Kernel with FWHM of 4 pixels where each pixel corresponds to an angular size of 1.5’. The resulting resolution of the map is 7′ . 2.3. Method 30
The individual positions for the different cores were converted into map pixel values. This was done to a) extract the required information (column density, dust intensity, etc) from maps and b) make small cutout images of individual cores which for the Taurus Molecular Cloud is 40 × 40 pixels (∼ 2.45 × 2.45 pc).
Figure 2.9: High pass filter applied on individual cores in the Taurus molecular cloud to filter out large-scale structures around the core. The first image shows the core along with the diffused emission from its filament, the second picture shows the image in the Fourier space, the third image shows the Fourier space mask being applied to filter out certain frequency scales, and finally the last image shows the filtered core.
Depending on the resolution of the maps and the pixel size, we made “postage stamp” maps of each of the individual cores for Stokes I, Q, and U using the cutout2d python function. Once the postage stamp maps of the cores were ready, we divided each of the cores by its brightest central Stokes I pixel value. This was done to ensure that the 2.3. Method 31
bright cores do not bias our stacking analysis.
Next, for the P lanck maps of the Taurus molecular cloud, we applied a high pass filter as shown in Figure 2.9 to enhance the area around the cores and filter out the large-scale structures like filaments. This was achieved by converting the image into Fourier space using the scipy.fft function, applying a circular mask to block all the low frequencies in the spectrum, and finally taking an inverse Fourier transform. The reason for choosing a circular mask over a rectangular one was to eliminate edge effects such as ringing in the Fourier space. We filtered structures ≥ 0.24 pc by this method.
Figure 2.10: 70 µm dust emission map of the individual cores of NGC1333. Cores with a bright 70 µm source at the center have been labeled in red. 2.3. Method 32
We then calculated the polarized fraction (p) and Polarized intensity (P ) using Eq. 2.9 and Eq. 2.10 before stacking them. Next, we stacked each of the maps of the individual cores from the I, Q, U , P , and p data taking the average and repeating for the different core types (starless, prestellar, protostellar).
Figure 2.11: 2D Gaussian fit of the total intensity (I) maps of the three core cat- egories in NGC1333: cores with an isolated 70 µm source (left), cores without 70 µm source (middle) and, “uncategorized” (right), extracted from the JCMT POL-2 data.
2.3.2 JCMT POL-2 data for NGC1333
For our POL-2 map of NGC1333, we did not use the Herschel Gould Belt Survey core catalog as the POL-2 data has a higher resolution than the long wavelength Herschel maps. We used the JCMT Gould Belt survey core catalog as discussed in Section 2.2.4. We expect cores with compact 70 µm emission to harbor a luminous internal source, which may mean that they are forming a protostar (Dunham et al., 2008). Therefore we binned them based on the 70 µm dust emission map (shown in Figure 2.10) into three categories: cores with an isolated 70 µm source (labeled 2.3. Method 33
in red in Figure 2.10), cores without 70 µm source and those where we can not be sure if there is 70 µm source associated with the core (from now “uncategorized”) after making the 12 × 12 pixel (∼ 0.24 × 0.24 pc) “postage stamp” maps of each core. NGC1333 has a compact HII region which means that it has extended 70 µm emission where the radiation from the young star has heated the surrounding dust. We fitted a 2-D Gaussian to our stacked total intensity (I) maps of these different core types to extract their respective pixel sizes as shown in Figure 2.11. Using the 2D Gaussian fit parameters we performed aperture photometry on the stacked polarized intensity maps (P ) maps of these different core types using the photutils python library. This was done to compare the three different core categories. We constructed a circular aperture with a 0.05 pc (∼ core size) radius around the center of the stacked polarized intensity maps (P ) and an annular aperture for the background ranging from an inner radius of 0.08 pc to an outer radius of 0.12 pc from the center. For each of the categories, aperture statistics were performed to find the mean, median, standard deviation, and sum for both the central aperture and the background aperture. Only the exact intersection of the aperture with each pixel was calculated. The total number of pixels within the central aperture was determined by dividing the mean by the sum of the aperture which resulted in an average of 13 pixels for each of the categories. Next, we calculated the background-subtracted flux f ′ which can be defined as. f ′ = ΣIa − I¯b Na , (2.13)
where ΣIa is the net flux from the aperture, Na is the number of pixels in the aperture, and I¯b is the mean flux in the background annular aperture. The error on 2.4. Results 34
this background-subtracted flux f ′ can be approximated as,
p σf ′ = Na σp,b , (2.14)
where σp,b is the standard deviation of each pixel in the background calculated from aperture statistics. The results of this aperture photometry are shown in Table 2.1.
2.4 Results
2.4.1 Taurus Molecular Cloud
The Taurus molecular cloud was analyzed using the P lanck 353 GHz Stokes I, Q, U maps with ∼0.2 pc resolution. We analyzed over 435 starless, 52 prestellar, 38 candidate prestellar, and 11 protostellar cores using our stacking analysis. On the stacked I, Q, and U maps (Figure 2.12), we see that the polarization intensity levels increase as we go from starless to protostellar cores. The excess intensity for the protostellar cores, however, could be due to stacking elongated filaments within which the cores are embedded as shown in Figure 2.8. We also made stacked maps of the polarized intensity (P ) and polarized fraction (p) (Figure 2.13). Again, we see a bright source at the center of the polarized intensity map indicating extra polarization at the source. However, at the ∼0.2 pc spatial scale, we cannot resolve individual cores, and these signals are likely coming from larger clumps or long filamentary structures within which the cores are embedded. For this reason, a higher-resolution map is necessary. To illustrate this further we plotted all the individual postage-sized maps of the 2.4. Results 35
Figure 2.12: Stacked I, Q, and U maps of the different types of cores in the Taurus molecular cloud.
prestellar cores in Figure 2.14. This clearly shows us multiple similar filamentary structures are being repeatedly added in our stacking analysis, thereby contributing largely to our signal. Also, the cores appear offset from the center in Figure 2.14 as the Herschel catalog has identified a faint core, which is not apparent in the low-resolution Planck maps, while the visible “core” is a bright adjacent clump.
2.4.2 NGC1333
The higher 14.4′′ resolution helps us resolve the cores down to ∼0.02 pc, which significantly improves the S/N for our analysis. As discussed in Section 2.3 we performed aperture photometry on the three categories of core types: those with an 2.4. Results 36
Figure 2.13: Stacked maps of polarized intensity, polarized fraction, and distribution of the polarized percentage for the different core types in the Taurus molecular cloud. The individual images have been normalized to the peak Stokes I value in the I map.
isolated 70 µm source, without a 70 µm source, and “uncategorized” based on the 70
µm dust intensity map as shown in Figure 2.17. The statistics from the photometry have been summarized in Table 2.1. Table 2.1 shows us that the cores with a bright 70 µm source have a higher overall polarized emission compared to the polarized emission in the background aperture, in contrast to the cores without a bright 70 µm source. We can also see that the background subtracted flux for the cores with a 70 µm source is higher than the ones without a source. This statistical comparison however doesn’t yet apply to the cores which are “uncategorized”. The “uncategorized” cores do not have a central bright 70 µm source but rather have diffused structures spread across our 12 × 12 pixels (∼ 0.24 pc ×0.24 pc) postage stamp maps. The overall higher polarization levels for 2.5. Discussion 37
Figure 2.14: Individual filtered Stokes I maps of Taurus prestellar cores extracted using the Herschel derived cores catalog showing the filamentary structures.
the various aperture statistics could indicate that the extended structures are likely heated dust from nearby protostars and HII regions.
2.5 Discussion
With a jump in resolution from P lanck (5′ ) to JCMT (14.4′′ ), we have improved the resolution of our stacking analysis. The photometry results from the POL-2 data (Table 2.1), show us that cores with 70 µm source in the stacked polarized intensity maps (P ) show a higher polarization emission when compared to the cores without 70 µm source. This indicates that cores with protostellar activities are more polarized and thus perhaps the dust grains are more efficiently aligned. The “uncategorized” cores have much higher background and source polarization, which could indicate that the polarization tends to be stronger in areas where there is external heating leading to extended 70-µm emission as seen in figure 2.15. However, we can’t yet 2.5. Discussion 38
Figure 2.15: Herschel 70 µm dust emission map of NGC1333 from the Herschel Gould Belt Survey Catalog and the POL-2 map of NGC1333. The green circle represents the cores.
say we find evidence of significant polarized emission coming from cores without a protostellar source. A potential problem with our analysis could be that the compact polarization signals that we see arise from stacked filaments with different orientations, and we aren’t detecting polarization from the cores at all. This needs further investigation. The future goal of this work is to apply this analysis at random locations with sim- ilar column densities around the cores and find an even better proxy to categorize cores with protostellar activities like stellar outflows. The “uncategorized” section of cores could be analyzed using Variable-Pixel Linear Reconstruction (or “drizzling”) algorithm (Fruchter and Hook, 2002), which could preserve the photometry and res- 2.5. Discussion 39
Figure 2.16: Stacked I, Q, U , P , and p maps of NGC1333 cores
Figure 2.17: Aperture photometry performed on the stacked polarization intensity
(P ) maps of the three core types: with 70 µm source, without 70 µm source, and “uncategorized”. The white circle represents the circular aperture with a 2.35-pixel radius (∼ 0.05 PC) and the red circles represent the annular aperture for the back- ground with an inner and outer radius of 4 pixels (0.08 pc) and 6 pixels (0.12 pc).
olution. “Drizzling” also weighs input images according to the statistical significance of each pixel, and removes the effects of geometric distortion on both image shape and photometry, thereby helping us analyze the off-centered cores. As of now, we have only applied this analysis over 57 cores in NGC1333 which is not a sufficiently large statistical sample. We would need to run this analysis over other molecular clouds to have a significant statistical sample of different types of cores. The TolTEC camera on the LMT with its upcoming Fields in Filament survey at 5′′ resolution will be a crucial data set for us to analyze as it will cover 7 dense 2.5. Discussion 40
Table 2.1: Summary of the statistics of aperture photometry performed on the
stacked polarized intensity maps of NGC1333 of the three core types based on 70 µm dust intensity map classification: with an isolated 70 µm source, without 70 µm source, and “uncategorized”. Note that the units are arbitrary since each map is normalized relative to the brightest central pixel in the I map.
filaments which contain many cores. Also, the upcoming CCAT observatory with its fast mapping speeds and >20 times the resolution of P lanck will be important for this study. Hence, in the upcoming chapter, we look into the polarization-sensitive Prime-Cam instrument of the CCAT Observatory. Chapter 3
Atmospheric Characterization of the CCAT
Observatory site
In Chapter 2 we saw that there is a clear need for more high-resolution data to establish how the polarization levels change through various evolutionary stages of star-forming cores. The upcoming CCAT Observatory along with its polarization- sensitive Prime-Cam instrument could provide the data needed to study a statisti- cally significant number of cores. Therefore in this chapter, we introduce the CCAT Observatory, discuss the Prime-Cam instrument along with the MKID detectors, and finally characterize the atmosphere above the CCAT site to establish the need for tone-tracking firmware to efficiently read out the detectors.
3.1 What is the CCAT?
The CCAT is an upcoming observatory having the 6-meter mm/sub-mm Fred Young Submillimeter Telescope (FYST) at 5600 meters above sea level on the Cerro Cha- jnantor plateau. The construction phase is expected to lead to first light in 2024.
41 3.1. What is the CCAT? 42
It uses a novel “Cross-Dragone” design which provides high throughput (end-to-end
effective area divided by the geometric area of an unobstructed 6-meter aperture) and a wide field of view that will be capable of illuminating more than 100,000 latest- generation MKID detectors enabling rapid mapping of large areas of the sky (Stacey et al., 2015). CCAT will be able to observe the sub-mm radiation coming from cold dark molec- ular clouds and dusty regions where stars and planets form. These regions are not visible at optical wavelengths as they have high dust extinction. As we saw in Chapter 2, star formation is a multiscale process where dense cores form inside filamentary structures in molecular clouds. The FYST telescope will have twice the angular resolution of the Herschel Space Observatory and four times that of SOFIA, an airborne telescope. This ability of CCAT will enable it to study molecular clouds in our Milky Way and the surrounding Large Magellanic Cloud (LMC) and Small Mag- ellanic Cloud (SMC) with incredible depth, revealing the star-forming dense cores as the spectrum of thermal dust emission from molecular cloud peaks in the sub-mm. Also, with thousands to tens of thousands of detectors per band, the Prime-Cam in- strument will have a much greater mapping sensitivity than any other ground-based sub-mm telescope. CCAT will have three main instruments within it to take advantage of its novel optics and location, namely:
1. Prime-Cam: A wide-field submillimeter to millimeter wavelength camera (Vav-
agiakis et al., 2018)
2. CHAI: The CCAT Heterodyne Array Instrument (Barrueto et al., 2022)
Table 3.1: Prime-Cam band properties such as central frequency, bandwidth, beam size, and total number of MKID detectors(Aravena et al., 2022).
3. P-Spec: An imaging Spectrometer (Nikola et al., 2022)
Figure 3.1: A cross-sectional view of the modular cameras on Prime-Cam centered
at different wavelengths, E.M. Vavagiakis (used with permission from the author).
This Chapter is focused on the Prime-Cam instrument. Prime-Cam will have
seven sub-camera modules as referred in Table 3.1 that will be mounted within a 1.8-meter diameter cryostat cooled by a dilution refrigerator as shown in Figure 3.1. Compared to the 10,000 bolometer detectors on the JCMT SCUBA-2/POL-2 instru- 3.2. What are MKIDs? 44
ment (Friberg et al., 2016) the 850 GHz module will have 20,808 MKIDs, provid- ing unprecedented broadband intensity and polarization sensitivity for ground-based sub-mm/mm camera. Prime-Cam has ∼ 1.5λ/D pixel spacing and is diffraction- limited at all wavelengths. Scaled to a 6m telescope, the Prime-Cam 850 GHz will initially deliver a 1◦ field of view for each sub-camera with an overall field of view equivalent to ∼3◦ in diameter with a resolution of 15′′ . The 850GHz band will have about the same resolution as the 350 GHz band of SCUBA-2, but it will make much larger maps and therefore cover entire molecular clouds which will include over a thousand cores (Aravena et al., 2022). The 850 GHz module is funded by the Canada Foundation for Innovation (CFI) and will be built by scientists from the University of British Columbia.
3.2 What are MKIDs?
MKIDs are a type of superconducting detector based on the phenomenon of kinetic
inductance, which refers to the change in the inductance of a superconductor when its electrons are set into motion (Mauskopf, 2018). In an MKID, the superconducting material is typically made from a thin film of a high-transition-temperature superconductor, such as niobium nitride (NbN) or titanium nitride (TiN). The film is patterned into a series of resonant elements, typically small loops as shown in Figure 3.2. These resonant elements are then coupled to a microwave transmission line, which their response to incoming radiation to be measured. The resonant frequency of an MKID can be defined as, 3.2. What are MKIDs? 45
Figure 3.2: a) Schematic of a single kinetic inductance detector capacitively coupled
where C capacitance, Lg , and Lk are the geometric and kinetic inductance of
the resonator. The kinetic inductance Lk depends on the number density of Cooper pairs. When the KID absorbs a photon, it breaks a Cooper pair of electrons in the superconductor thereby changing its number density. These free electrons then move through the material, creating a cloud of charge carriers with a net motion that changes the inductance of the resonant element, thereby shifting the resonance frequency (f0 ) of the detector (Figure 3.3). The change in f0 is proportional to the power absorbed (Pabs ) by the MKID (Mauskopf, 2018),
1 df0 ∝ −Pabs . (3.2) f0 dPabs
This change in frequency can be measured using just one transmission line, pro- 3.2. What are MKIDs? 46
viding a highly sensitive measurement of the absorbed photon. MKIDs are a rel- atively new technology for mm and sub-mm ground-based cameras (e.g. Balloon- borne Large Aperture Submillimeter Telescope (BLAST)) and so far no operational ground-based sub-mm telescope uses MKIDs.
Figure 3.3: Basic operation of an MKID (a) Photons are absorbed, breaking Cooper pairs in the superconducting film and changing the inductance. (b) The film is placed in a resonant circuit. (c) The amplitude and (d) phase change of the resonance due to inductance change are measured, indicating the photon energy, Day et al. (2003)(used with permission from the author).
To read out MKID arrays, a comb of tones of resonant frequencies is transmitted
along a single transmission line connected to every detector. By tuning each element of the array to a slightly different resonant frequency, the phase and magnitude of 3.2. What are MKIDs? 47
each detector can be individually read out from the same signal line. This is a huge advantage compared to traditional bolometers which require multiple wires for each detector, and Transition Edge Sensor (TES) which use time-multiplexing rather than frequency multiplexing. Therefore, we can read out tens of thousands of these detectors with compact and low-power readout electronics. The major issue with MKIDs arises when they are subjected to large changes to the incident power. As the incident power loading on the detector changes due to the change in the column of atmosphere above, it saturates the MKIDs (Figure 3.4), making it difficult to detect the small changes in the probe tones.
Figure 3.4: Blackbody power sweep on a test MKID resonator with a temperature that varies from 3 K to 35 K corresponding to an optical power of 0 to 48 pW. The plot on the left illustrates the S21 sweep of a resonator at various optical loads. The plot on the right shows the decrease in the quality factor (Qi) as the incident power on the detector increases, Chapman et al. (2022) (used with permission from the author).
As we can see in Figure 3.4 there is a diminishing return from the MKIDs as incident power increases, making them difficult to read out. This could be a problem for Prime-Cam. Despite the excellent altitude and the site’s dry conditions, the 3.3. Atmospheric simulations 48
MKID arrays’ high responsivity is highly susceptible to saturating the MKIDs during its scans, especially at the 850 GHz band. In this chapter, we will describe a thorough atmospheric characterization of the CCAT site using the am Atmospheric Model Paine (2012) as the atmosphere will be the dominant noise source for the MKIDs.
3.3 Atmospheric simulations
Atmospheric modeling is important when designing instruments and planning as-
tronomical surveys because it helps us understand how the atmosphere affects our observations. Based on the simulation results, we would want to test whether changes in atmospheric loading will push the MKID detectors out of their linear frequency response regime and therefore require more complicated instrumentation or software to track the frequency shifts from changing power levels (tone tracking). Atmospheric Model (am) (Paine, 2012) is a radiative transfer tool for computa- tions of propagation paths at the microwave to submillimeter wavelengths. Am models the atmosphere as a sequence of planar layers as shown in Figure 3.5 and calculates parameters like thermal emissions, absorption, optical depth (τ ), brightness (Tb ), and Rayleigh-Jeans temperatures (TR−J ). The radiation, initially with a background temperature of T0 , moves through a succession of layers. These layers correspond to the atmosphere’s horizontal strata and generally denote a section of a propagation path. A single layer, for instance, might represent a horizontal path between two sites in the atmosphere at a constant temperature and pressure. Each layer contains a mixture of various absorbing species, such as dry air, and has a defined column density. 3.3. Atmospheric simulations 49
Figure 3.5: Planar layers of the atmospheric model showing how the radiation propa- gates through the individual layers with defined height (h) and midpoint parameters like pressure (P ), temperature (T ), and column density(N ) (Paine, 2012).
The incoming radiation propagates through the layers at a fixed zenith angle θ with negligible refraction. The zenith angle can be varied between simulations, which will further decide the path length of the incoming radiation. For models in am where the layer temperature is set by the layer midpoint tem- perature Ti , the Planck radiation function,
2hν 3 B(ν, T ) = , (3.3) c2 (ehν/kT − 1)
is computed for constant T = Ti , throughout the layer. In this case, the radiative 3.3. Atmospheric simulations 50
transfer across a layer takes the form
Ii (ν) = B(ν, Ti )(1 − e−τi (ν) ) + Ii−1 (ν)e−τi (ν) , (3.4)
where Ii (ν) is the radiance at the lower boundary of the ith layer, and τi (ν) is the optical depth of the layer. On the other hand, for models with temperatures defined at the layer boundaries, the Planck function is assumed to vary linearly with optical depth across the layer,
The optically thick limit of these expressions is simply,
pi (τ ) → 1. (3.7)
An iterative method is applied to compute the radiative transfer for a complete
stack of layers, starting at the top of the stack where the initial radiance spectrum is
I0 (ν) = B(ν, T0 ). (3.8)
3.3. Atmospheric simulations 51
Eq. 3.4 or Eq. 3.5 is applied iteratively to produce the model output radiance:
I(ν) = In (ν). (3.9)
The total opacity,
X τ (ν) = τi (ν), (3.10) i
through the stack of layers is computed along with the radiance.
The default direction of propagation in am is from the top layer (1) to the base layer (n). The radiance I(ν) and the optical depth τ (ν) are the fundamental spectral quantities computed by the radiative transfer model. The Brightness Temperature, Tb (ν), is obtained by inverting Planck’s radiation function as stated in Equation 3.3 and solving for T . Substituting the radiance I(ν) for the Planck radiance B(ν, T ) we get,
hν Tb (ν) = 3 . (3.11) k ln(1 + c2hν 2 I(ν) )
The Rayleigh-Jeans temperature, on the other hand, is calculated by taking the limit of Eq. 3.3 for hν << kT and inverting to obtain
c2 TR−J (ν) = I(ν). (3.12) 2kν 2
Hence, the Rayleigh-Jeans temperature is the low-frequency limit of the bright-
ness temperature Tb (ν). TR−J has an important physical significance; when it is multiplied by the Boltzmann constant kTR−J (ν), we get the power absorbed by an 3.4. Procedure 52
antenna per unit bandwidth coupled to a single transmission line. To calculate the power loading (in Watts) on the detector we multiply this by the bandwidth ∆ν :
∆P = kTR−J (ν)∆ν. (3.13)
3.4 Procedure
Every computation with am is defined by a configuration file containing a description
of the model atmosphere or propagation path and other parameters that control the model computation or fit. This file contains replaceable parameters labeled %1, %2, etc., which can be replaced at the am command line based on the required amount of water vapor profile specific to the site. The am configuration files can be passed to am by their respective name at the command line, like,
$ am model.amc
am computes the files and writes two output streams:
1. Standard output stream: which stores the output spectra
2. Standard error stream: which stores the model configuration data
The standard output stream is then written out to a “.out” file which can be named in the command line,
$ am model.amc > model.out
The structure of a single layer am configuration is shown in Figure 3.6, This
3.4. Procedure 53
Figure 3.6: Typical am file format.
example here shows a laboratory setup in which a spectrometer views a 77K liquid nitrogen cold load, through a meter of air, at 1 atm ambient pressure and 20◦ C ambient temperature. This file is configured to compute transmission from 0 GHz to 2500 GHz. Depending on the site characteristics, we can construct configuration files and modify existing ones based on our needs. We can also specify the parameters we want to compute and stream in the output file. Since we are interested in characterizing the atmosphere above the CCAT observatory, we modified the ALMA site files because of its proximity to the CCAT site at the peak of the Cerro Chajnantor plateau. The ALMA site file has been taken from the am “cookbook” which contains configuration files for high-altitude ground-based mm/sub-mm telescopes. To match the conditions at the Cerro Chajnantor peak, the base layer from the ALMA configuration file at 554 mbar was removed. The base layer for the CCAT now starts at 550 mbar, which is the mean atmospheric pressure at the CCAT Observatory taken from the observatory’s meteorological archival data (2006-2014)1 . 1 CCAT Observatory meteorological archival data (2006-2014) 3.4. Procedure 54
This pressure difference makes sense as the CCAT site is 600 m above the ALMA site.
Figure 3.7: Modified am ALMA site file where the base ALMA layer was removed; thus, the modified site file starts from 550 mbar to 0.1 mbar.
Each site file for high-altitude ground-based mm/sub-mm telescopes can be found under the am cookbook after installing the software2 . There are individual site files for each month along with their annual average for different H2 O percentile profiles. H2 O percentiles represent the PWV content in the atmospheric layers. Firstly, we compared our model files against the unmodified Atacama Cosmology Telescope (ACT) site file, broadly used within the rest of the CCAT collaboration (Choi et al., 2020). ACT is a 6m CMB mapping telescope on Cerro Toco (5100m) 2 am software 3.4. Procedure 55
and is within close proximity to the CCAT site. Our modified ALMA files matched closely with the ACT files (Figure 3.8).
Figure 3.8: Comparison between the ACT and our modified ALMA files for at- mospheric transmission and Rayleigh-Jeans temperature as a function of elevation angle.
Secondly, we simulated the overall site characteristics for a frequency range of 0
to 900 GHz, covering all three CCAT bands (Figure 3.1). This includes finding the frequency vs. atmospheric opacity for given relative humidity and Liquid Water Pro- file (LWP) value, which tells us if our bands have any bright H2 O or O2 atmospheric absorption lines. These parameters are passed during the simulation. Next, we compute the annual and monthly atmospheric transmission, Rayleigh- Jeans temperature, and brightness temperature as a function of frequency for various H2 O percentiles. These results will show us which months/time of the year is best suited for observations, and also help us spot how the sky brightness and Rayleigh- Jeans temperature varies throughout the year. Next, we perform a detector band-wise (Figure 3.1) “skydip” analysis showing 3.4. Procedure 56
how the various parameters like atmospheric transmission, brightness, and Rayleigh- Jeans temperature and power loading on the detector vary with elevation angle. Finally, we calculate the actual frequency shifts and power loading on the detector as a function of the elevation angle. This means we can find the change in incident power on the detector for a certain change in elevation angle. This is important as it allows us to test what detector responsivities are best suited to the range of power loads that will be incident on the detectors. The power loading on the detector can be estimated as,
∆P = TR−J k ∆ν ηoptical , (3.14)
where TR−J is the Rayleigh-Jeans Temperature (in K), k is the Boltzmann Constant (1.380649 × 10−23 m2 kg s−2 K −1 ), ∆ν is the bandwidth of the array, and ηoptical is the optical efficiency of the telescope which refers to the amount of net energy which reaches the detectors after passing through the optics of the telescope. For a given change in power load ∆P the actual frequency shift in Hz is,
∆f = R fres ∆P, (3.15)
where R is the responsivity of the detector array and fres is the resonant frequency of the MKID detector readout tone (assumed to be 800 MHz (Lourie et al., 2018)). 3.5. Results 57
3.5 Results
In this section, we will analyze the results of our simulations. For the full set of simulation results over all possible observing months and many different observing conditions see Appendix A.
3.5.1 Overall CCAT Observatory site characteristics
Figure 3.9: Frequency vs opacity for different H2 O profile percentiles. The shaded bands represent the Prime-Cam frequency bands.
Figure 3.9 shows the frequency vs. opacity for 33.67% Relative Humidity (RH) 3.5. Results 58
and 100 g m−2 LWP at the CCAT site for different H2 0 profiles. The RH and the LWP values are obtained from the 2004-2016 meteorological data. This plot shows us if the frequency bands of the Prime-Cam detector contain any bright molecular line; the H2 O rotational lines show up at 570 GHz and 750 GHz. There is no atmospheric contamination at our frequency bands by any molecular species.
Figure 3.10: Atmospheric transmission vs. Frequency for different months at the CCAT site. The vertical lines, along with the shading, represent the band centers along with the respective bandwidth.
Figure 3.10 shows annual atmospheric transmission vs. frequency for the 50th percentile water profile at the CCAT site. A value of 1 represents 100 % transmission of the radio waves or a transparent atmosphere, whereas 0 represents no transmission hence an opaque atmosphere. As we can see, the atmosphere is highly transparent for the 280 and 350 GHz bands during the winter months (July to November). However, the transmission value falls off with increasing frequency. The 850 GHz array occupies the last atmospheric window where the transmission is at most ∼45%. It is also clear that the transmission is significantly better for the colder and drier Chilean winter 3.5. Results 59
months from July to November. Hence, these months are most suitable for high- frequency observations, especially for high-frequency 850 GHz observations. The atmospheric transmission deteriorates with higher H2 O percentiles. For transmission simulations of the different H2 O profile percentiles see Appendix A.
3.5.2 Detector band-wise skydip analysis
Next, we look at the am simulations for detector band-wise skydip analysis, where we simulate the change in parameters like atmospheric transmission, opacity, brightness, and Rayleigh-Jeans temperature as a function of elevation angle of the Prime-Cam instrument for various H2 O profiles and months. As mentioned before, we look at a skydip of 10◦ to 85◦ with 5◦ steps in elevation angle. Figure 3.11 shows us that the atmospheric opacity decreases as elevation in- creases, which is understandable as the path length through each of the layers gets shorter as compared to a lower elevation angle with a longer path length. Overall the atmospheric opacity also increases as a function of frequency as we can see the 850 GHz band is almost a magnitude higher than the 280 GHz band. Next, we look at the Rayleigh-Jeans temperature (TR−J ) for the various frequency bands as a function of elevation. Rayleigh-Jeans temperatures (Figure 3.12) are significantly lower for the 280 GHz frequency band and it increases as we go higher up in frequency. However, it seems to decrease as the elevation angle increases. Also, for the 850 GHz frequency band, a lower H2 O profile results in much lower Rayleigh-Jeans temperature values. There is almost an order of magnitude increase in Rayleigh-Jeans temperature from 250 3.5. Results 60
Figure 3.11: Skydip analysis of atmospheric opacity for different elevation angles. Each panel shows the annual average of the three different frequency bands of CCAT
GHz to 850 GHz.
3.5.3 Power loading on the detector and actual frequency shifts
The total incident power on a detector is a sum of a static component (ptele (static)) from the thermal radiation of the telescope’s optics and a dynamic component (Patm (dynamic)) due to the atmosphere, which changes as a function of elevation angle, weather conditions, and gas pockets with different humidity drifting past the telescope, 3.5. Results 61
Figure 3.12: Skydip analysis of Rayleigh-Jeans temperature (TR−J ) for different ele- vation angles. Each panel shows the three different frequency bands of CCAT.
Ptotal = Ptele (static) + Patm (dynamic). (3.16)
For our simulations, we assumed that ptele (static) = 0, as we are only concerned about how the change in power will affect the detector readout. As seen in Eq. 3.14 the atmospheric power directly depends on the Rayleigh-Jeans temperature (TR−J ), frequency bandwidth (∆ν), and the optical efficiency (ηoptical ) which we assume to be 10%. We expect the optical efficiency to be somewhere between 10-30%, so 10% is a conservative assumption. 3.5. Results 62
Figure 3.13: Skydip analysis of loading power of the atmosphere in (pW ) for different elevation angles. Each panel shows the annual average of three different frequency bands of CCAT.
As we can see in Figure 3.13, the power loading on the individual MKID detectors shows us that for a higher frequency 850 GHz array, the atmospheric power loading is about a magnitude higher than the lower frequency 280 GHz array. Again, the power loading values fall off as we look higher up in elevation as the telescope is observing through less of the atmosphere. Winter months result in slightly lower loading power and, thus, better observations. Now to calculate the actual shifts in the frequency of these MKID arrays due to changing loading power, we would need the responsivity characteristics of the 3.5. Results 63
detectors. Responsivity is the measure of the input-output gain of a detector system.
In other words, for a typical photodetector, it is a measure of electric power output per optical input. In simple terms,
q R=η , (3.17) hν
where η is the quantum efficiency of the detector for a given frequency, q is the electron charge, ν is the frequency of the incident photon, and h is Planck’s constant. Using Eq. 3.15, we can estimate the actual shift in frequency. The MKID resonator parameters have been taken from the BLAST-TNG 850 GHz (Lourie et al., 2018) array where Qr = 36380 (Quality factor), Qc = 58290 (Coupling factor), and, fres = 800 M Hz (Resonant Frequency).
Figure 3.14: Responsivity of the 280 GHz MKID detector, Jordan Wheeler (used with permission from the author).
We used the responsivity values as shown in Figure 3.14 for our analysis. For 3.5. Results 64
powers higher than 8.5 pW, we assumed a constant R = −1.5 × 108 .
Figure 3.15: Skydip analysis of the expected shift in frequency in MHz for different elevation angles. Each panel shows the annual average of three different frequency bands of CCAT
As we can see in Figure 3.15, the actual frequency shifts are significantly higher for the 850 GHz array. The frequency shifts are important because a large enough fre- quency shift can throw the MKID arrays from their linear regime, thereby requiring us to constantly re-tune the arrays. The linear regime is an intrinsic property of a resonator wherein the phase and frequency are linearly related as shown in Figure 3.16. This is crucial for tone- tracking as this conversion factor allows for an easy change from phase feedback to 3.5. Results 65
Figure 3.16: Linear regime of a resonator showing how the phase and frequency are linearly related. The solid lines show the phase of the resonator in both rest and modulated state and the dashed lines represent the slope of the two curves, Hoh et al. (2020)(used with permission from the author).
frequency feedback. The linear regime, especially for the sensitive 850 GHz band, can be estimated by looking into the complex forward transmission of the resonator or its S21 characteristic, Qr 1 S21 = 1 − . (3.18) Qc (1 + j2Qr x)
Where Qr and Qc are the total and coupling quality factors, respectively, x is the fractional frequency shift and j is the imaginary unit. The FWHM of a detector with the adapted properties of the CCAT MKIDs is:
fres 800 F W HM = = = 22, 000 Hz, (3.19) Qr 36380
and (1/10)th of line width (Figure 3.16) or FWHM is estimated as the linear regime 3.5. Results 66
for the MKID detector array, beyond which the phase and the frequency are not linearly related. With all of these values, we can calculate how much change in elevation would throw off the MKIDs from their linear regime, when the frequency changes by more than 2,200Hz.
Figure 3.17: Possible changes in elevation before being thrown off the 2.2 KHz linear regime. Each panel shows the annual average of three different frequency bands of CCAT
Figure 3.17 shows us that for the 280 and 350 GHz we can move up in elevation to a few degrees before getting thrown off the linear regime and retuning the MKIDs. However, for the 850 GHz array, with the 280 GHz responsivity values, it is pretty much impossible to stay within the linear regime even for 0.5◦ changes in elevation. 3.5. Results 67
In Eq. 3.15 we saw that the actual frequency shifts are related to the responsivity of the MKID detector array. Therefore, using a less responsive MKID array can significantly push down the power loading on the detector array and maintain the linear regime as shown in Figure 3.18. For detector band-wise “skydip” simulations varying over different months see Appendix A.
Figure 3.18: Possible changes in elevation before being thrown off the 2.2 KHz linear regime. Panel a shows the results from using the 280 GHz responsivity values, and panel b shows the results from using a MKID detector with lower responsivity (a factor of 10 less).
Overall, it will be much harder to keep the high-frequency detectors in the linear regime for scans at low elevation, compared to high elevation. For the 850 GHz array, we would need a tone-tracking algorithm to ensure that the MKIDs are always within their linear operation regime. The plan is to develop a Proportional-Integral- Derivative (PID) based FPGA tone-tracking firmware to achieve this in the near future. In this chapter, we introduced the CCAT observatory along with the polarisation- sensitive Prime-Cam instrument, discussed the functioning of the latest generation 3.5. Results 68
microwave kinetic inductance detectors (MKIDs), and finally characterized the at- mosphere above the CCAT site to better plan observing strategies. We found that dry conditions of the Chilean winter months are suitable for observation, especially for the high-frequency 850 GHz band. We also found that low-elevation scans would require a tone-tracking algorithm as it will be hard to keep the high-frequency de- tectors (850 GHz) in their linear regime of operation. Chapter 4
RFSoC-based backend for Balloon-borne VLBI
Experiment (BVEX)
4.1 Introduction
The Balloon-borne VLBI Experiment or BVEX is an upcoming pathfinder experi-
ment set to launch in the summer of 2025 from Timmins, Ontario, which has the goal of demonstrating the feasibility of a high-frequency stratospheric interferometry station. VLBI, is a special interferometric technique used in radio astronomy to im- age distant objects with very high resolution using continental baselines. Using this technique, EHT was able to image the shadow of black holes M87* and SAG A*, achieving resolutions up to 25 µarcsecond (Ball et al., 2019). SMBHs with masses ranging from ∼ 106 to ∼ 108 M⊙ are found at the centers of almost all galaxies (Cattaneo et al., 2009). They play a crucial role in the overall evolution and star for- mation rate of the galaxy. Hence, it is very important to characterize more of these SMBHs and directly imaging their event horizon gives us a wealth of crucial infor- mation. However, huge distances to these black holes mean we need higher angular
69 4.2. What is VLBI? 70
resolution. As discussed in Chapter 1 the maximum resolution for radio interferom-
λ etry is θ = B where λ is the wavelength of the radiation and B (called baseline) is the separation between the individual telescopes. Therefore, we would either need baselines larger than VLBI already achieves, i.e., the Earth’s diameter, or push to even higher frequencies. The issue is that the Earth’s atmosphere severely impacts the incoming radio waves at higher frequencies (see Chapter 3). Hence, ground- based telescopes operate up to a certain frequency range before getting limited by atmospheric noise. Stratospheric balloons, on the other hand, operate above 99.5% outside the Earth’s atmosphere, thereby allowing us to push our observations to a higher frequency. Since balloons can be launched virtually anywhere on the Earth, these mobile interferometry stations will increase the coverage of the u − v (refer to Section 4.2.2) plane by filling in where we do not have a ground-based telescope allowing us to get a better-reconstructed image.
4.2 What is VLBI?
VLBI, or Very Long Baseline Interferometry, is a special radio interferometric tech-
nique consisting of baselines that extend up to several thousands of kilometers. Single-dish radio telescopes usually have a relatively low angular resolution as they are limited by their diameter. We would need sub-arcsec resolution to observe dis- tant astronomical radio sources, resulting in impractically large dish sizes. Large dishes also deal with issues like steerability, compromised source tracking accuracy, gravitational sagging, and maintaining required surface accuracy. To overcome all these issues and increase the angular resolution, we use multiple smaller dishes sep- 4.2. What is VLBI? 71
arated from each other by some distance to simulate one large dish effectively. This is called interferometry.
4.2.1 Radio interferometry: brief overview
Figure 4.1: A simple two-element multiplying interferometer observing a radio source
s centered on ν = ω/(2π)s.
Figure 4.1 illustrates the simplest form of a radio interferometer (N = 2) which
is pair of two radio telescopes separated by some distance, whose output voltages are correlated (multiplied and then averaged). Even for N ≫ 2, there will always be N (N − 1)/2 independent pairs of such interferometers. ⃗b is the baseline vector representing the distance between the two dishes, ŝ is a unit vector that represents the pointing direction of both the dishes, and θ is the angle such that plane waves from a distant point source travel an extra distance ⃗b.ŝ = b cos θ to reach antenna 1. 4.2. What is VLBI? 72
So, the output of antenna 1 lags in time by a geometric delay τg .
⃗b.ŝ τg = . (4.1) c
Assuming a quasi-monochromatic interferometer responding to a narrow frequency
band ∆ν ≫ 2π/τg centered on ν = ω/(2π). The raw output voltages from antennas 1 and 2 after being mixed down to a suitable bandwidth at time t can be written as,
V1 = V cos[ω(t − τg )], (4.2)
V2 = V cos(ωt), (4.3)
The correlator then multiplies these two voltages,
It then takes an average to remove the high-frequency term cos(2ωt − ωτg ) while keeping the slowly varying term. Finally, we have,
V2
R = ⟨V1 V2 ⟩ = cos(ωτg ) (4.5) 2
Where R is the time-averaged correlator response and the voltages V1 , and V2 are proportional to the electric field produced by the source multiplied by the voltage gain of the two antennae and their receivers. Therefore the correlator output amplitude V 2 /2 is proportional to the flux density s of the point source multiplied by (A1 A2 )1/2 , 4.2. What is VLBI? 73
where A1 and A2 are the effective collecting areas of the two antennas. The correlator response R varies sinusoidally with the Earth’s rotation; which changes the source direction relative to the baseline vector. These sinusoids are called fringes, and the fringe phase is,
ω ϕ = ωτg = b cos θ, (4.6) c
where ϕ depends on θ as follows,
dϕ ω b sin θ = − b sin θ = −2π . (4.7) dθ c λ
Detecting these fringes is one of BVEX’s primary science goals. The fringe period ∆ϕ = 2π corresponds to an angular shift ∆θ = λ/(b sin θ). The phase of the fringe is an extremely sensitive measure of the source position if the projected baseline b sin θ is a multiple of many wavelengths. The point-source response of an isotropic interferometer is a sinusoid covering the sky. An isotropic interferometer is sensitive to only one Fourier component of the brightness distribution with an angular period of λ/(b sin θ). To improve the point source response of an interferometer, we need more Fourier components: i.e., more baselines. For an interferometer with N antennas, we will have N (N − 1)/2 pairs of antennas, each of which would be a two-element interferometer. Hence, the effective point-source response will be an average of the outputs of all the two-element interferometers rapidly approaching a Gaussian as N increases. 4.2. What is VLBI? 74
4.2.2 Earth-Rotation Aperture Synthesis
The anticipated baseline coverage of an interferometer whose components are fixed
to the ground changes due to Earth’s rotation. As the Earth rotates, all baselines of an interferometer constrained to an east-west line will remain in a single plane per- pendicular to the North-South rotation axis. The computational benefit of limiting all baselines to two dimensions is that the brightness distribution of a source may be computed as the two-dimensional Fourier transform of the observed visibilities. This cartoon illustrates Earth’s rotation aperture synthesis by a two-element inter-
Figure 4.2: Cartoon illustrating Earth’s rotation aperture synthesis by a two-element
interferometer and its corresponding u − v plane.
ferometer aligned with the east-west direction at a latitude of +40◦ as viewed from a source at declination δ = +30◦ . Let us assume u to be an east-west component of the projected baseline and v as the north-south component of the projected baseline in wavelengths, so, during 12 hours centered on the source transit, the interferometer 4.3. Key objectives of BVEX 75
will trace out an ellipse on the (u, v) plane. The maximum value of u equals the separation of the antenna in wavelengths. The maximum value of v is smaller by the source declination factor sin δ. If the interferometer has more than two elements, the u − v coverage will become a series of concentric ellipses having the same shape. Therefore, the synthesized beam in the east-west Earth rotation aperture synthesis can approach an elliptical Gaussian with enough dishes.
4.3 Key objectives of BVEX
As mentioned earlier, BVEX will be a pathfinder experiment demonstrating that we
can track the changes in the position of our balloon telescope well within a fraction of the observing wavelength to perform VLBI with a ground-based telescope as shown in Figure 4.3. Also, since stratospheric balloons operate over 99.5% of the atmosphere, it will allow us to access higher frequencies and help us push the baselines and cover more of the u − v plane. The primary objectives of BVEX are:
1. To characterize the positional stability of balloon-borne telescope to a) at least
40 µm precision over a 0.1-second timescale, and b) at least 100 µm precision over a 1-second timescale;
2. To track the 3D position of the balloon-borne telescope to better than 1m
accuracy on 1-sec timescales;
3. To characterize the phase stability of the stratosphere and the vibrational en- vironment of a stratospheric balloon gondola; 4.3. Key objectives of BVEX 76
4. To simultaneously observe a bright radio source along with a large ground-
based dish to detect fringes and demonstrate VLBI.
To achieve all these objectives, we are building a balloon-borne 1-m K-band (21- 23 GHz) radio telescope that will fly on the Centre National d’Études Spatiales (CNES) gondola through the Canadian Space Agency (CSA) STRATOS program in the summer of 2025 from Timmins, Ontario. During its short ∼12 hour flight, BVEX will record the data from all its sensors like GPS, accelerometer, etc., along with raw radio signals from the radio dish, re-quantize them to 2 bits, find the spectrum, and store all of the data into disks so that we can reconstruct the position and correlate the signals with a ground-based dish. Ultimately BVEX will demonstrate if we can track position and phase well enough to design a higher-frequency (∼350 GHz) VLBI telescope that could partner with the ngEHT (Tiede et al., 2022). Balloon- borne stratospheric VLBI stations like BVEX and its upcoming high frequency (ν >200 GHz) successor will have the following advantages over ground-based stations:
1. Ability to push to higher frequencies: Figure 4.4 shows that pushing to 690 GHz will help us resolve twice as many black holes as with the current EHT standard.
2. Increased sensitivity of dish: As we know, for a two-element interferometer,
the sensitivity is defined as,
21/2 kTs σs = , (4.8) Ae (∆ντ )1/2
where Ts is the noise temperature of the system and Ae is the effective area 4.3. Key objectives of BVEX 77
Figure 4.3: Cartoon showing simultaneous VLBI observation between a balloon and a ground-based telescope.
of the dish. At high frequencies, Atmospheric noise is the biggest contributor
to Ts . Since balloon telescopes observe through < 1% of the atmosphere, the signals from the astrophysical sources are less attenuated, and the overall systemic noise is reduced by a factor of ∼100.
3. Good baseline sensitivity: Since balloons have limitations in terms of pay-
load volumes, therefore, we can only launch relatively small telescopes (<2.5 meters). The loss in the baseline sensitivity is offset by pairing with a large ground-based dish, as the geometric mean of the system-equivalent flux den- sity of each telescope sets the baseline sensitivity. Therefore, even a small dish (d∼2m) high-frequency balloon-borne telescope paired with a large ground- based station like the Green Bank Telescope (GBT), which has a high System 4.3. Key objectives of BVEX 78
Figure 4.4: 230-690 GHz flux vs. supermassive black hole shadow sizes with known masses. The dashed vertical lines show the fringe resolution of an Earth-sized inter- ferometer at frequencies between 86 GHz and 690 GHz, Johnson et al. (2019)(used with permission from the author).
Equivalent Flux Density (SEFD), will provide the necessary sensitivity on long baselines. The Root Mean Square (RMS) noise is given by,
where SEF D1 and SEF D2 are the system equivalent flux densities of the two telescopes and ηQ is the quantization efficiency. SEFD is generally defined as the flux density of a radio source that doubles the system temperature. Lower values of the SEFD indicate more sensitive performance.
4. Mobility: One of the biggest advantages of a balloon-borne VLBI station is the
4.4. Stages of BVEX 79
ability to launch the telescope virtually anywhere. As we saw, the ability to
detect interference fringes depends on the Fourier transform of the VLBI station locations from the point of view of the source, with the longest baseline setting the maximum resolution. Therefore, we can pick complementary locations compared to the fixed ground site with a balloon telescope, thereby giving us better coverage in the u − v plane. The balloon also moves giving even more baselines (although it doesn’t really move on a highly predictable trajectory).
4.4 Stages of BVEX
In this section, we will discuss and design the timing reference chain, stable operation of the OCXO, RFSoC-based backend readout, and a spectrometer.
Figure 4.5: Block diagram showing the stages of BVEX, which includes the dish, receiver box, mixing chain, and readout electronics (Credit: Felix Thiel). 4.4. Stages of BVEX 80
4.4.1 BVEX front end
The front end of BVEX consists of a 1-m Cassegrain dish (223K-36/.396/595), the feedhorn, receiver antenna, polarizer, and a coax adapter leading to the Low-Noise Amplifier (LNA) (ERZ-LNA-2100-2700-25-2) (low noise amplifier). The entire radio telescope and the LNA box sit on an elevation mount. The 1-m Cassegrain dish is
Figure 4.6: The BVEX front end includes the 1-m radio telescope, elevation mount, and the LNA box. The bottom panel shows the LNA box design. (Credits: Felix Thiel)
fiberglass reinforced with a lightweight aluminum frame weighing around 12 kg. The dish has a 45 dB peak gain and an effective collecting area of 0.445m2 at 22.5 GHz (A = Gλ2 /4π). The angular size of the beam is about 1◦ . The LNA box consists of two LNAs, each having a 25dB gain, a 2dB noise figure, 4.4. Stages of BVEX 81
and a 169.6 K noise temperature. A high-pass filter is placed between the two LNAs to let RF frequencies between 21-23 GHz pass while filtering out the rest. The RF signal from the LNA box between 21-23 GHz is too high to be digitized by an ADC; this signal needs to be mixed down to a suitable frequency to be processed further. We typically use an Radio Frequency (RF) mixer, a local oscillator, and a bandpass filter to achieve this.
Figure 4.7: Mixing stage of BVEX.
The BVEX mixing stage consists of an RF mixer (FMMX1030) running from 16 to
30 GHz, a 19 GHz 5.5 V DC local oscillator (PLDRO-19000-10), a 2-4 GHz bandpass filter (ZBSS-3G-S+), and an Intermediate Frequency (IF) low noise amplifier (ZX60- 63GLN+). The mixer has an IF frequency range from DC to 8 GHz and a net conversion loss of 9.5 dB. Our RF signals from 21-23 GHz are mixed along with a 19 GHz reference LO signal to get an intermediate frequency band from 2-4 GHz. The higher frequency RF+LO component running from 40-42 GHz and the RF image (lower sideband) component running from 15-17 GHz is rejected using the 2-4 GHz band pass filter (refer to Felix Thiel’s MSc thesis for more details on the receiver design). 4.5. Timing reference design and challenges 82
4.5 Timing reference design and challenges
In VLBI, signals emitted by distant radio sources are received and recorded at several antennas. At each antenna (VLBI station), a very stable frequency standard is required to provide a reference signal that enables time stamping of the signals as they are being recorded. For each VLBI observation, the correlation of the timestamped signal from the participating station yields the differences among the arrival times of the specific radio wave at the antennas. These time differences are used to calculate the locations of the antennas respect to each other. This is extremely important as it allows us to precisely calculate the different baselines, correlating the signals without degradation. VLBI also measures the phase differences between the signals received by different telescopes. The stability of timing is directly related to the stability of the phase measurements. Any fluctuations or instabilities in timing can introduce phase errors, degrading the quality of the interferometric data and affecting the final image quality. Phase Noise (L(f )) is a method of describing the stability of an oscillator in the frequency domain. It is defined as noise power in a 1 Hz bandwidth, at a particular frequency removed from the carrier,
P ower density in a single 1 Hz side band
L(f ) = . (4.10) T otal signal power
A quantity of special interest in VLBI is the coherence time (τc ), which for a 4.5. Timing reference design and challenges 83
phase error of 1 radian is defined as,
2πν0 τc σy (τc ) ≈ 1, (4.11)
where ν0 is the observing frequency and σy is the Allan deviation. Therefore, a
stable timing reference is essential with sufficiently low Allan Deviation (ADEV) to maintain accurate phase measurements and reduce phase errors. ADEV is a measure of frequency stability in clocks, oscillators, and amplifiers. For BVEX considering observations at 23 GHz (ν0 ) for 1 second integration time, the minimum required ADEV for a 1-radian phase error would be,
1 σy ≈ ≈ 7 × 10−12 . (4.12) 2πν0 τc
Usually, for ground-based VLBI stations, we use hydrogen masers to provide
stable timestamping, however, due to volume and power restrictions on a balloon, we cannot fly such a system. For BVEX, we therefore need to use something of a much smaller form factor that can give us the required stability. Oven-controlled crystal oscillators, or OCXOs, are a great way of achieving the required stability while being just under 250 gm. OCXOs are temperature-controlled devices that maintain a constant operating temperature for a quartz crystal. This prevents changes in the specified frequency due to variations in ambient temperature. The frequency output of a crystal oven depends on the quartz crystal’s temperature. Inside an OCXO is a proportionally controlled heater circuit in addition to the oscillator circuit that uses a heating 4.5. Timing reference design and challenges 84
Figure 4.8: The BVEX timing reference chain which consists of the ultra-precise oven-controlled 10 MHz crystal oscillator (OCXO), distribution amplifiers, and a time interval counter(TICC) to compare the 10 MHz reference with a PPS from the GPS.
element to heat the crystal and a thermistor as a temperature sensor. Overall, these elements are part of a bridge network that acts as a closed-loop feedback network to hold the temperature of the crystal oscillator at a constant temperature, thereby giving us a stable frequency output. For BVEX, we are planning on flying the RK409AVNS Ultra-stable OCXO, which is Rakon’s new range of crystal oscillators specifically designed for near space, low 4.5. Timing reference design and challenges 85
Figure 4.9: a) Frequency stability as a function of temperature for a quartz-based
crystal oscillator (Walls and Vig, 1995) b) inside closed loop of an oven-controlled crystal oscillator c) Rakon RK409AVNS OCXO.
earth orbit environments. These oscillators can withstand extreme shock during rocket launches and the radiation environment of space. A balloon environment, although not as close as a near-space environment, serves as a great proxy to test technology for space. The RK409 has an ADEV of 6 × 10−13 thereby giving us a phase error of 0.09 radians at 23 GHz with 1 sec of integration time according to Eq 4.11 and 4.12. The timing chain that we designed for BVEX, as shown in Figure 4.8, consists of the OCXO as its primary source of reference signal. This 10 MHz signal is attenuated by a 6dB attenuator and fed into a distribution amplifier which splits the signal 3 ways leading into the readout electronics, local oscillator, and a time to digital 4.5. Timing reference design and challenges 86
Components Vendor Part number
OCXO RAKON RK409AVNS 6 dB attenuator Mini-Circuits VAT-6A+ Distribution Amplifier Quantic-Wenzel LNDA-3 Programmable Integrated Circuit (PIC) Microchip Technology PIC16F8DIL Time interval counter TAPR TAPR-TICC Global Positioning System (GPS) receiver Novatel OEM729 RFSoC 4X2 AMD-AMD-Xilinx RFSOC4X2-ZU48DR
Table 4.1: Components used in the timing reference chain along with vendor and part number.
comparator. The local oscillator uses this 10 MHz as its reference to generate the 19 GHz signal for mixing down our RF frequencies. The readout electronics, which is the RFSoC-based FPGA, use the 10 MHz reference from the OCXO for timestamping the raw signals. Finally, the 10 MHz signal is compared against a PPS (pulse per second) from a GPS to record the phase drifts. As discussed earlier, it is extremely crucial to record the phase drifts to account for any phase offset while correlating the signals. Since we cannot compare a 10 MHz signal with a 1 Pulse-Per-Second (PPS) signal, we would be using a PIC microcontroller-based frequency divider (PICDIV) that functions like a series of synchronous decade counters. A simple software program toggles the output pin after every 5 million clock cycles. Thus a 10 MHz input clock signal becomes a 1 Hz square wave, just like a PPS signal. The output accuracy of such a system is equal to the input accuracy, and the jitter is extremely low (< 2 ps) since the PIC is a synchronous device. Finally, to compare and record the phase offsets between the two PPS signals, we plan to use the TAPR Timestamping/Time Interval Counter (TICC) board as 4.5. Timing reference design and challenges 87
referred to in Table 4.1, a time interval counter for comparing PPS signals. The TICC can output timestamp data for each channel, or the time interval between the two channels. It can perform more than 100 measurements per second (both channels) with a resolution of about 60 picoseconds and RMS jitter of fewer than 100 picoseconds that yields a one-second noise Allan Deviation of about 7x10−11 . Finally, this output is fed into the RFSoC and recorded onto a disk along with the raw voltage signals.
4.5.1 OCXO operational environment and Pressure Vessel
Our oven-controlled crystal oscillator (OCXO) has strict operational requirements
to ensure a 10 MHz stable output. It requires:
1. 25◦ C operating temperature;
2. 1 bar of pressure;
3. Sinusoidal vibrations (for each axis, tested over 2 minutes) not exceeding 11mm peak for 5 to 21 Hz and 20mm for 21 to 100 Hz frequency ranges;
4. Random vibrations (for each axis, tested over 2 minutes) not exceeding 0.7 g 2 /Hz for 100-400 Hz frequency ranges.
The phase noise performance of crystal oscillators is generally susceptible to random
vibrations, especially low-frequency sine vibrations. Sinusoidal vibrations produce spectral lines at ±fν from the carrier, where fν is the vibrational frequency:
ΓAf0 L = 20 log , (4.13) 2fν 4.6. BVEX: backend and digitization design 88
where Γ is the g-sensitivity for the RK409AVNS crystal oscillator (0.8ppb/g),
A is the vibration level (0.001 g), and f0 is the frequency of the crystal oscillator output signal (10 MHz). These low-frequency sine vibrations will increase our phase noise as shown in Eq. 4.10, limiting the time we can integrate. During launch, the balloon payload gondola is lifted off by an auxiliary balloon and is released after the balloon takes off. To end the flight, the “cutter”(cutting device) is activated, causing the balloon to separate from the rest of the flight train and rip open. There is a parachute system that brings down everything below the cutting device. In both of these events, large accelerations are expected. Therefore, we need a system to provide a stable operational environment to the OCXO. To address all of the operational requirements of the OCXO and the rest of the sensitive RF components like the mixers and the local oscillator concerning a balloon environment, I plan on designing a temperature-controlled pressure vessel that also provides vibrational isolation by absorbing some of the post-launch and the termination shocks. Our strategy is to have the OCXO powered off during the launch and termination phase. We will obtain more information on the vibrational environment of the gondola by flying a suite of test sensors as a piggyback experiment in August 2023.
4.6 BVEX: backend and digitization design
The raw voltage RF signals from the front end need to be digitized first to store them for correlation later on. Our mixing stage as mentioned above mixes down the RF signals (21-23 GHz) to intermediate frequencies (2-4 GHz). To digitize a signal 4.6. BVEX: backend and digitization design 89
with a bandwidth of 2 GHz (2 to 4 GHz) we need to sample it at a minimum of 4
GHz according to the Nyquist theorem,
fs B< , (4.14) 2
where B is the bandwidth of the signal and fs is the sampling frequency.
Traditionally power-hungry bulky systems consisting of separate ADCs and FP- GAs were used to digitize the radio signals. With the advent of the latest RFSoC (Radio frequency system on chip) by AMD-Xilinx, GHz-range ADCs, and DACs are now integrated into one device along with the FPGA hardware. This direct-RF sub- system implementation allows for fewer components and increased system flexibility. The external RF front-end signal chain as a result can be streamlined, reducing the net power consumption and system complexity. For BVEX we are using the RFSoC 4x2 shown in Figure 4.10 based on AMD- Xilinx’s latest RFSoC Gen3 XCZU48DR FPGA (Table 4.1) with 4 14-bit ADCs capable of sampling at 5 GSPS and 2 DACs at 9.85 GSPS. The ethernet ports support up to 100 Gbps ethernet for high-speed data transfer. The entire BVEX IF stage is summarized in Figure 4.11. The primary goal for the RFSoC 4x2 is:
1. Digitize the data at 5 GSPS;
2. Re-quantize and timestamp the data to 2 bits according to the typical VLBI standard;
3. Calculate the amplitude spectrum of the source;
4.6. BVEX: backend and digitization design 90
Figure 4.10: The latest generation AMD-Xilinx RFSoC 4x2 board with 4 5 GSPS ADCs and 2 9.8 GSPS DACs.
4. Packetize and store all of the data via the high-speed 100 Gbps ethernet port onto a recording unit.
The RFSoC 4x2 will record the output from the time interval counter (TICC) as described in our timing reference chain along with the 10 MHz reference signal and the raw IF voltages. The 10 MHz reference signal from the OCXO will be used to timestamp all our data. Before going into the details of the FPGA designs, we would like to go over the basics of ADCs and all their related terminologies specific to radio astronomy.
1. Quantization of ADC is the process where the Least Significant Bit (LSB) is determined if the input analog voltage lies in the lowest sub-range of the input 4.6. BVEX: backend and digitization design 91
Figure 4.11: Block diagram showing the IF stage design on the RFSoC FPGA and the preliminary data storage plans.
voltage range. For example, let us consider an ADC with Vref = 2V and the resolution is 3-bit. The 2 V reference signal is then divided into 23 or 8 sub- divisions, so the LSB voltage is within 250 mV. Now an input voltage of 0 V, as well as 250 mV, is assigned to the same output digital code 000. The input analog voltage range from 251 mV to 500 mV will be assigned to the same digital code 001 and so on. Due to the digital nature of the ADC, continuous output values are not possible. The perfect ADC performs the quantization process during conversion. This results in a staircase transfer function where the step represents one LSB. The RFSoC 4x2 board with 14-bit (214 = 16,384) would have an LSB voltage of 61 µV . In our case since we will be observing from the stratosphere the antenna temperature will be around 6 K which gives 4.6. BVEX: backend and digitization design 92
us a power of 0.1656 pW assuming a 2 GHz bandwidth,
PT A = kTA ∆ν, (4.15)
where k is Boltzmann’s constant (1.38 × 10−23 J/K), TA is the antenna tem-
perature and ∆ν is the bandwidth. For a 50 Ohm system, this corresponds to an RMS voltage of 2.9 µV for which we would need a gain of at least 38 dB. Our receiver gain for all of the components combined is 54 dB which leaves us with plenty of room for sampling our data efficiently. Usually, since VLBI deals with a lot of data, a smaller number of bits is preferable for reducing the data size and increasing the processing speed. Modern VLBI uses 1-bit,2-bit, or 4-bit quantization for most campaigns (Whitney et al., 2010). Since the ADCs in the RFSoC samples at 14 bits we would need to re-quantize that down to 2 bits before feeding it in the correlator. 2-bit is the default input for the DiFX correlator (Deller et al., 2011) which will be used to correlate our signals post-flight.
2. The noise floor of the ADC is a parameter that determines the lowest signal strength of an input signal that the receiver can recover successfully. If the receiver noise floor within the demodulated bandwidth is larger than the re- ceived signal itself, the input signal cannot be demodulated. The noise floor of an ADC is the difference between the ADC Noise Spectral Density (NSD) in dBm/Hz to the thermal noise density (kT in dBm/Hz) normalized to a 1-Hz 4.6. BVEX: backend and digitization design 93
bandwidth as shown,
(Vpp )2
Pf ull−scale [dBm] = 10 log 1000 , (4.16) 8ZIN
where Vpp is the peak-to-peak voltage. The noise floor for the ADCs in the RFSoC 4x2 was checked with zero input power as shown in Figure 4.12. The input to all the ADCs was terminated by a 50 Ω terminator. The ADCs were configured to sample at 4 GSPS and the ADC samples were captured in chunks of 4096 samples and the corresponding FFTs were calculated using the numpy.fft python library.
Figure 4.12: The noise floor (flat-response) of the 2 ADCs of the RFSoC 4x2, with a tone injected at 100 MHz. The other spurs in the spectrum are due to Radio Frequency Interference (RFI) in the lab. 4.6. BVEX: backend and digitization design 94
The four RF ADC tiles numbered 224 to 227 consisting of two ADCs, ADC01 and ADC23 both sampling at a maximum of 5 GSPS have an interleaved structure. Therefore, each ADC slice in this test is sampled at one-eighth the overall ADC sample rate, i.e. 512 MHz which shows up in the Fast Fourier Transform (FFT) as expected.
3. Spurious-free dynamic range (SFDR) is used to calculate the dynamic perfor-
mance of an ADC, and it is defined as the ratio of the power in the fundamental tone to the power of the worst spur (harmonic distortions/ spikes) in the fre- quency range under investigation.
power in the f undamental tone
SF DR = (4.17) power in the worst spur in the f requency range
SFDR is a crucial specification for radio astronomy receivers as it determines
the smallest signal that can be observed over the spectrum of interest. The RFSoC ADCs have SFDR around 83 to 90 dB across their effective frequency range (∼4 GHz) which is adequate for most radio astronomy applications. In order to calculate the SFDR of the ADCs we performed a simple loop-back frequency sweep test. A loopback test involves sending a signal from a device’s output, routing it back to the input, and analyzing the returned signal to assess the device’s performance and functionality. The DACs were connected into the ADCs and a range of frequencies as listed in Table 4.2 were injected. The received frequency was measured and the SFDR was calculated.
4. Dynamic Range is the range between the noise floor and its specified maximum 4.6. BVEX: backend and digitization design 95
Table 4.2: Result from the sweep test on the ADC where we send in a comb of tones through the DAC and measure the SFDR of the individual ADCs as it receives them.
output level. An ADC’s dynamic range is the range of signal amplitudes that the ADC can resolve. RFI is the main driver of sampling bit-width for modern- day radio environments. RFI is usually orders of magnitude stronger than the actual astronomical signal of interest; therefore, it requires a large dynamic range of the digitized waveform. The theoretical maximum dynamic range of an ADC in dB is given by,
DR = 20 log10 (2nbits ) ≈ 6.02nbits (4.18)
ADCs are practically imperfect analog devices, and their ENOBs (effective- number-of-bits) are lower than the number produced by the ADC.
5. Signal-to-noise and distortions (SINAD) and ENOBs: are a comprehensive
indicator of the overall dynamic range of an ADC, as it takes both noise and distortions into account.
RM S of f undamental tone SIN AD = . mean value of root sum square of all other components (4.19) 4.6. BVEX: backend and digitization design 96
Using SIN AD the ENOB (Kester, 2009) can be calculated from,
SIN AD − 1.76dB EN OB = . (4.20) 6.02
The ENOB of an ADC decreases with an increase in frequency as at higher
frequencies, the quantization error becomes more significant as the input sig- nal changes rapidly, reducing effective resolution. To measure the ENOB, we measure a perfect full-scale sine wave into our ADC, subtract the sine wave, and measure the RMS value of the resulting noise. Then we determine the bit width of an ideal ADC with the same RMS noise, which results in the ENOB. The ENOB, shown on the left axis of Figure 4.13, decreases from 11.87 to 10.49
Figure 4.13: SINAD and ENOB of ADC, with tones generated by signal generator from 300 to 1500 MHz (Liu et al., 2021). 4.7. CASPER-based RFSoC spectrometer 97
bits as the input frequency increases. However, this shouldn’t be an issue for our purposes as we plan to re-quantize the data down to 2 bits.
4.7 CASPER-based RFSoC spectrometer
I designed an FPGA-based spectrometer for BVEX using the Collaboration for As- tronomy Signal Processing and Electronics Research (CASPER) toolflow (Werthimer, 2011). CASPER is an international collaboration of scientists and engineers produc- ing resources to aid astronomy research. The primary goal of CASPER is to simplify the process of radio astronomy instrumentation by developing open-source hardware and software. The CASPER library is designed to be easily integrated with FPGAs, GPUs, and CPU processors. Their well-maintained repo and documentation make it extremely easy to deploy designs on FPGAs. The CASPER toolflow as shown in Figure 4.14 is also called the mlib devel repos- itory containing FPGA-based DSP libraries and programming tools maintained by the collaboration. This toolflow can be accessed via MATLAB’s GUI-based design- ing tool SIMULINK, which lets us generate and implement signal processing designs onto the FPGA by converting them into bitstreams. A bitstream is a binary se- quence comprising a sequence of bits used to program the FPGA. CASPER also provides a Python-based software library called the casperfgpa for interacting with the running designs. For BVEX, we are running the latest version of the CASPER toolflow, Ubuntu 20.04, and MATLAB 2021.a, AMD-Xilinx Vivado 2021.1, casperfpga 0.4.4, and the latest RFSoC 4X2 image. The board boots up from the RFSoC 4x2 CASPER 4.7. CASPER-based RFSoC spectrometer 98
image flashed on a micro-SD card. The casperfpga and its necessary supporting Python libraries have been installed onto a separate virtual environment for easy version control. The SIMULINK environment and all the libraries can be launched from the terminal by activating the virtual environment first and executing the ./mlib devel bash script. Once a design is completed, it can be compiled by passing the jasper command in the MATLAB terminal. The jasper command first uses AMD-Xilinx’s System Generator to compile AMD-Xilinx blocks in the design into a circuit via bitstreams. Next, it runs Vivado’s synthesis, implementation, and route tools to create the physical hardware for the FPGA. After a successful compilation, we get .fpg and .dtbo files as output. These files are then uploaded onto the FPGA and can be later accessed by the casperfpga python library. At its core, a spectrometer is an analysis filter bank that takes in a time domain input signal and transforms it into a frequency domain representation to see what frequency components are in our signal. Since this is an FPGA-based digital system, spectrometers are typically designed utilizing the FFT algorithm. However, we typ- ically use a Polyphase Filterbank (PFB) approach to improve compute speeds and better spectral bin performance. The first digital filter spectrometer used for astronomy was developed by Weinreb (1963). It was a 1-bit system used to observe an 18 cm OH absorption line in the spectrum of Cassiopeia A. This was the first evidence of OH in the interstellar medium. 4.7. CASPER-based RFSoC spectrometer 99
Figure 4.14: Block diagram of the CASPER Toolflow setup.
4.7.1 Polyphase FIR filters and filterbanks
An Finite Impulse Response (FIR) filter shown in Figure 4.15 is a windowed moving average of an input sequence x(n). An FIR filter computes the sum,
K−1 X y(n) = h(k)x(n − k), (4.21) k=0
where y(n) is the output sequence, and h(k) is a set of K weighting coefficients. The upper summation bound, K, is called the number of taps. Polyphase decomposition is a common technique used in Digital Signal Processor (DSP) to decompose an input 4.7. CASPER-based RFSoC spectrometer 100
Figure 4.15: N-tap FIR filter.
sequence x(n) into a set of P sub-sequences, xp (n′ ),
xp (n′ ) = (↓ P )(z −p )x(n). (4.22)
A simple even and odd decomposition of signal x(n) is achieved when P = 2,
xo (n′ ) = {x(0), x(2), x(4), ...}, (4.23)
x1 (n′ ) = {x(1), x(3), x(5), ...}. (4.24)
A signal may be more generally decomposed into P different “phases.”
A PFB is a computationally efficient filter bank implementation constructed from an FFT preceded by a polyphase FIR filter at the front. Polyphase filter bank- based spectrometers offer lower spectral leakage than other architectures like the Autocorrelation spectrometer (ACS) or the Fourier transform spectrometer (FTS). The polyphase filter bank takes advantage of the fact that a lowpass filter with coefficients h(k) can be converted into a quadrature bandpass filter with a central frequency ν by simply multiplying the coefficients by ei2πν . For a decimating lowpass 4.7. CASPER-based RFSoC spectrometer 101
filter, the output from each branch will be,
M X −1 yp (n′ ) = hp (m)xp (n′ − m) (4.25) m=0
where hp (m) are coefficients from the lowpass filter.
To get the output y(n’), we would typically add across the P branches (i.e., across the sub-filters yp ). If we feed the branches (sub−filters) into a DFT with P inputs rather than just adding up the P branches, we then obtain,
P X −1 ′ Y (k, n ) = yp (n′ )e(−2πikp)/k , p=0 P −1 M −1 (4.26) X X = [hp (m)e(−2πikp)/k ]xp (n′ − m). p=0 m=0
This is equivalent to a set of decimating polyphase filters, where the center frequency of each filter is shifted by an amount p/P ; this is a PFB. The output is squared and time-averaged from here to obtain the signal’s power spectral density (PSD). Figure 4.16 shows a simplified diagram of a PFB.
Figure 4.16: Simplified block diagram of a polyphase filter bank.
4.8. Spectrometer design 102
4.8 Spectrometer design
When designing a spectrometer, there are three key parameters we would need to keep in mind: bandwidth, frequency resolution, and time resolution. Bandwidth is defined as the width of the frequency spectrum in Hz. Bandwidth depends on the sampling rate; for complex sampled data, this is equivalent to:
1 ∆ν = sampling rate = (4.27) sampling period
In contrast, the rate is half for real or Nyquist sampled data, as two samples are required to reconstruct a given waveform.
For BVEX since we are Nyquist sampling our data at 3932.16 Msps our effective bandwidth for the spectrometer is 3932.16/2 = 1966.08 MHz. Frequency Resolution is defined as the width of each frequency bin. It is given by ∆ν ∆f = (4.29) no. of channels
This parameter is a measure of how precisely one can measure a frequency. The frequency resolution for our spectrometer with an effective bandwidth of 1966.08 MHz and 2048 frequency channels would be 1966.08/2048 = 0.96 MHz. Time Resolution is defined as the spectral dump rate of the instrument. Usually, we accumulate multiple spectra to average out noise; i.e., the more accumulation 4.8. Spectrometer design 103
you do, the lower the time resolution. For looking at short timescale events, such as pulsar bursts, a higher time resolution is necessary; conversely, if we want to look into a weak HI signal, a long-time accumulation is required. This spectrometer design utilizes the CASPER Yellow (XPS), Green (DSP), and AMD-Xilinx blocks. We will first go over an overview of the signal flow in the design before going into the details of each block.
1. Firstly, the AMD-Xilinx token is paced to allow the System Generator to be
called to compile the overall design.
2. Next, the platform yellow block, i.e., the hardware platform, is set to the target RFSoC platform, which in our case is the RFSoC 4x2. The User IP clock rate (245.76 MHz) is specified in the RFSoC datasheet.
3. Input signals are first digitized by the RFDC block. The RFDC block allows access to the ADCs on the RFSoC, which is configured to sample in real → I/Q mode at a sample rate of 3932.16 Msps, with a decimation factor of 2, and the Numerically Controlled Oscillator (NCO) set to 983.04 MHz. This results in an effective 0 to 1966.08 MHz of bandwidth.
4. The output samples from the Radio Frequency Digital Downconverter (RFDC) block are presented in parallel after each clock cycle and expressed as fixed- point integers.
5. The parallel time samples are then passed through the polyphase FIR filter 4.8. Spectrometer design 104
and FFT block, collectively called the PFB. When the output data from the RFDC are complex, the FFT size is selected to be 211 = 2048.
6. The AMD-Xilinx delay blocks are placed all over the design to help meet the timing standards of the FPGA. This consumes more resources but eases signal timing-induced restrictions.
7. The frequency bin outputs of the PFB are passed through the “power” blocks, which convert from complex-valued outputs to real-valued outputs.
8. The bin power enters the vector accumulators, which are Block Random Acess Memory (BRAM) on the FPGA fabric, which is used to accumulate the signals to a specific length.
9. Finally, the accumulated signal is fed into another BRAM-based yellow block which stores them before dumping them out via the ethernet onto the host machine.
4.8.1 The RFDC Block
The ADCs are used to digitize the incoming signals, which is crucial to the spectrom- eter. For the RFSoC, the ADC is represented by the RFDC (RF Data Converter) yellow block. Every clock cycle, the inputs are sampled and digitized to a 14-bit, 2’s complement binary number representation. These samples are packed, and the Most Significant Bit (MSB) is aligned to a 16-bit word representation. The samples are then presented in a parallel sequence on the output interface. The 16-bit word 4.8. Spectrometer design 105
represents numbers from -32768 to 32767, including 0. This gives us plenty of range to represent our signals for BVEX. The RFSoC 4X2 is a dual-tile platform meaning the 4 ADCs are present and embedded across two tiles on the board. For dual-tile platforms in I/Q digital output modes, the in-phase and quadrature data are produced from different ports. The first digit of the signal name in this mode corresponds to the tile index. The second digit, on the other hand, is 0 for inphase and 1 for quadrature data from the ADC. For example, with four samples per clock, this is four complex samples with two complex components coming from different ports, m00 axis tdata for in-phase data on tile 0 ADC 0 ordered {I3, I2, I1, I0} and m01 axis tdata for tile 0 ADC 1 with quadrature data ordered {Q3, Q2, Q1, Q0}. Each tile tab on the RFSoC4x2 is configured with the same clock settings to ensure the signal syncs throughout the tiles. The output data are now I/Q with a decimation factor of 2 and 8 per sample clock. The decimation factor is simply the input-to-output rate ratio. Decimation filters help us reduce the computational resources required to process the signal. To digitally downconvert the IF signals, the mixer type within the RFDC is set to “Fine” with the Real → I/Q mode enabled, and the NCO set to 983.04 MHz. We finally get digitally down−converted signals which are first shifted down by the fine mixer and then passed through the decimation filters. 4.8. Spectrometer design 106
4.8.2 Munge Blocks
The Munge blocks are used to reorder the output data from the RFDC. We do this because the RFDC outputs the time samples with the newest time samples in the MSB(most significant bit), and the CASPER blocks want the oldest sample in the MSB. The Munge block reverses the word order because the RFDC is configured to output eight samples.
4.8.3 Polyphase FIR
The Polyphase filter bank consists of two main blocks. The first is a FIR block, and the second is the FFT block. For complex data types, we used the “pfb fir” block. A polyphase FIR divides the input signal into parallel “taps” and then applies finite impulse response filters (FIR). FIR filters are used to implement almost every sort of frequency response digitally. The output of this block is still time domain signals. The signals are mapped out in the frequency domain when combined with the FFT block.
4.8.4 CASPER FFT block
The green CASPER FFT block is the heart of the spectrometer, which runs on the FFT butterfly biplex algorithm for improved efficiency and reduced computational complexity. A butterfly operation is a fundamental computational step, involving the combination of two complex numbers through specific arithmetic operations to com- pute discrete Fourier transform values. In a traditional butterfly structure, two com- plex numbers are multiplied by twiddle factors and added or subtracted to produce 4.8. Spectrometer design 107
the output values. A twiddle factor, in FFT algorithms, is any of the trigonometric constant coefficients that are multiplied by the data in the course of the algorithm (Gentleman and Sande, 1966). In the biplex design, the butterfly structure is split into two parallel paths, one for addition and one for subtraction. This separation enables more efficient hardware implementation and reduces the required complex multiplications. The quantity of ports determines how many simultaneous bin out- puts are generated. The left and right halves of the FFT spectrum are mirror images for a real input signal. The block does not generate the imaginary (negative channel) indices. Hence, for a 4096-point FFT, we get about 2048 channels. As a result, there are only half as many parallel outputs. These parallel FFT outputs will produce sequential channels on every clock cycle.
4.8.5 Power Block
The power block is used in the design as it computes the power of a complex number. At its core, the block takes real and imaginary numbers from the FFT output, squares, and sums it up.
4.8.6 Sync Gen Block
The CASPER DSP blocks require a synchronized pulse to function. In the case of a spectrometer, sync pulses define the spectrum frames. For BVEX, we would use a 1 PPS signal from our OCXO to provide our design with a pulse. However, since we do not have our OCXO yet, we used the internal voltage-controlled oscillator in the RFSoC using the sync gen block to build this design. If a sync comes out of 4.9. Results of the spectrometer 108
alignment, it automatically resets itself.
4.8.7 Vector Accumulator Block
This block is used in our design to accumulate the output from the FFT before sending it out to the BRAM of the FPGA memory to prevent it from overflowing. This block accumulates approximately 28 bits each second. Each accumulator will capture a fraction of the total output bandwidth in parallel, and the frequency bins will have to be interleaved to reconstruct the output spectrum while reading back.
4.8.8 Shared BRAM blocks
This is the final block in our spectrometer design from which we read all our data via a Python script. The BRAM is a dual-port RAM module built into an FPGA fabric to provide temporary on-chip storage for a relatively large set of data (up to 36 Kb), the raw spectra. While setting the parameters of this block, we had to ensure that it matched the data width of the vector accumulator block.
4.9 Results of the spectrometer
To read out the spectrum data from the RFSoC 4x2, I wrote a sample Python script that fetches the data from the BRAM within the FPGA fabric, analyses them, and finally plots the spectra. We defined a function that uses the casperfpga library to access the BRAM and fetch data from all eight channels in chunks. Each chunk 4.9. Results of the spectrometer 109
length has been defined by,
no. of F F T points 2048
lchunk = = = 256. no. of channel 8
The function then reconstructs the output by interleaving all parallel data from the eight channels. It finally returns a reconstructed raw data stream that is plotted as shown in Figure 4.18. The main goal of this spectrometer for BVEX would be to observe the spectra of bright water maser lines at 22 GHz. Water maser emissions from Orion KL observed between 1997 to 1998 have channel widths from 7.9 km/s to 10.8 km/s (Omodaka et al., 1999). These channel widths can be converted in Hz using the Doppler shift,
∆v ∆ν = × ν0 . (4.30) c
Where ∆ν is the channel width in Hz, ∆v is the channel width in km/s, c is the speed of light in km/s and ν0 is the rest frequency of observation in Hz. Applying Eq. 4.30 we get line widths between 0.58 to 0.74 MHz for the water maser emission line in Orion KL. Since we are sampling our data at 1966 MHz and using 2048 FFT points, we get a frequency resolution of df = f s/nF F T = 0.96 MHz. This should accommodate our signal leaving plenty of room. The Python code can be executed by passing the file name followed by the IP address of the board. It plots out the spectrum in real-time. For testing out our design we used a 580 MHz Voltage-Controlled Crystal Oscillator (VCXO) and con- nected it to the ADC of RFSoC 4x2 running the spectrometer. Figure 4.18 shows us 4.9. Results of the spectrometer 110
the result of our test as it successfully plots the spectra of the 580 MHz source.
Figure 4.18: RFSoC 4X2 based spectrometer reading out the 580 MHz sine wave spectrum generated from a VCXO.
Overall, in this chapter, we introduced the Balloon-borne VLBI Experiment
(BVEX) a pathfinder experiment that has the goal of demonstrating the feasibil- ity of a high-frequency stratospheric interferometry station. We briefly discussed the key objectives of BVEX and its front-end design. Finally, we designed the timing reference responsible for providing a stable 10 MHz output to the local oscillator and timestamp our signals, and also designed an RFSoC-based spectrometer which will be a key part of the backend digitization system. 4.9. Results of the spectrometer
Figure 4.17: Design of the CASPER-based wideband spectrometer for BVEX.
111 Chapter 5
Conclusion
5.1 Summary
The motivation for this work was to develop next-generation readout techniques and analysis to improve the resolution of microwave and sub-mm astronomical instru- ments in order to better understand the early stages of star formation. We performed a statistical analysis of the polarization properties of cores, in order to understand whether we can use thermal dust polarization to study the mag- netic fields of cores. We started this study by establishing how dust could be an excellent proxy for probing the plane-of-sky magnetic field orientation within these clouds. Analyzing the polarized thermal emission from various molecular clouds, we found that protostellar systems with a compact 70 µm source showed higher polar- ized intensity compared to the ones that did not have a 70 µm source indicating that perhaps the dust grains are more efficiently aligned in cores with protostellar activities. However, to truly establish this and map the local magnetic field morphol- ogy at the core scale, we need higher resolution maps of many other star-forming to
112 5.1. Summary 113
have a significant statistic. Hence to put some effort towards this, we looked at the upcoming Prime-Cam instrument on the 6-m sub-mm CCAT observatory. Next we carried out a thorough atmospheric modeling of the CCAT Observatory site in order to characterize the varying atmospheric power loading on the detectors used for the Prime-Cam instrument. We found that for the higher 850 GHz band, we would need an active tone tracking system even for low-elevation scans. Also, scans at low elevations (30◦ -40◦ ) with the 350 GHz array would also require a tone- tracking system. Overall, the Chilean winter months from May to September would be ideal for observing across all the frequency bands as the atmosphere is significantly drier and the temperatures are constant. CCAT, with its high spatial resolution, sensitivity, and wide area surveys, will be a crucial instrument for studying the polarization of young cores and understanding the role of large-scale magnetic fields in shaping their formation. To take this even further, we investigated a mobile stratospheric sub-mm VLBI station that can observe simultaneously with other ground-based sub-mm telescopes to produce images with micro-arcsec resolution. BVEX, a pathfinder experiment, will demonstrate this at K-band (21-23 GHz) pairing with a large ground-based dish like the Green Bank Telescope. In this work, we made an initial design towards RF- SoC FPGA-based backend readout and timing reference chain. Since VLBI requires precise time stamping of data, we are flying an ultra-stable oven-controlled crystal oscillator (OCXO) with an Allan deviation (ADEV) requirement of ∼ 7×10−12 at 23 GHz for 1-second timescales. Using the latest generation AMD-Xilinx RFSoC 4x2 FPGA, we designed a fully functioning PFB based 2048 channel spectrometer with 5.2. Next Steps 114
the CASPER toolflow. This spectrometer would take the spectra of bright water masers at 22 GHz and start recording the raw digitized and re-quantized data to the disks. The 5 GSPS sampling speed of the RFSoC 4x2 helped us eliminate the need for a dual-stage down conversion system, simplifying the entire telescope receiver design. Overall, BVEX will be a critical balloon-borne instrument to demonstrate the technologies required for high-frequency stratospheric VLBI stations, which can potentially be a part of ngEHT or the global VLBI network.
5.2 Next Steps
The immediate future work towards this study can be outlined as follows:
1. Investigating data from higher resolution maps of star-forming molecular clouds
like the upcoming fields in filament survey using the TolTEC camera on the LMT. The 5” resolution will allow us to probe the individual cores within the molecular clouds and help us map the magnetic field orientation within them.
2. Developing an active RFSoC-based tone tracking firmware to ensure that the
MKIDs are always on resonance and efficiently read them out. This will in- clude efforts toward designing an optimized PID-based algorithm and extensive testing on the FPGA.
3. Further work is required for analyzing the data from a piggyback sensor package experiment towards characterizing the vibrational environment of the BVEX gondola and designing an active vibration-damping system for minimizing the phase noise of the OCXO. 5.2. Next Steps 115
4. The complete digitization, re-quantization, and transportation of the UDP
(User Datagram Protocol) data packets to a NAS (Network-attached storage) based storage system need to be designed on the RFSoC.
5. The next step of BVEX will be a high-frequency mobile VLBI station which will complement ALMA and ngVLA and help us resolve more black hole shad- ows and distant star-forming cores. The ultimate version of BVEX will be a terahertz space-based VLBI network with cooled optics sensitive to far IR wavelengths. Such an instrument will be capable of resolving >100 black hole shadows and extragalactic star-forming cores. Bibliography
Yashar Akrami, Kristian Joten Andersen, Mark Ashdown, Carlo Baccigalupi, M Bal- lardini, Anthony J Banday, RB Barreiro, Nicola Bartolo, S Basak, K Benabed, et al. Planck intermediate results-lvii. joint planck lfi and hfi data processing. Astronomy & Astrophysics, 643:A42, 2020.
FO Alves, Pau Frau, Josep Miquel Girart, Gabriel AP Franco, Fábio P Santos, and Helmut Wiesemeyer. On the radiation driven alignment of dust grains: Detection of the polarization hole in a starless core. Astronomy & Astrophysics, 569:L1, 2014.
BG Andersson, A Lazarian, and John E Vaillancourt. Interstellar dust grain align-
ment. Annual Review of Astronomy and Astrophysics, 53:501–539, 2015.
Ph André, Alexander Men’shchikov, Sylvain Bontemps, V Könyves, F Motte,
N Schneider, P Didelon, V Minier, P Saraceno, Derek Ward-Thompson, et al. From filamentary clouds to prestellar cores to the stellar imf: Initial highlights from the herschel gould belt survey. Astronomy & Astrophysics, 518:L102, 2010.
Philippe André, Shantanu Basu, and Shu-ichiro Inutsuka. The formation and evolu- tion of prestellar cores. Structure Formation in Astrophysics, page 254, 2009.
116 Bibliography 117
Manuel Aravena, Jason E Austermann, Kaustuv Basu, Nicholas Battaglia, Ben-
jamin Beringue, Frank Bertoldi, Frank Bigiel, J Richard Bond, Patrick C Breysse, Colton Broughton, et al. Ccat-prime collaboration: Science goals and forecasts with prime-cam on the fred young submillimeter telescope. The Astrophysical Journal Supplement Series, 264(1):7, 2022.
David Ball, Chi-kwan Chan, Pierre Christian, Buell T Jannuzi, Junhan Kim, Daniel P Marrone, Lia Medeiros, Feryal Ozel, Dimitrios Psaltis, Mel Rose, et al. First m87 event horizon telescope results. i. the shadow of the supermassive black hole. 2019.
I Barrueto, UU Graf, CE Honingh, K Jacobs, M Justen, H Krüger, M Schultz,
K Vynokurova, L Weikert, S Wulff, et al. Ccat-prime heterodyne instrument (chai) advances. In 2022 47th International Conference on Infrared, Millimeter and Terahertz Waves (IRMMW-THz), pages 1–2. IEEE, 2022.
Rebecca A Bernstein, Wendy L Freedman, and Barry F Madore. The first detections of the extragalactic background light at 3000, 5500, and 8000 åi results. The Astrophysical Journal, 571(1):56, 2002.
Paola Caselli, Priscilla J Benson, Philip C Myers, and Mario Tafalla. Dense cores in dark clouds. xiv. n2h+ (1-0) maps of dense cloud cores. The Astrophysical Journal, 572(1):238, 2002.
A Cattaneo, SM Faber, J Binney, A Dekel, J Kormendy, R Mushotzky, A Babul,
PN Best, M Brüggen, AC Fabian, et al. The role of black holes in galaxy formation and evolution. Nature, 460(7252):213–219, 2009. Bibliography 118
Scott C Chapman, Anthony I Huber, Adrian K Sinclair, Jordan D Wheeler, Jason E
Austermann, James Beall, James Burgoyne, Steve K Choi, Abigail Crites, Cody J Duell, et al. Ccat-prime: the 850 ghz camera for prime-cam on fyst. In Millimeter, Submillimeter, and Far-Infrared Detectors and Instrumentation for Astronomy XI, volume 12190, pages 77–91. SPIE, 2022.
Steve K Choi, Jason Austermann, Kaustuv Basu, Nicholas Battaglia, Frank Bertoldi, Dongwoo T Chung, Nicholas F Cothard, Shannon Duff, Cody J Duell, Patricio A Gallardo, et al. Sensitivity of the prime-cam instrument on the ccat-prime tele- scope. Journal of Low Temperature Physics, 199:1089–1097, 2020.
Event Horizon Telescope Collaboration et al. First m87 event horizon tele- scope results. iv. imaging the central supermassive black hole. arXiv preprint arXiv:1906.11241, 2019.
Thomas M Dame, Dap Hartmann, and P Thaddeus. The milky way in molecular clouds: a new complete co survey. The Astrophysical Journal, 547(2):792, 2001.
Peter K Day, Henry G LeDuc, Benjamin A Mazin, Anastasios Vayonakis, and Jonas Zmuidzinas. A broadband superconducting detector suitable for use in large arrays. Nature, 425(6960):817–821, 2003.
AT Deller, WF Brisken, CJ Phillips, John Morgan, W Alef, R Cappallo, E Middel-
berg, J Romney, H Rottmann, SJ Tingay, et al. Difx-2: a more flexible, efficient, robust, and powerful software correlator. Publications of the Astronomical Society of the Pacific, 123(901):275, 2011. Bibliography 119
Yasuo Doi, Tetsuo Hasegawa, Ray S Furuya, Simon Coudé, Charles LH Hull, Doris Arzoumanian, Pierre Bastien, Michael Chun-Yuan Chen, James Di Francesco, Rachel Friesen, et al. The jcmt bistro survey: Magnetic fields associated with a network of filaments in ngc 1333. The Astrophysical Journal, 899(1):28, 2020.
Bruce T Draine. Interstellar dust grains. Annual Review of Astronomy and Astro- physics, 41(1):241–289, 2003.
BT Draine and Aigen Li. Infrared emission from interstellar dust. i. stochastic heating of small grains. The Astrophysical Journal, 551(2):807, 2001.
BT Draine and Joseph C Weingartner. Radiative torques on interstellar grains: I.
Michael M Dunham, Antonio Crapsi, Neal J Evans II, Tyler L Bourke, Tracy L Huard, Philip C Myers, and Jens Kauffmann. Identifying the low-luminosity pop- ulation of embedded protostars in the c2d observations of clouds and cores. The Astrophysical Journal Supplement Series, 179(1):249, 2008.
Felix Thiel. Improving the resolution and sensitivity of microwave astronomical
Laura M Fissel, Peter AR Ade, Francesco E Angilè, Peter Ashton, Steven J Ben- ton, Mark J Devlin, Bradley Dober, Yasuo Fukui, Nicholas Galitzki, Natalie N Gandilo, et al. Balloon-borne submillimeter polarimetry of the vela c molecular cloud: systematic dependence of polarization fraction on column density and local polarization-angle dispersion. The Astrophysical Journal, 824(2):134, 2016. Bibliography 120
Per Friberg, Pierre Bastien, David Berry, Giorgio Savini, Sarah F Graves, and Kate Pattle. Pol-2: a polarimeter for the james-clerk-maxwell telescope. In Millime- ter, Submillimeter, and Far-Infrared Detectors and Instrumentation for Astronomy VIII, volume 9914, page 991403. SPIE, 2016.
AS Fruchter and RN Hook. Drizzle: A method for the linear reconstruction of
undersampled images. Publications of the Astronomical Society of the Pacific, 114 (792):144, 2002.
W Morven Gentleman and Gordon Sande. Fast fourier transforms: for fun and profit. In Proceedings of the November 7-10, 1966, fall joint computer conference, pages 563–578, 1966.
Reinhard Genzel, Frank Eisenhauer, and Stefan Gillessen. The galactic center mas- sive black hole and nuclear star cluster. Reviews of Modern Physics, 82(4):3121, 2010.
Jo R Herrnstein, James M Moran, Lincoln Jared Greenhill, PJ Diamond, M Inoue,
N Nakai, M Miyoshi, C Henkel, and A Riess. A geometric distance to the galaxy ngc4258 from orbital motions in a nuclear gas disk. Nature, 400(6744):539–541, 1999.
Roger H Hildebrand. Magnetic fields and stardust. Royal Astronomical Society,
Jonathan R Hoh, Adrian Sinclair, and Ryan Stephenson. Maximizing dynamic range Bibliography 121
of microwave kinetic inductance detectors through high-speed tone tracking. In
Millimeter, Submillimeter, and Far-Infrared Detectors and Instrumentation for Astronomy X, volume 11453, page 114530B. SPIE, 2020.
WS Holland, D Bintley, EL Chapin, Antonio Chrysostomou, GR Davis, JT Dempsey,
WD Duncan, M Fich, P Friberg, M Halpern, et al. Scuba-2: the 10 000 pixel bolometer camera on the james clerk maxwell telescope. Monthly Notices of the Royal Astronomical Society, 430(4):2513–2533, 2013.
Charles LH Hull, Richard L Plambeck, Woojin Kwon, Geoffrey C Bower, John M
Carpenter, Richard M Crutcher, Jason D Fiege, Erica Franzmann, Nicholas S Hakobian, Carl Heiles, et al. Tadpol: A 1.3 mm survey of dust polarization in star-forming cores and regions. The Astrophysical Journal Supplement Series, 213 (1):13, 2014.
Tsuyoshi Inoue, Ryo Yamazaki, and Shu-ichiro Inutsuka. Turbulence and magnetic field amplification in supernova remnants: Interactions between a strong shock wave and multiphase interstellar medium. The Astrophysical Journal, 695(2):825, 2009.
Shu-ichiro Inutsuka, Tsuyoshi Inoue, Kazunari Iwasaki, and Takashi Hosokawa. The formation and destruction of molecular clouds and galactic star formation-an ori- gin for the cloud mass function and star formation efficiency. Astronomy & Astro- physics, 580:A49, 2015.
Michael Johnson, Kari Haworth, Dominic W Pesce, Daniel Palumbo, Lindy Black- burn, Kazunori Akiyama, Don Boroson, Katherine L Bouman, Joseph R Farah, Bibliography 122
Vincent L Fish, et al. Studying black holes on horizon scales with space-vlbi. Bulletin of the American Astronomical Society, 51(7):235, 2019.
Terry Jay Jones, M Bagley, M Krejny, B-G Andersson, and P Bastien. Grain align- ment in starless cores. The Astronomical Journal, 149(1):31, 2014.
F Kemper, C Jäger, LBFM Waters, Th Henning, FJ Molster, MJ Barlow, T Lim,
and A De Koter. Detection of carbonates in dust shells around evolved stars. Nature, 415(6869):295–297, 2002.
Walt Kester. Understand sinad, enob, snr, thd, thd+ n, and sfdr so you don’t get lost in the noise floor. MT-003 Tutorial, www. analog. com/static/importedfiles/tutorials/MT-003. pdf, 2009.
Helen Kirk, Doug Johnstone, and Mario Tafalla. Dynamics of dense cores in the perseus molecular cloud. The Astrophysical Journal, 668(2):1042, 2007.
A Lazarian and Thiem Hoang. Alignment of dust by radiative torque: Recent de- velopments. arXiv preprint arXiv:0901.0146, 2008.
Alex Lazarian and Thiem Hoang. Radiative torques: analytical model and basic properties. Monthly Notices of the Royal Astronomical Society, 378(3):910–946, 2007.
Aigen Li and BT Draine. Are silicon nanoparticles an interstellar dust component?
The Astrophysical Journal, 564(2):803, 2002. Bibliography 123
Chao Liu, Michael E Jones, and Angela C Taylor. Characterizing the performance of high-speed data converters for rfsoc-based radio astronomy receivers. Monthly Notices of the Royal Astronomical Society, 501(4):5096–5104, 2021.
Nathan P Lourie, Peter AR Ade, Francisco E Angile, Peter C Ashton, Jason E
Austermann, Mark J Devlin, Bradley Dober, Nicholas Galitzki, Jiansong Gao, Sam Gordon, et al. Preflight characterization of the blast-tng receiver and detector arrays. In Millimeter, Submillimeter, and Far-Infrared Detectors and Instrumen- tation for Astronomy IX, volume 10708, pages 52–66. SPIE, 2018.
KA Marsh, Jason Matthew Kirk, Ph André, MJ Griffin, V Könyves, P Palmeirim,
A Men’Shchikov, Derek Ward-Thompson, Milena Benedettini, DW Bresnahan, et al. A census of dense cores in the taurus l1495 cloud from the herschel gould belt survey. Monthly Notices of the Royal Astronomical Society, 459(1):342–356, 2016.
PD Mauskopf. Transition edge sensors and kinetic inductance detectors in astro-
nomical instruments. Publications of the Astronomical Society of the Pacific, 130 (990):082001, 2018.
A Men’Shchikov, Ph André, P Didelon, F Motte, M Hennemann, and N Schneider.
Thomas Nikola, Steve Choi, Cody Duell, Rodrigo G Freundt, Zachary Huber, Yaqiong Li, Kshama Malavalli, Mike Niemack, Kayla Rossi, Gordon Stacey, et al. Bibliography 124
Ccat-prime: the epoch reionization spectrometer for primce-cam on fyst. In Mil-
limeter, Submillimeter, and Far-Infrared Detectors and Instrumentation for As- tronomy XI, volume 12190. SPIE, 2022.
Nishio, Takeshi Miyaji, Nanako Motiduki, Masaki Morimoto, Hideyuki Kobayashi, and Tetsuo Sasao. The enormous outburst of the 7.9 km s-1 water-maser feature in orion kl. Publications of the Astronomical Society of Japan, 51(3):333–336, 1999.
Scott Paine. The am atmospheric model. Zenodo, 2012.
HAL Parsons, DS Berry, MG Rawlings, and SF Graves. The pol-2 data reduction cookbook, 2018.
Kate Pattle, Derek Ward-Thompson, David Berry, Jennifer Hatchell, Huei-Ru Chen, Andy Pon, Patrick M Koch, Woojin Kwon, Jongsoo Kim, Pierre Bastien, et al. The jcmt bistro survey: the magnetic field strength in the orion a filament. The Astrophysical Journal, 846(2):122, 2017.
Stefano Pezzuto, Milena Benedettini, J Di Francesco, P Palmeirim, S Sadavoy, EU-
GENIO Schisano, G Li Causi, Ph André, D Arzoumanian, J-Ph Bernard, et al. Physical properties of the ambient medium and of dense cores in the perseus star- forming region derived from herschel gould belt survey observations. Astronomy & Astrophysics, 645:A55, 2021.
Jaime E Pineda, Doris Arzoumanian, Philippe André, Rachel K Friesen, Annie Zav- agno, Seamus D Clarke, Tsuyoshi Inoue, Che-Yu Chen, Yueh-Ning Lee, Juan D Bibliography 125
Soler, et al. From bubbles and filaments to cores and disks: gas gathering and growth of structure leading to the formation of stellar systems. arXiv preprint arXiv:2205.03935, 2022.
Frédérick Poidevin, Diego Falceta-Gonçalves, Grzegorz Kowal, Elisabete de Gouveia
Dal Pino, and Antonio Mário Magalhães. Magnetic field components analysis of the scupol 850 µm polarization data catalog. The Astrophysical Journal, 777(2): 112, 2013.
Riwaj Pokhrel, Philip C Myers, Michael M Dunham, Ian W Stephens, Sarah I Sa- davoy, Qizhou Zhang, Tyler L Bourke, John J Tobin, Katherine I Lee, Robert A Gutermuth, et al. Hierarchical fragmentation in the perseus molecular cloud: from the cloud scale to protostellar objects. The Astrophysical Journal, 853(1):5, 2018.
SI Sadavoy, J Di Francesco, Ph André, S Pezzuto, J-P Bernard, Sylvain Bontemps,
E Bressert, S Chitsazzadeh, C Fallscheer, M Hennemann, et al. Herschel observa- tions of a potential core-forming clump: Perseus b1-e. Astronomy & Astrophysics, 540:A10, 2012.
SI Sadavoy, J Di Francesco, Ph André, Stefano Pezzuto, J-P Bernard, A Maury,
A Men’Shchikov, F Motte, Q Nguy, N Schneider, et al. Class 0 protostars in the perseus molecular cloud: a correlation between the youngest protostars and the dense gas distribution. The Astrophysical journal letters, 787(2):L18, 2014.
Juan D Soler and Patrick Hennebelle. What are we learning from the relative ori- entation between density structures and the magnetic field in molecular clouds? Astronomy & Astrophysics, 607:A2, 2017. Bibliography 126
PM Solomon, AR Rivolo, J Barrett, and A Yahil. Mass, luminosity, and line width relations of galactic molecular clouds. The Astrophysical Journal, 319:730–741, 1987.
GJ Stacey, F Bertoldi, B Magnelli, J Stutzki, P Schilke, J Glenn, J Bally, R Gio-
vanelli, D Riechers, A Karim, et al. The ccat observatory: Science and facility. EAS Publications Series, 75:419–425, 2015.
Paul Tiede, Michael D Johnson, Dominic W Pesce, Daniel CM Palumbo, Dominic O
Chang, and Peter Galison. Measuring photon rings with the ngeht. Galaxies, 10 (6):111, 2022.
Eve M Vavagiakis, Zeeshan Ahmed, Aamir Ali, Kaustuv Basu, Nicholas Battaglia, Frank Bertoldi, Richard Bond, Ricardo Bustos, Scott C Chapman, Dongwoo Chung, et al. Prime-cam: a first-light instrument for the ccat-prime telescope. In Millimeter, Submillimeter, and Far-Infrared Detectors and Instrumentation for Astronomy IX, volume 10708, pages 187–202. SPIE, 2018.
F.L. Walls and J.R. Vig. Fundamental limits on the frequency stabilities of crys- tal oscillators. IEEE Transactions on Ultrasonics, Ferroelectrics, and Frequency Control, 42(4):576–589, 1995. doi: 10.1109/58.393101.
Sander Weinreb. A digital spectral analysis technique and its application to radio astronomy. PhD thesis, Massachusetts Institute of Technology, 1963.
Dan Werthimer. The casper collaboration for high-performance open source digital Bibliography 127
radio astronomy instrumentation. In 2011 XXXth URSI general assembly and
scientific symposium, pages 1–4. IEEE, 2011.
Alan Whitney, Mark Kettenis, Chris Phillips, and Mamoru Sekido. Vlbi data inter- change format (vdif). In Proceedings of the Sixth General Meeting of the Interna- tional VLBI Service for Geodesy and Astrometry, 2010. Appendix A
Atmospheric simulations for the CCAT-prime site
In this section, we will review the results from am simulations for characterizing the atmosphere above the CCAT-prime site at the summit of the Cerro-Chajnantor plateau 5600 m above sea level.
A.1 Overall site characteristics
Figure A.1 shows us the annual average of atmospheric transmission as a function of
frequency for various H2 O percentiles. A value of 1 represents 100 % transmission of the radio waves or a transparent atmosphere, whereas 0 represents no transmission hence an opaque atmosphere. As we can see, the atmosphere is highly transparent for the 280 and 350 GHz bands; however, the transmission value falls off as we go higher in frequency. The 850 GHz array occupies the last atmospheric window where the transmission drastically reduces on increasing H2 O percentiles. Next, in we look into the monthly variations of atmospheric transmission and brightness temperature for different H2 O profiles to determine which time of the
128 A.1. Overall site characteristics 129
Figure A.1: Annual average of frequency vs. transmission for various H2 O profiles and months. The vertical lines represent the three different frequency bands of CCAT-prime.
year is ideal for observations.
The 5th percentile H2 O profile A.2 is the best-case scenario, usually the dissolved H2 O percentile in the atmosphere will be higher than this. As we can see in terms of atmospheric transmission and brightness temperature, June to November performs better than the other months as these months tend to be drier and colder. These months are suitable, especially for high-frequency 850 GHz observation. The 25th (A.3) and 50th (A.4) percentile H2 O values closely represent the ground A.1. Overall site characteristics 130
Figure A.2: 5th percentile H2 O where the first panel shows atmospheric transmission and the second panel shows brightness temperature (Tb ) as a function of frequency.
condition at the CCAT-prime site. Again, the winter months from May to November perform much better than the summer months. Average atmospheric transmission hovers near 50% for the 25th percentile and touches the 25% mark from March to May for the 50th percentile H2 O profile. The 850 GHz sky brightness temperatures settle at the high 200 K almost throughout the band for the 50th percentile H2 0. The 75th (A.5) and 95th (A.6) percentile H2 O profiles are the worst-case scenario for the CCAT-prime site. The am simulation results show that observations would be difficult even for the 350 GHz band. The atmospheric transmission for the 850 GHz frequency band is close to 0, i.e., the atmosphere is almost opaque. The overall conclusion is higher frequency 850 GHz observations should preferably be carried out from May to November, i.e., during the cold, dry months. Observing conditions across all the frequency bands deteriorate with increasing H2 0 percentile values. A.1. Overall site characteristics 131
Figure A.3: 25th percentile H2 O where the first panel shows atmospheric transmission and the second panel shows brightness temperature (Tb ) as a function of frequency.
A.1.1 Detector-band wise analysis
In this section, we will see how the various atmospheric parameters like opacity, transmission, and detector parameters like power loading and frequency shifts vary as a function of elevation angle for various months. These results will help plan scanning strategies for the Prime-Cam instrument. Figure A.7 shows the variation of atmospheric opacity as a function of elevation angles for different months. The atmosphere is on average almost a magnitude more opaque for the 850 GHz band. Also, the opacity drops with an increase in elevation angle as we observe through less of the atmospheric column above the telescope. Next, Figure A.8 shows us the variation of Rayleigh-Jeans temperature as a func- tion of elevation angles for different months. Rayleigh-Jeans temperature increases A.1. Overall site characteristics 132
Figure A.4: 50th percentile H2 O where the first panel shows atmospheric transmission and the second panel shows brightness temperature (Tb ) as a function of frequency.
drastically with frequency as the 850 GHz Prime-Cam band has a significantly higher temperature compared to the 280 GHz and 350 GHz bands. A similar trend is seen with a decrease in the elevation angle as we observe through less of the atmospheric column above the telescope. As discussed in Chapter 3, the loading power on the detectors in our case is predominantly from the atmosphere and depends on Rayleigh-Jeans temperature and detector bandwidth (3.14). Figure A.9 shows how the loading power on the detectors increases by almost an order of magnitude from the 280 to 850 GHz band considering a 10% optical efficiency. Again, the loading power decreases with an increase in elevation angle. Finally, Figure A.10 shows us the possible changes in elevation before being thrown off the linear regime as a function for different months. This shows us that A.1. Overall site characteristics 133
Figure A.5: 75th percentile H2 O where the first panel shows atmospheric transmission and the second panel shows brightness temperature (Tb ) as a function of frequency.
it will be much harder for the 850 GHz band to maintain its linear regime before it needs to be re-tuned especially for low elevation angles. Therefore, for the 850 GHz array, we would need a tone-tracking algorithm to ensure that the MKIDs are always within their linear operation regime. A.1. Overall site characteristics 134
Figure A.6: 95th percentile H2 O where the first panel shows atmospheric transmission and the second panel shows brightness temperature (Tb ) as a function of frequency. A.1. Overall site characteristics 135
Figure A.7: Skydip analysis of atmospheric opacity as a function of elevation angle
for different months. Each panel shows the annual average of the three different frequency bands of CCAT-prime. A.1. Overall site characteristics 136
Figure A.8: Skydip analysis of Rayleigh-Jeans temperature as a function of elevation
angle for different months. Each panel shows the annual average of the three different frequency bands of CCAT-Prime. A.1. Overall site characteristics 137
Figure A.9: Skydip analysis of detector loading power as a function of elevation angle for different months. Each panel shows the annual average of the three different frequency bands of CCAT-prime. A.1. Overall site characteristics 138
Figure A.10: Skydip analysis of possible changes in elevation before being thrown off the linear regime as a function for different months. Each panel shows the annual average of the three different frequency bands of CCAT-prime.
