Double-lens Scintillometry: The variable scintillation of pulsar B1508+55

We report on observations of PSR B1508+55's scintillation at the Effelsberg 100-m telescope spanning from early 2020 to early 2022. In the autumn of 2020, close to the time the pulsar was predicted to cross echoes in its pulse profile, a sudden transition in the scintillation arcs from peculiar stripe-like features to parabolic arclets was observed. To infer a geometric model of the scattering we measure the effects of the annual velocity curve of Earth, of the relative movement of the line of sight, and of the projection of points on a second scattering screen and develop novel methods to do so. The latter phenomenon was discovered by this study and strongly indicates a two-screen scattering geometry. We derive an analytical two-screen model and demonstrate in a Markov Chain Monte Carlo analysis as well as simulations that it can be successfully applied to explain the observations by interpreting the transition as a change of relative amplitudes of images as well as a shift in the orientation of anisotropy. The collection of methods we demonstrate here is transferable to other pulsars with the potential to strongly improve constraints on scattering models.


INTRODUCTION
Compact radio sources like pulsars, fast radio bursts (FRBs), and a few quasars emit spatially coherent radiation which can be scattered by the ionized interstellar medium (IISM). Interference of the scattered paths of propagation leads to a variation of the intensity of the measured signal that is called scintillation. The characteristic time and frequency scales of this variation depend on the geometrical situation. In particular, they depend on the proper motion of the source perpendicular to the line of sight. The science of using scintillation to determine astrometric properties like orbits of pulsars and the distances of scattering structures within the IISM has become known as scintillometry, while imaging these scattering structures can be understood as a form of in-line holography (Walker et al. 2008).
Interstellar scintillation (ISS) can be observed as intensity variations on many different scales. The variations studied here manifest in intensity variations within the timespans and bandwidths of single observations that tend not to alter the mean flux of the source. This kind of scintillation is caused by the Earth moving through an interference pattern. There are also long-time intensity variations. These are an incoherent effect that is caused by relative movement such that different parts of the scattering material cross the line of sight. Thus, the amplification of the flux changes.
Scintillation resembles an interference pattern in the dynamic spectrum over time and frequency. Its structure is often more clear in the power spectrum of the dynamic spectrum, which is called the seconday spectrum. Secondary spectra are dominated by parabolic ★ E-mail: tsprenger@mpifr-bonn.mpg.de structures refered to as scintillation arcs. These arcs were first described in Stinebring et al. (2001). Walker et al. (2004); Cordes et al. (2006) successfully explained them by invoking thin and highly anisotropic scattering screens whose interstellar optical properties were analysed by Gwinn et al. (1998).
The physical nature of these screens, that are solely observable through the phase delay they impose on passing radiation, still needs to be determined on the small angular scales that cause scattering. Supernova remnants (SNR) are very likely physical associations, as shown for the Crab pulsar (Cordes et al. 2004), J0538+2817 (Yao et al. 2021), and likely for J0540−6919 (Geyer et al. 2021). Other likely associations are hot stars (e.g. Walker et al. 2017) and extended HII regions (e.g. Mall et al. 2022). Scattering is caused by variable densities of dispersive material that delays the arrival of pulses. Thus, scattering screens also cause variations in dispersion which is quantified by the dispersion measure (DM) of a pulsar. This correlation has been observed (e.g. McKee et al. 2018).
On the theoretical side, Pen & Levin (2014); Simard & Pen (2018) proposed folded sheets of scattering material where anisotropy is a result of inclination and that could resemble corrugated reconnection sheets of plasma at the boundary between magnetic field configurations. Gwinn (2019); Gwinn & Sosenko (2019) proposed magnetic noodles of plasma whose anisotropy is maintained by parallel magnetic fields. This paper is organized as follows. An introduction to scintillation and PSR B1508+55 is given in Section 1. The data is described in Section 2 and Section 3 gives an overview of the observed phenomena. Sections 4 to 6 contain the three methods used to extract information from the data, and whose physical interpretation is visu- dynamic spectrum ( D , ) secondary spectrum ( ) complex amplitude distribution of images Table 1. Glossary of used variables. All vectors are two-dimensional because the dimension along the line of sight has a negligible impact on scintillation. An extended notation is used for the two-screen theory in Sections 7 to 9 and defined in Appendix A and Figs. 12 and A1. alized in Fig. 1. To explain these, a two-screen model is introduced in Section 7. Computational methods are described in Section 8 and the results based on them are presented in Section 9. We discuss ramifications of our work, and future avenues for research in Section 10.

Definitions and basic Theory
A glossary of the variables used is given in Table 1.
Scintillation inflicts a complex amplification on the observed electric field ( , ). Since the incoming radiation is not temporally coherent, its phase contains no usable information. Hence, its square modulus ( , ) = | | 2 ( , ) -the dynamic spectrum -is used. The dynamic spectrum still contains information that is not originating from scintillation in form of the intrinsic amplitude of the signal. As a result, the mean amplitude cannot be used for this analysis. Furthermore, intrinsic variations, such as in particular pulse-to-pulse variations, present a source of noise on the scintillation signal.
The secondary spectrum ( D , ) = |˜( D , )| 2 is a function of the Fourier conjugates of time and frequency, which are called Doppler rate D and delay respectively. These names are due to their geometric interpretation as differential Doppler shift and temporal delay of a signal. For a given lensing structure at angular coordinate -called image henceforth -from the line of sight, these coordinates obey a parabolic relation to each other. This is the reason for the appearance of scintillation arcs. One way to parameterize these relations is Here, is the speed of light, eff is the effective distance and eff is the effective velocity. The effective picture of these quantities is an observer moving with eff and rays originating on a plane with distance eff (see Sprenger et al. 2021, Figure 3). They can be translated to the distances and velocities of the Earth ( ⊕ ), the scattering screen ( scr ), and the pulsar (  Techniques not used in this work include direct measurement of the impulse response function using giant pulses (Main et al. 2017) or cyclic spectroscopy (Demorest 2011;Walker et al. 2013;Dolch et al. 2021). Instead, we will focus on secondary spectra, which map out the power belonging to structures on the screen (Walker et al. 2004, Eqs. (1.1) and (1.2)).

PSR B1508+55
The pulsar B1508+55 (J1509+5531) is located in the constellation Draco and circumpolar at the Effelsberg telescope and the Low Frequency Array (LOFAR). It was discovered by Huguenin et al. (1968) with a spin period of 0.739681922904(4) seconds (Hobbs et al. 2004) and Stovall et al. (2015) measured a dispersion measure (DM) of 19.6191(3) pc cm −3 . No orbital companion is known. Chatterjee et al. (2005) found that B1508+55's perpendicular proper motion is one of the largest known with the most recent measured value 963 +61 −64 km s −1 by Chatterjee et al. (2009). Origins of such high velocities were explored by Gvaramadze et al. (2008). Chatterjee et al. (2009) determined this velocity by measuring the pulsar's distance psr and its angular proper motions RA and Dec in right ascension and declination respectively: psr = 2.10 +0.13 −0.14 kpc , Before the establishment of secondary spectra as the primary tool to study scintillation, the characteristic time and frequency scales of this pulsar's scintillation have been measured in a number of studies by Rickett (1970), Armstrong & Rickett (1981), Lyne & Smith (1982), Smith & Wright (1985), Bhat et al. (1999a,b,c), and Kondratiev et al. (2001). These studies showed that the scintillation scales of B1508+55 are unstable and strongly varying between observations. More recently, this result was summarized and confirmed by Wang et al. (2018b,a) who also published dynamic spectra at high frequencies of 2250 MHz. Kaspi & Stinebring (1992) and Esamdin et al. (2004) studied the long-term flux-density variation of the pulsar, where the former found an unusually long refractive time scale that was not present anymore in the data of the latter.
The first secondary spectra of this pulsar were then published by Stinebring (2007b). They already show the same visual characteristics as the earliest observations done for this work. In particular, parallel stripes seem to be convolved with the scintillation arc. This feature does not belong to the expected phenomenology produced by a single thin screen and is thus the first major motivation for this study.
The second major motivation for this study are the results of Wucknitz (2018) who studied echoes in the pulse profile at 150 MHz observed with LOFAR. These echoes could be directly imaged via very long baseline interferometry (VLBI) revealing linearly aligned structures of one arcsecond in total length. The pulsar is moving along these structures and was expected to cross two brighter ones at some point in 2021 and 2023 respectively. Bansal et al. (2020) independently observed the echoes below 100 MHz using the LWA1 and also concluded scattering in the ISM as their origin.
In the most recent study of B1508+55's scintillation, Marthi et al. (2021) observed at 550-750 MHz simultaneously at the Giant Metrewave Radio Telescope (GMRT) and at the Algonquin Radio Observatory (ARO). Their secondary spectra also show the peculiar parallel stripes and they confirmed consistent movement of features between observations. Importantly, they found a low scintillation correlation between the sites as well as a travel time of the scintillation pattern that was much lower than the scintillation time scale and concluded the presence of a second screen in addition to the one causing the echoes.

OBSERVATIONS AND DATA
Our observations are from the Effelsberg 100-m telescope. While the main observing campaign started in March 2020 using the Lband receiver P217mm 7-Beam (18-21 cm prime focus receiver), the settings were chosen under the impression of five early observations in November 2019 and January 2020 of which some used a different receiver -the ultra broadband receiver (UBB) -that also covered higher frequencies of up to 2600 MHz. This study includes 76 distinct observations of 247 hours in total.
The settings used for all following observations and the reduction of baseband data were optimized such that both axes of the secondary spectrum covered approximately twice the average region of signal, asserting sensitivity to a wide range of curvatures and possible outliers of higher delay. The usable band was 1.27-1.45 GHz. The data were folded with 128 phase bins, 12.5 kHz channels, and 10 second subintegrations with DSPSR (van Straten & Bailes 2011).
The folded data were then read to a 4D data cube of phase, polarization, frequency and time using PSRCHIVE (van Straten et al. 2012), of which only the total intensity was used to continue. After dedispersing, the bandpass is determined as the median over time of the phase-averaged intensity and corrected for by division of the data by the bandpass.
The following recipe for noise removal was used: The off-pulse region was selected as the phase bins of the averaged pulse profile whose value is below the median. The standard deviation of that off-pulse region as a function of time and frequency is then used as a measure of radio-frequency interference (RFI). If this value is more than 10 above the mean corrected for outliers, it was masked. After that, the remaining noise was estimated by the mean over the off-pulse region for each bin of time, frequency, and polarization and subtracted accordingly.
Dynamic spectra were then created by multiplying along the phase axis by a new mean pulse profile created from the corrected data and subsequently summing over that axis. Since the overall intensity is irrelevant for studies of scintillation on the time scales regarded here, the dynamic spectra were finally divided by their standard deviation, which creates typical intensity values in the range of single-figure numbers.
As can be seen in Eq. (1.1), the Doppler rate of a fixed angle depends also on frequency, which leads to broadening of scintillation arcs for wider bandwidths. To avoid this, the Fourier transform over time is rescaled to a common reference frequency. This method was introduced as the NuT transform in Sprenger et al. (2021). It has the additional advantage that all spectra could be scaled explicitly to a frequency of 1400 MHz, avoiding errors in curvature measurement due to varying relative signal strength of different frequencies. Furthermore, to compute secondary spectra, the mean of dynamic spectra was subtracted, which reduces the influence of power in the centre. Figure 2 shows the secondary spectra of all observations from Effelsberg used for this study. To highlight the prominent phenomena and their evolution, characteristic examples of dynamic and secondary spectra are shown in Fig. 3.

Weak and Strong Scattering
The terminology used here when referring to the strength of scattering needs to be distinguished from the strength of scintillation. Both terms have been used synonymously in the past (e.g. Rickett 1990). However, there are two measures of strength in the context of scintillation such that we chose to explicitely discriminate between them as follows. Weak scintillation manifests in intensity variations which do not modulate the intensity significantly compared to the mean intensity and results from weaker phase variations on the screen. Strong scintillation results in strong interference patterns that can completely extinguish the signal. PSR B1508+55 clearly exhibits strong scintillation. The further division into weak and strong scattering is essentially defined over the existence or non-existence of arclets in addition to the main arc within secondary spectra. Weak scattering is weak because the scattered images are much dimmer than the central direct line of sight, such that only interferences between this central image and the outer images are observable. In the case of strong scattering, interferences between all images are observable, without interference involving central images being orders of magnitude brighter.

Arclets and Stripes
As described above and first explained by Walker et al. (2004), arclets arise from the interference between images that are not at the centre Figure 2. Secondary spectra obtained from observations of B1508+55 from the Effelsberg 100-m telescope for this study revealing a significant change in phenomenology in autumn of 2020. The dates are given in day, month, and year in this order. For this plot only, all secondary spectra have been rescaled to the same level of noise. In all plots of secondary spectra in this work, pixels below the resolution of the plot have been combined by taking the mean. The color scale is logarithmic and spans roughly four orders of magnitude. The spectra have been zoomed in -the full axes for most data extend to 50 mHz and 40 s. of the distribution. Phenomenologically, they manifest as downward parabolae whose apexes align on the main arc which has the same curvature as the arclets. In this case, Eqs. (1.1) and (1.2) for the locations of power need to be updated to include the two angular coordinates in question: The features we refer to as stripes are elongated parallel distributions of power in the secondary spectrum whose positions are aligned along a parabolic arc. Stripes were already observed for B1508+55 by Stinebring (2007b). They are an unusual feature which has rarely been observed. Stinebring (2007a) found similar stripes for B0355+54 and Smirnova et al. (2020) observed B0834+06 showing a similar state of scintillation for a limited amount of time. They are interesting for the general understanding of pulsar scintillation because they are highly ordered structures that are not predicted naturally by the single thin screen theory without finetuning of two-dimensional structures that would not be expected to be stable over such a long time. On the other hand, they are regular and stable enough such that the mechanism Figure 3. Examples of dynamic spectra (top row) and secondary spectra (bottom row) that have been chosen to highlight the different states of scintillation. The dynamic spectra have been cut only for the plots to produce comparable scales and the Doppler rate axes have been scaled such that the widths of the arcs are the same. The color scale is linear for the dynamic spectra and logarithmic in powers of ten for the secondary spectra. behind them should involve few enough degrees of freedom to be accessible to measurements.
Additional interesting aspects of the stripes are that their width is not observed to depend on the timespan of the dynamic spectrum as might be expected from a moving central bright image. Also, although the stripes of a single spectrum look very similar to each other, they cannot be simply deconvolved into a thin parabola and a single stripe. The stripes are thin and completely parallel in this observing band. Their brightness distribution can be shifted from the main arc but shows only minor scattering compared to the offsets of the other stripes. Their length has no clear cutoff and depends on the brightness of the individual feature and the data quality. However, beyond the overall brightness, their length does not depend on the location in the spectrum.
Although the appearance of stripes correlates with (almost) achromatic temporal variations in the dynamic spectra, these are unlikely to be a mere effect of pulse-to-pulse variations on larger time scales because they are missing in dynamic spectra of later observations where stripes and arclets coexist. Other reasons will be discussed below by showing that different scattered images of the pulsar exhibit different intensity variations, which contradicts an origin at the level of the pulsar's emission.

Transition to new State of Scintillation
The secondary spectra of B1508+55 observed until September 22 at Effelsberg differ strongly from those observed on November 10 and later at the same telescope (see Fig. 3).
Before the transition, the secondary spectra showed stripes. As there are no arclets, the initial state of B1508+55 can be described as a case of weak scattering. Many of the observations even show a secondary arc. It becomes more visible the narrower the parabola of the primary arc. Rarely, hints of a third and fourth arc can be seen that are even wider than the first two ones.
Intriguingly, the first observation after the transition already shows clear and strong arclets which means the pulsar now experiences strong scattering. The stripes are completely gone at first. In later observations however, flat features reappear at the apexes of arclets. Nevertheless, the pulsar remains in a state of strong scattering. These additional stripes, which coincide with an overall blurring of all other structures, are only visible during times of a lower width of the arc. Due to the annual cycle of arc curvatures, this is happening during the summer months of each year (see Figs. 2 and 3).

SCINTILLATION ARC EVOLUTION
The Earth's velocity vector changes its relative angle to the line of sight to the pulsar during the course of a year. Depending on the scattering geometry, this changes the time scale of scintillation (see Eq. (1.1)) and hence alters the width of scintillation arcs. Measuring scintillation arcs in multiple observations during different times of the year can thus be used to infer a geometric model of the scattering material. Variations can also be caused by orbital motions of the pulsar. In this case they provide constraints on orbital parameters. Recently

Modeling
Following Eqs. (1.1) and (1.2), a one-dimensional (1D) distribution of scattered paths results in a parabolic distribution of power = 2 D . The curvature parameter of scintillation arcs is given by (4.1) eff, is the effective velocity vector defined in Eq. (1.4) projected onto the 1D distribution of images, i.e. the screen's axis. The distribution of power can be more complicated, e.g. a filled parabolic envelope (isotropic screen), downward arclets (interference between outer images, see Eqs. (3.1) and (3.2)), and multiple parabolae (multiple screens).
The curvature has become the standard parameter to characterize scintillation arcs although it has some disadvantages. First, it is frequency dependent. To get a single number characterizing the scintillation, it needs to be scaled to some reference frequency, which is typically chosen to be 1400 MHz. Furthermore, is suboptimal for fitting since the region covered by the corresponding parabola within the secondary spectrum does not shift linearly with this parameter but rather with 1/ √ ∝ | eff, |. In particular, diverges in the important region around zero effective velocity. Finally, the parameter has no recognizable physical meaning by itself, which would help characterize a measured scintillation arc by labeling it with this number. For all these reasons, we use a different parameter than the curvature in this work: The physical meaning described by the temporal derivative of the square root of delay -which we refer to as the effective drift rate -is explained in Section 5. The parabolic equation now becomes (4.3) This notation naturally leads to a linear dependence in D -√ space which will be utilized throughout this work.

Weak Scattering with Stripes
The stripes complicate the arc measurement significantly because their power distribution varies between observations and is not necessarily centered on the main arc which creates an asymmetry between the left and right part of the arc. Methods that rely on summing the intensity along hypothetical parabolae like employed in Reardon et al. (2020) or on the identification of brightest points like in Main et al. (2020) are not optimal in this case.
To avoid systematic errors caused by the lack of proper modeling of the stripes, we chose to measure the arcs manually, which is feasible with regard to the number of observations. In order to do this, the square roots of secondary spectra are transformed from D -space to D -√ space for two reasons: First, this transformation visibly increases the signal to noise for the important regions of higher delay. By taking the square root, this transformation also achieves a uniform distribution of signal per area on the screen when pixels of the data before transformation are summed accordingly into pixels after the transformation (compare similar transformations in Sprenger et al. 2021). Second, parabolae get transformed to straight lines which are  Example of a wavefield obtained by phase retrieval as described in Baker et al. (2022). On the left side the corresponding secondary spectrum is shown. This data was taken on December 6 in 2021. The color scale is logarithmic and different for both plots. more efficiently identified by the human brain. The same method was applied to the secondary arc in the cases it was sufficiently visible. Figure 4 shows an example for such a manual measurement. The delay axis is resampled by summing over all pixels in that fall into a new pixel of √ while correctly accounting for fractional coverage. The color scale is logarithmic. The linear relation searched for becomes (4.4) Since the uncertainty of this method is dominated by unmodelled systematic deviations, a formal measure would underestimate the uncertainty. Instead, we chose the range in which the arc could reasonably be located by eye and used this as the standard error to be conservative.

Strong Scattering
While the presence of arclets introduces similar problems as the stripes for algorithms based on the distribution of intensity, it enables the option of using the eigenvector decomposition method introduced by Baker et al. (2022) in -space (Sprenger et al. 2021) which in principal can lead to a much higher precision than previously mentioned methods.
They showed how it is possible to utilize also the phases of the conjugate spectrum that are usually dismissed by using its squarethe secondary spectrum. As a result, sufficiently stable eigenvectors can be found even from small subsets of the dynamic spectrum, such that arc curvature estimation and phase retrieval is possible. We applied this technique of phase retrieval to our data of B1508+55 and found convincing results for subset sizes in pixels of The change of scintle width was taken into account by correcting according to the manually measured scintillation arc width. Other parameters were the range of delays < 15 s up to which corresponding angles are regarded for the eigenvectors and the number = 151 of discrete steps over which to evaluate the eigenvector. These numbers where the minimum found such that artefacts and offarc power were minimized in the wavefield that was computed from the square root of the dynamic spectrum with the retrieved phases applied. An example for an obtained wavefield is shown in Fig. 5.
However, we found that the method of finding the curvature belonging to the strongest possible eigenvector as described in Baker et al. (2022) deviated systematically from measurements done with the manual method described above. This deviation is strongest when the secondary spectra show reappearing stripes in addition to the arclets. We thus conclude that this method is also compromised by the peculiarity of this pulsar. Instead of resorting fully to manual measurements, we modified the eigenvector method to be more robust as described in the following. The original method relies on finding the largest eigenvalue over a distribution of curvatures evaluated on a small subset of data and then measuring the curvature by fitting those curvature measurements over frequency. In this case, individual eigenvalue distributions mostly did not follow a single peak distribution. Thus, we changed the order, sampled over the effective drift rate √ instead of the curvature to remove the frequency dependence and summed all eigenvalues before searching for the peak. Examples of the results are shown in Fig. 6. Still, this method remained biased because of the presence of a second peak that has no physical counterpart in the form of a secondary arc. This second peak is particularly prominent at lower numbers of than optimal for phase retrieval but also less wide such that its effect can be decreased slightly. Still, the main peak usually lies within the 1 range of the manually determined result. Hence, the mean of both results was taken and the standard deviation of both measurements taken as the uncertainty of the combined measurement. Whenever the eigenvalue peak was not compatible with the manual measurement only the latter with its error was used. In cases where the deviation of both methods was smaller than the uncertainty of a parabola fitted to the peak, the uncertainty obtained from the latter was used.

Additional Public Data
Two recent publications contained scintillation arc data close to the time of the Effelsberg observations of this study and at close enough frequency bands such that they could be added to the analysis. Marthi et al. (2021) report data taken from the GMRT over 550-750 MHz on August 6 and August 12 of 2017. They measured an arc curvature of ≈ 0.76 s 3 for both observations. For both dates, we added this value to our analysis after applying the correct frequency scaling and estimating the uncertainty as the mean uncertainty of the Effelsberg weak scattering arc measurements since the data looked qualitatively similar. Wang et al. (2018b) report data taken from the Jiamusi 66 m telescope over 2180-2320 MHz on December 13 of 2015 and October 31 of 2017. The corresponding dynamic spectra are publicly available. We took them from Wang et al. (2018a) 1 and the secondary spectra we formed from these dynamic spectra revealed measurable scintillation arcs. The latter observation even shows a secondary arc as bright as the primary arc. We applied the same arc measurement technique as was used for the weak scattering Effelsberg data.

FEATURE ALIGNMENT
The features of a scintillation arc that manifest in distinct blobs, arclets, or stripes can be translated to individual images of the pulsar. If the same features appear in several subsequent observations, they can be aligned, which provides a measurement of the movements of these features along the arc with respect to its origin. Hill et al. (2005) showed that the angles corresponding to different individual arclets o B0834+06 change by the same angular velocity between observations. They concluded that the scattering structures in the screen are approximately constant. Simard & Pen (2018) computed the angular movement of images for a specific model of the screen and predicted non-linear movement and even truncation of the images around the origin of the scintillation arc.

Modeling
The time derivative of the crossing point 0 with respect to some fixed location within the screen can be derived from geometrical arguments (see Fig. 7). As an angular motion it is the physical velocity ISS of This term equals the relation used by Brisken et al. (2010) and was branded scintillation speed by Cordes & Rickett (1998). Conveniently, the scintillation speed can be multiplied by psr /( psr − scr ) to obtain the effective velocity and the remaining distance factors combine such that the result can be written in terms of effective quantities as defined in Eqs. where 0 is the point in time where is defined as the current position of the image relative to the line of sight. However, if screens are indeed close to one-dimensional, images will never stay at a fixed position on the two-dimensional screen but rather move such that the line of images is never offset from the direct line of sight to the pulsar. This means that in practice only the angular coordinate relative to the screen matters and eff has to be replaced by its projection eff, onto the screen's axis. This result can be shown to be identical to the formula that has been used by Main et al. (2020); Marthi et al. (2021). As mentioned earlier, the effective drift rate √ is used in this study. The reason for using the temporal derivative of √ instead of the temporal derivative of D is to avoid a frequency dependence. Both quantities are linearly proportional to the angular position on the screen. Thus, the effective drift rate can be derived from Eqs. (1.2) Figure 8. Results of feature alignment using the peaks of template correlation. From the two highest peaks within a fixed range, the one closer to the mean of the results from arc curvature measurement of the two compared observations was chosen. This method is described in detail in Section 5.2.1. and (5.1): (5.7)

Weak Scattering
In the case of weak scattering with stripes, aligning features across different observations is eased by the fact that the stripes are horizontal in delay and thus a feature can easily be identified with a single angular coordinate. To be model independent, angular coordinates should be expressed in √ . Since the features move along the arc, it is beneficial to shift the zero axis in Doppler rate to the centre of the arc at each delay. This step is only done to ease the comparability between observations. Although the centre of the arc depends on a measurement of that arc, this is not introducing a bias to the result of the alignment because the angular coordinate of a feature is not affected by this shift. The square root scaling of the delay can further be used to increase the signal to noise for high delay feature by summing the pixels such that a uniform pixel width is reached in √ -D space is obtained.
To measure the effective drift rate √ between two observations, at least one common feature needs to be identified. In the case of weak scattering we employed a template correlation algorithm to find the best fitting shift. We chose a subset of the data with high enough cadence and then created √ -D diagrams as described above. All diagrams are created in the same way: (i) Take the amplitude of the conjugate spectrum, which equals the square root of the secondary spectrum.
(ii) Compute coordinates in 800 equal steps in √ ranging from 8 to 0 s as − √ 8 to 0 √ s following the arc downwards for negative Doppler rates and from 0 to 8 s as 0 to √ 8 √ s following the arc upwards for positive Doppler rates.
(iii) Using the measured arc curvature, compute coordinates in Figure 9. Aligned spectra according to matching features. Shown are the coordinates along the main arc against the modified Julian dates of the observations. Top: During the weak scattering regime it sufficed to only transform one axis, such that the unlabelled coordinate along the y-axis is D . Bottom: In the strong scattering regime the data had to be transformed to √ -√ space (model-independent -) to transform arclets to lines. The contraction of the y-axes increases the remaining distortions due to reappearing stripes. D ranging from −2 mHz to 2 mHz in 40 equal steps after shifting the zero point of D to √ / √ .
(iv) Compute each pixel's value by taking it's range in D and and summing over those values counting partially covered pixels by their fractional area covered.
(v) Only for the measurement of the the effective drift rate √ , divide each diagram by the sum of all diagrams to correct for common artifacts and the brightening of features once they approach the origin.
The shift between two observations is measured as follows: (i) In the first observation, take a region centered on the arc with 10 pixels width in D and 90 pixels width in √ as a template to correlate with.
(ii) The correlation with each region of the same size in the second region is computed. For each value of shift in √ the shift with the highest correlation in D is chosen to correct for slightly wrong arc curvature values.
(iii) The above steps are repeated for the next template shifted by 10 pixels in √ . The new correlations are summed onto those already obtained.
(iv) For the next step the search range of shifts is reduced to shifts corresponding to curvatures that do not deviate from the mean curvature of the two observations by more than 30%.
(v) Within this region the highest two values are identified with the additional constraint of the second one being at least 5 steps apart from the highest one.
(vi) If there is only one peak, it is chosen as the best fit value. If there are two peaks, the one closer to the mean curvature is chosen.
The method of aligning features works for any pair of observations that is close enough in time to have detectable overlapping features. Here, we included every pair of observations that was separated by three weeks or less and obtain a point of data for each combination. Figure 8 shows the results of this method compared to the allowed range and the values from the arc measurement. To use this data, its uncertainty has to be determined. This is not straightforward as the probability of two spectra being a match is distributed in a non-Gaussian way over all possible shifts. To estimate the uncertainty a posteriori instead, we subtracted a sinusoidal fit with a period of one year from the data and computed the uncertainty such that the variance is 1 2 which corresponds to a normal distribution. The downside of this procedure is that the uncertainty estimate is modeldependent, which has to be kept in mind if a model is considered which is not linearly dependent on a velocity vector. This concern does not apply to this study though. Moreover, the data points are not completely independent because there cannot be more independent information than the number of observations because each of them resembles a fixed shift relative to a previous point in time. Fig. 9 shows the aligned √ -D diagrams according to the measured shifts of adjacent observations. For the purpose of this plot the measured shifts from non-adjacent observations were ignored. In addition to the obtained effective drift rates, the alignment allows for two additional conclusions: First, the images fluctuate between observations, which means their amplitude is neither constant nor completely determined by their distance to the origin. Second, images can be identified over a long time, even after crossing the origin. This places constraints on plasma screen models because it requires quite different values of phase gradients to be present in close vicinity on the screen.

Strong Scattering
The parabolic arclets present in strong scattering pose a challenge for feature alignment because each feature is spread over a region in Doppler rate and delay that cannot be identified with a single angular position without using a measured curvature. Furthermore,in some observations stripes reappeared in addition to the arclets, which further distorted the position of features.
To make an alignment possible, the secondary spectra are transformed to incoherent -diagrams as described in Sprenger et al. (2021). Correlation peaks obtained by the technique used for the weak scattering data scattered randomly in the allowed range and thus allowed for no measurement. Therefore, we aligned the diagrams manually by plotting them next to each other and determining the best shift.
The obtained alignment is shown in Fig. 9. The uncertainty is even less constrained than in the weak scattering case because of the significant differences of the quality of alignments of pairs. Thus, we chose only to regard the results from the alignment of adjacent observations. The uncertainty depends on two factors: Firstly, there are noisy and perhaps also physical variations of the location of features relative to each other over time. This is estimated to be a constant number. Secondly, the uncertainty grows with diminishing overlap of observations. This is estimated as a linear dependence on the time between observations. The total error is the quadratic sum of these two errors. To fix the two free parameters, the constant error is estimated to be equal to the error obtained for the correlation method on the weak scattering data. The linear error is then fixed by again 2 = 1 scattering around a sinusoidal fit.

EIGENVECTOR DECOMPOSITION ANALYSIS
Eigenvector decomposition is a powerful tool to obtain the brightness distribution of images on the screen. This technique is based on the -transformation that was introduced by Sprenger et al. (2021) and transforms the secondary spectrum into a space where it, in the case of a single and one-dimensional screen, resembles a matrix equal to the outer product of the complex amplitude distribution ( ) as a function of the angular coordinate along the screen's axis. It has been employed in Section 4.2.2 to measure the curvature of scintillation arcs.
This method can also be used to study the obtained brightness models for different parts of the data, as demonstrated by Rickett et al. (2021). In our case, the small required subset size allowed us  Fig. 10) of data taken on April 12 in 2020 before the transition. Top right: Corresponding secondary spectrum. Below are shown in this order: The corresponding dynamic spectrum and the contributions to it obtained from inverting the Fourier transform after masking around the regions A and B in the conjugate spectrum. While the slope of modulation cannot be measured, a variation of the magnification of different images is clearly visible and shows that the almost achromatic structures in early dynamic spectra are caused by this modulation.
to obtain a huge set of eigenvectors for every observation. They can best be visualized by taking the mean of their amplitudes over time or frequency, in order to suppress noise. An example of such is shown in Fig. 10. In order to increase resolution, < 15 s and = 601 was used for this case. To assert statistical independence, the data subsets were chosen not to overlap here.
The eigenvector analysis reveals another special phenomenon of B1508+55: The brightness distribution is affected by a modulation that travels along the screen over time. The velocity of this modulation can be measured. To do this, we divided the frequency-averaged eigenvectors by their temporal mean. Then, the Hough transformation is applied to the eigenvectors over time to obtain a diagram over the -intersect at = 0 and the slope of the modulation. Since not only lines of power but also lines of absent power exist, the most efficient way to find the best fitting slope is to take the peak of the standard deviation over the -intersect as a function of slope. The result is quite stable for all observations after the transition, allowing only for annual variations of less than 0.05 √ /h. If a constant value is assumed, the result for the modulation speed of 0.6560(56) √ /h can be obtained by taking the mean and using the standard deviation divided by √ as its uncertainty. An alternative way of measuring this eigenvector modulation is to mask the arclet or stripe belonging to a single image, invert the Fourier transform, and average the amplitude of the result over fre-quency. This method yields the temporal evolution of eigenvectors without regard for the phases but slightly more stable in the case of weak scattering. Still, before the transition, the quality of eigenvectors is too poor to measure the modulation speed. Nevertheless, there are indications that the modulation is already present. An example is shown in Fig. 11.
Since the modulation is traveling along the full visible screen within hours and is multiplicative and thus separable from the structure on the screen, its physical source is most likely located behind the screen causing the main arc. A modulation like this would naturally occur from a second screen if the extra delay is negligible such that there is no variation in frequency.
The vertical structures in the dynamic spectra of weak scattering could then be interpreted as the modulation of the bright central image of weak scattering. These structures correspond to the stripes in Fourier space. Thus, the stripes and the eigenvector modulation result from the same mechanism in this picture. Putney & Stinebring (2006) already reported on many pulsars with multiple scintillation arcs suggesting that scattering situations with more than one relevant screen are not unlikely. Two screens have been proposed as explanations for special scintillation phenomenafor example by Liu et al. (2016). Here, the evidence for two screens is very strong because such a theory consistently explains the special phenomena of stripes and eigenvector modulations (see Section 6).

TWO-SCREEN MODELING
Systems of multiple lenses have been studied theoretically by Er et al. (2022);Feldbrugge (2020). However, the stability of image positions we observed and the high number of images in the brightness distribution both imply that the lensing region of a single image is much smaller than the observable part of the screen. Thus, phenomena predicted by those studies on the scale of lenses that are resolved enough to be described analytically might not be accessible to our observations. On larger scales and without knowledge on the exact distribution of phases imposed by the screen, we can only rely on a number of data-driven approximations instead: (i) Anisotropy: The phase screen is very anisotropic, i.e. the phase varies exclusively in one direction on the screen.
(ii) Stability: Along this axis of anisotropy, there are only certain points of stationary phase that allow propagation through the screen. On the large scales of interest, these points are stable in position over time and frequency.
The second point can be interpreted as follows: An image exists in a region where the mean phase gradient is zero. On smaller scales within this region the distribution contains many slightly different phase gradients such that the combined emission of the whole region remains stable even if small changes in the geometry happen because of movements of observer and pulsar or slightly different frequencies used.
In Appendix A, we derive a double screen theory based on these approximations in order to model the phenomena we observed. Fig. 12 shows the definition of distances and velocities in a two-screen geometry.
On each of the screens, propagation is limited to lines of scattering material with a suitable phase gradient. In the presence of a single screen, the dimension along these parallel filaments can be ignored because the radiation will always follow the direct line of sight in dimensions without variation in phase delay. In the case of two screens observer first screen second screen pulsar Figure 12. The geometric situation of two scattering screens. The green dots show the points where the screens were intitially crossed while the blue circles show the current crossing point of the ray (red line) after some movement indicated by arrows. The green and yellow line represent a line of constant dispersive phase on the screens at initial and current time. In addition to the physical movement of the scattering material, the crossing point moves along these lines to the minimum of the total phase gradient. whose axes of anisotropy are tilted to each other, both dimensions have to be taken into account. As a result, scintillation arcs and feature movements are the product of mutual interaction of the geometries of both screens. In addition, the screen closer to the Earth is illuminated not by a single point source but by multiple paths of propagation through the farther screen. Thus, a scintillation pattern is imprinted on the closer screen which can be observed by a variation of the amplitudes of images on the closer screen. If the variation in frequency can be neglected, the result is an eigenvector modulation like described in Section 6, whose modulation speed can be derived geometrically.

Modelling
We use coordinates in the equatorial coordinate system and define angles of velocities and screen orientations with respect to the axis of positive right ascension in mathematically positive direction on the sky. For each scattering screen in the model there are three free parameters: The distance , the angle of anisotropy and the local ISM velocity , along this axis of anisotropy. The direction of positive angles on the screen is arbitrary such that the same model is described when changing by ±180 • and inverting the sign of . We use Astropy (Astropy Collaboration et al. 2018Collaboration et al. , 2013 to compute the velocity of the Earth projected on the sky towards PSR B1508+55 within this coordinate system. Accounting also for the additional velocity component from the rotation of the Effelsberg 100-m telescope around the Earth was found to be slightly below the measurement uncertainty and thus dismissed.
The secondary arc may relate to another screen than the one responsible for the eigenvector modulation. However, we will begin with the simplifying assumption that there are only two screens, and show that all of the observed phenomena in the data are fully consistent with this model.
The annual variation of the primary arc apparently changed consistently in phase at the point of the transition from weak to strong scattering (see Section 9). There is no reason to assume that both screens changed properties simultaneously such that this shift of phase is also very likely originating from the first screen. The natural candidate for such a shift is the angle of anisotropy. A direct result of a shifted angle of anisotropy is a sensitivity to the second component  Table 2. Used formulae to model different data types using models of one or two screens. The relations for two screens are derived in Appendix A.
of the velocity vector because the scintillation depends on the projection of the local ISM velocity onto the axis of anisotropy. Thus, the minimal model has eight free parameters: The distance of the first screen, the distance of the second screen, the angle ,WS of the first screen before the transition, its angle ,SS after the transition, the angle of the second screen, the amplitude of the ISM velocity at the first screen, its angle , and the projected ISM velocity , at the second screen. Additional three parameters that will be marginalized over are the pulsar's distance and its velocity vector on the sky.
In addition, a one-screen model is explored for comparison. Again, only the angle of the screen is allowed to change with the transition. As this model does not give any prediction on the secondary arc or the modulation, these will be modelled by simple constants in order to compare the model to the full data. Thus, this model has seven free parameters.
In order to judge how well these models fare in general to describe the data, a simple heuristic model is introduced as a control. Here, three simple annual sine waves are used for the primary arc and its features before the transition, and separate ones for arc and features after the transition because there is a visible difference between those results. The secondary arc and the modulation are again modelled by constants. This model has eleven free parameters since each sine wave has three free parameters: √ = sin [2 ( − 0 )/(1year)] + . (8.1)

Parameter Inference
We use the Markov Chain Monte Carlo (MCMC) ensemble sampler emcee (Foreman-Mackey et al. 2013) for Bayesian inference of model parameters. The likelihood is constructed as independent Gaussian probability densities for each measurement. The modulation speed is lacking a reliable method of uncertainty estimation (see Section 6). Thus, a conservative assumption of using the standard deviation of all its measurements is used. This assumption is compatible with an approximately constant value, which is expected for a distant screen. In order to avoid double counting of information in the case of the pair alignments during weak scattering, the likelihood contributions of this data were weighted by the number of observations divided by the number of pairs. This is an approximate way to account for dependent measurements. For the parameters psr , , and Gaussian priors are created according to the larger uncertainty of the two-sided uncertainties in Chatterjee et al. (2009) (see Eqs. (1.5) to (1.7)). A Gaussian prior with a standard error of 100 km/s is applied to all ISM velocities in order to keep values within a reasonable range of nearby star velocities as measured by Gaia Collaboration et al. (2018). For all other parameters flat priors are used. Table 2 summarizes the formulae used for different types of data for models of one or two screens.

Simulations
We perform simulations with the purpose of demonstrating that the theoretical framework described above can reproduce the observed phenomenology of B1508+55. Hence, dynamic spectra of same range and resolution as real observations are simulated. The parameters used equal the best fit values from the two screen model defined above.
The measurement of scintillation arcs, feature movement and modulation speed place no constraints on the distribution of images and their amplifications . To ensure comparability between the observations, the same image positions were used for all simulations.
Only on the first screen the amplifications of images were changed to replicate the change between weak and strong scattering. We found that the brightest central image could be kept in position and thus be treated in the same way as any other image to replicate the data.
Methods to infer information from image distributions are left for future work. Here only a proof of concept is desired. Thus, no attempt was made to exactly match the observed images.
Numerically, simulations were performed by computing the electric field as on a grid of time and frequency. For the phase Φ, Eqs. (A8) and (A15) were used. The dynamic spectrum is obtained by taking the square modulus, after which the secondary spectrum can be obtained as usual.
Using the same methods as for the real data, the eigenvector modulation and its speed can be measured for the cases where the first screen is strongly scintillating. Here, the measured modulation speed can be directly compared to the theoretical expectation. Figure 13 shows the accumulated results of all measurements of scintillation arcs, feature movements, and modulation speed. All observables can be expressed in effective drift rates. In addition, model predictions drawn from the posterior of the two-screen model are shown. The modulation speed and the secondary arc only show small annual variation. The primary arc before the transition clearly does not agree with the primary arc after the transition. While the feature alignment is compatible with the scintillation arc before the transition, there is a slight offset between them after the transition. The two-screen model is able to reproduce all of the measurements. Table 3 shows the best-fitting values for the parameters of the two-screen and one-screen model with the standard deviation of the posterior distribution in parenthesis. The 2 per number of degrees of freedom is also shown to compare the models. To provide a baseline in case of wrongly estimated uncertainties, the 2 /ndf is also shown for the heuristic model. The pulsar's distance and velocity are not given because they could not be constrained any further than their priors.

Model Fits
The two-screen model is the only model that gives a physical prediction for all observables. At the same time it fits the data better than the one-screen model and the comparison to the heuristic model does not suggest the existence of badly modelled parts of the data. Even the small deviation between the effective drift rates of scintillation arc measurements and feature movement measurements during the regime of strong scattering is covered by the two-screen model. The only other explanation for this shift would be an undiscovered systematic error as such a deviation is impossible in a one-screen model.
These numbers can be used to infer physical scales from the measurements. According to the best-fitting two-screen model an arc containing power up to 10 s of delay corresponds to structures within a range of two astronomical units and the modulation speed translates to a velocity of about 8600 km/s. The illumination pattern needs roughly 10 hours to travel along a screen like this.
If the data for the secondary arc is omitted in the likelihood, the best-fitting model still predicts values for the secondary arc that are compatible with the data, although less constrained. The same is true if the modulation is omitted in the likelihood. This finding provides further indication that the two phenomena originate from the same second screen, which can also be reproduced in simulations.

Echo Crossing and Screen Orientation
A useful side effect of a change of the direction of anisotropy is sensitivity to the 2D vector of ISM velocity at the first screen. This means that the point where the direct line of sight to the pulsar crossed the screen in the past can be reconstructed unambiguously. The result is shown in Fig. 14. After integrating the movement backwards to the time Wucknitz (2018) observed B1508+55 with LOFAR, the localized echo components can be compared with the observations at Effelsberg.
The result is an agreement with our model in three aspects. The transition coincides with the beginning of bright outer components. Furthermore, since an image forms where the phase is stationary, the orientation of constant scattering material has to be perpendicular to the line to the origin of the observation and can thus be computed. Right at the transition, this orientation changes gradually from the initial angle to the one after the transition. Finally, the distance of the screen as measured by Wucknitz (2018) is compatible with our result for the first screen. Hence, the echoes can be explained by an effect of the same scattering screen on a larger scale.
Noteworthy properties of the best-fitting model are that the two screens are nearly perpendicular to each other while the second screen is close to parallel to the pulsar's velocity vector. This might be a reason for the images on the second screen being stable enough over time such that two-screen phenomena are so prominent in the scintillation of B1508+55. Short-lived images would produce an effectively comoving central image like described at the end of Section 7 which would reduce the first arc's observables to being compatible with a one-screen theory. However, the modulation of eigenvectors would still be visible. Such cases may be realized in many other pulsars or in the faint additional arcs of B1508+55. Figure 15 shows simulations according to the best-fitting two-screen model. This is another important check if the two-screen model can potentially explain all the phenomena because the shift in appearences of stripes and arclets was not used for the measurements the model was fitted to. By changing the relative amplitudes of images on the first screen, it is indeed possible to reproduce stripes, secondary arcs and a shift from weak to strong scattering. Figure 13. Measured values for the effective drift rate for both arcs and the modulation as well as models according to parameters drawn from the posterior. All times are shifted to the reference year 2020. (WS) and (SS) indicate data taken before the transition (weak scattering) and after the transition (strong scattering). 400 models are shown to illustrate the uncertainty whereas each line has the same color as the data it belongs to, except for the Jiamusi and GMRT data which belong to (WS) arc data. On the right axis, the effective drift rates are translated to the effective velocity over square root of effective distance.  Table 3. Best-fitting parameters and standard deviation of posteriors of the two-screen and one-screen model as well as a heuristic reference model. The last column lists the reduced 2 . These fit parameters are defined in Section 8.1.

Comparison of Data and Simulations
Only for the observation on September 22 in 2020 one bright image had to be added to the second screen to reproduce the bright observed feature on the secondary arc. This is indicative of significant evolution of structures on the second screen, which is unfortunately impossible to quantify due to the secondary arc being fuzzy in all other observations. Even the reappearance of stripes on top of arclets during times of narrower scintillation arcs naturally occurs. This last effect does not even need a change of any of the image amplitudes.
One interesting observation is that to reproduce the last observation before the transition (September 22, 2020) it was necessary to make images on the incoming side much brighter than on the outgoing side. It is possible to interpret this as the incoming of a region of stronger scattering which later is responsible for the strong scattering and may even be identified with the onset of one of the echoes studied by Wucknitz (2018); Bansal et al. (2020).
As shown in Fig. 16, the simulations also reproduce the eigen-vector modulation. Since the dynamic spectra were simulated, and not the eigenvector modulation directly, the method of measuring the modulation speed can be validated.

CONCLUSIONS
We reported a series of observations of PSR B1508+55 using the Effelsberg 100-m telescope. These observations covered roughly two and a half years and showed a large variation of scintillation arcs. Furthermore, in the secondary spectra we observed a transition from a weak scattering phase with multiple arcs and hitherto unexplained stripe-like features to a strong scattering phase with clear arclets.
The aim of this paper is to develop a model of the scattering geometry leading to the peculiar scintillation of the pulsar B1508+55, and thereby to introduce a number of new techniques for data analysis. We leave the physical interpretation of the scattering material as well  Examples of simulated dynamic spectra (top row) and corresponding secondary spectra (bottom row) using the model predictions at the corresponding dates. Noise and pulse-to-pulse variation have been added to make the results compatible to Fig. 3. as possible implications on properties of the pulsar itself for future studies. It shall be mentioned, though, that there are candidates: An A2 star is located 1.37 pc from the line of sight to the pulsar at a distance of 120 ± 8 pc (Wucknitz 2018). Also, the edge of the local hot bubble lies at about 100-150 pc away from us in the direction of the pulsar (Liu et al. 2017). Both are consistent with the distance of the closer screen. For the more distant screen, there is little known material along this line of sight that could be associated with it.
In addition to applying recently developed methods to measure scintillation arcs, we introduced methods to infer also independent information from the movement of features within the arc. The alignment of features over months showed consistent movement through the full arc, which is in disagreement with the model by Simard & Pen (2018) and indicates the presence of very small regions on the screen that contain a huge range of different dispersive phase gradients.
For the first time, we found a third observable by measuring the temporal evolution of the amplitude of images. This modulation speed could be measured for all observations that included arclets and provided a strong indication to interpret the scattering as a system of two screens.
To test this theory, we developed an analytic two-screen model. This model was then used to predict the observables as well as to simulate dynamic spectra. For comparison, we modelled all applicable phenomena also in the well-tested single-screen framework. Finally, we performed Bayesian inference on the obtained data and compared the models as well as the ramifications of the obtained parameters. We found that the best-fitting model agrees well with the VLBI results of echo components as observed by Wucknitz (2018). Marthi et al. (2021) had previously suggested a two-screen model to explain the dynamic cross-spectrum seen on a 10000km baseline; our results strongly support this claim, and have allowed to measure the distance and orientation of both screens.
Although the two-screen model is able to reproduce the observations we reported on as well as observations by others, much of the underlying physics remains to be uncovered. In particular, the coinciding occurrence of the echo crossing, the change of the axis of anisotropy to being parallel to the pulsar's motion, and the onset of brighter scattered images suggest a common physical origin. Furthermore, the existence of strong two-screen effects in a system whose second screen is perpendicular to the first screen and aligned with the combined velocity of the pulsar and the first screen might be more than a coincidence.
Future study of this pulsar will include the analysis of data taken at LOFAR at a much lower frequency in order to study the validity of the model under these circumstances and to further explore the connection of the transition to the echoes that were discovered at these frequencies.
The transition underwent by B1508+55 allowed us to observe its scintillation through the lens of two different sets of phenomena which enabled a robust test of a two-screen scattering model. Since this explanation worked here, it is reasonable to assume similar models for other scintillating pulsars. Furthermore, we discovered the phenomenon of eigenvector modulation, which is a consequence of two-screen scintillation. In this case it was so strong that it also affected the shape of the secondary spectra by introducing the stripe feature, but it may reveal itself also in sources which otherwise appear dominated by one strongly scintillating screen. We have shown how to combine measurements of scintillation arc curvatures, fea-ture movement along these arcs, and the speed of modulation along the features on these arcs and that the latter two can add a lot of information to the well-tested study of scintillation arc curvatures.
Scattering screens seem to be ubiquitous such that a situation of more than one scattering screen is not expected to be rare. In particular, double-lens scintillometry should be a useful tool for sources that have a local environment likely hosting a screen in addition to the screens expected closer to Earth, namely pulsars in Supernova remnants like the Crab and B0540-69, the galactic centre magnetar, and extragalactic sources like pulsars in the Magellanic clouds or fast radio bursts (FRBs).
where is the elementary charge and is the mass of the electron. Even with the above approximations being valid, the screen is still likely a very large structure. Hence it can be approximated as being frozen in time, which makes it convenient to separate the frame velocity of the screens from coordinates within those screens, i.e. = + . For the remaining part of this section we will use upper case letters for time dependent variables and lower case letters for coordinates on the frozen comoving planes. A sketch illustrating these notations is given in Fig. A1.
In its most general form, the phase inflicted on a path of propagation with respect to the unperturbed line of sight is given by and the electric field at the observer can be computed via the integral Here, is the total number of diffraction screens and , , and are the coordinate on the nth screen, its velocity, and its inflicted phase respectively. −1, is the distance between the n'th screen and the screen before it from the view of the observer. Screen 0 is the observers plane denoted with in the following while screen N+1 is the source plane denoted with .
The model we are interested in has only two screens. The one closer to the observer will be denoted with and the other one with . Furthermore, we assume these screens to be one-dimensional, such that they have the same high anisotropy that was assumed in observer first screen second screen pulsar Figure A1. The notation used in each of the four planes of a double-screen situation. Upper case symbols are coordinates fixed within the planes with respect to the initial line of sight, while lower case symbols are coordinates that are fixed to the moving object. The coordinate systems are aligned to the axes of anisotropy on the screens. In the case of Earth and pulsar, they are aligned to the closer respective screen, which eases derivations. Compare Fig. 12 for the notation of distances between the planes. the single screen model because the data shows no strong deviations of features from their main scintillation arc. This means the inflicted dispersive phase only depends on the coordinate parallel to the axis of anisotropy of the screen. Another argument for this model is that the apex of the parabola is always centered at the origin which means -as has been demonstrated by Shi (2021)-that the phase structure extends quasi infinitely perpendicular to the anisotropy's axis. Thus, this model is determined by the following phase: This expression can be reduced to the relevant coordinates and by performing the integrals over the perpendicular coordinates ⊥ and ⊥ . This amounts to integrals over a quadratic phase. Hence the stationary phase approximation yields the correct remaining phase up to constants and a prefactor, both of which we will again drop. Thus, we insert the following conditions: The result is a quite lengthy expression, so we will employ some abbreviations here. First, we concentrate on the partΦ that depends on the screen coordinates, such that