spedas / bleeding_edge

IDL-based Space Physics Environment Data Analysis Software (bleeding edge)
http://www.spedas.org
Other
7 stars 0 forks source link

Develop code to infer THEMIS-E Bz component from Bx, By and particle data #160

Open jameswilburlewis opened 5 months ago

jameswilburlewis commented 5 months ago

With the latest Bz offset changes on THEMIS-E, we need a way to estimate the Bz component using particle data. From Vassilis:

As we discussed we would like to determine the field line direction from the ESA electrons in order to recover the Bz DSL component from the other two components. The goal is to obtain the angles theta and phi (latitude, azimuth) of that direction based on a principal axis analysis of electrons of a given energy, and do this as a function of energy. Each spin period, some energies may respond well, some may have noise and not respond well. We can only keep the energies that do respond well, and average the components of the unit vector to obtain the best direction.

Once we have theta, then we can find Bz = tan(theta)*sqrt(Bx^2+By^2) in DSL every spin period and create a corrected B_dsl product at spin resolution.

To obtain theta at a given energy, we compute the partial pressure tensor at that energy Pij. This can be done for each channel or a combination of channels, 3 at a time centered at one. So, 15 energies will provide 15-2= 13 channels, each with 3 energies. You may have to throw away the lowest few channels due to photoelectrons, so lets say 10 channels. Then use principal component analysis to determine 3 eigenvectors and eigenvalues for each channel. Two eigenvalues would be very similar (the electrons are gyrotropic) representing the P_perp, while the other is most dissimilar. The dissimilar eigenvalue would be P_par, and could be either below P_perp (the electron distribution would be like a pancake) or above P_perp (a cigar distribution) or equal (isotropic distribution, carries no information on B-direction). If the most dissimilar eigenvalue is (for example) < 20% away from the average of the other two most similar ones, then we can throw that point away. At the same time,nnnnn if the most similar eigenvalues are > 10% different from each other (value is an example), then again we throw away that point. These are 2 criteria for selecting potentially good B-direction solutions.

The collection of N good points at different energies results in N eigenvectors in DSL coord’s each time (1<N<10). You can average the few eigenvectors that satisfy these criteria to obtain a single eigenvector. You then determine theta and phi from that.

From the theta, you can determine Bz = tan(theta)*sqrt(Bx^2+By^2) in DSL as discussed above. A third criterion can be applied to check the validity of the solution. That is the eigenvector phi angle compared to that from atan(Bx, By). You can use that to determine how good the 2 criteria above are, and how much error to tolerate and still get N> at each time, in enough cases.

Let me now if this unclear. IDL has a built- in function for principal axis analysis, which is same as is used in minvar.pro (I think) based on Bij. You could see how it is implemented tjere.

jameswilburlewis commented 3 months ago

From Adrian P: Description of TH-E Bz recovery procedure

THE - Spin Axis Recovery Only use data from the SSL-Z component as the angle to the spin axis is smaller. This could change in the future due to spacecraft maneuvers that lead to a change in the spin axis.

  1. Download data and convert raw data to nT
  2. Rotate offset from calibration file into SSL system and subtract from the data 1
  3. Fit cosine with the following form to the data for a window length of 2.5 periods: Fit_B = (A + Bt) cos(ω (t - ϴ)) + C + Dt (1) Use the spin frequency ω from the state file and fit all other parameters
  4. Build running average with weights w with a window length of 1 minute w = 1 / ( f min / A^2) with f min as the mean squared difference of the fit from the data f min = Σ(Fit – Raw_Data)^2 / n_elements(Fit)
  5. Calculate angle α between sensor axis and spin axis by comparing the fit values with the IGRF model near perigee
  6. Divide retrieved offset (C + D*t) by the cosine of the angle α

1 We may have to determine the offset in a different way in the future as we have only the values for the spin plane and not the spin axis. One possibility: Compare the fit values with the IGRF model in the outer magnetosphere (30 nT - 100 nT).

jameswilburlewis commented 1 month ago

Some sample sav files with estimated Bz data are now available in /disks/themisdata/fgm_bz.