Completed 1.0 documentation except KS2005 functions

PlanetMag.m

@@ -1,6 +1,6 @@
 function [T_h, B0vec, B1vec1, B1vec2, B1vec3, outFname, header] = PlanetMag(moonName, era, ...
-    coordType, CALC_NEW, ALL_MODELS, DO_FFT, DO_MPAUSE, LIVE_PLOTS, specificModel, specificMPmodel, ...
-    outData, fPatternFT, nptsApprox, magPhase_deg)
+    coordType, CALC_NEW, ALL_MODELS, DO_FFT, DO_MPAUSE, LIVE_PLOTS, specificModel, ...
+    specificMPmodel, outData, coeffPath, fPatternFT, nptsApprox, magPhase_deg)
 % Evaluates planetary magnetic field for a time series at the specified moon location and inverts
 % for the complex amplitudes of oscillation in that moon's frame.
 %
@@ -45,6 +45,8 @@
 %   model See MpauseFld for a description of the options.
 % outData : char, 1xF, default='out'
 %   Directory to use for output of complex spectrum amplitudes.
+% coeffPath : char, 1xE, default='modelCoeffs'
+%   Directory containing model coefficients files.
 % fPatternFT : char, 1xG, default='FTdata'
 %   Pattern for file names of FFT spectrum data saved to disk.
 % nptsApprox : int, default=12*365.25*3*24
@@ -87,6 +89,7 @@
     if ~exist('specificModel', 'var'); specificModel = 0; end
     if ~exist('specificMPmodel', 'var'); specificMPmodel = 0; end
     if ~exist('outData', 'var'); outData = 'out'; end
+    if ~exist('coeffpath', 'var'); coeffPath = 'modelCoeffs'; end
     if ~exist('fPatternFT', 'var'); fPatternFT = 'FTdata'; end
     if ~exist('nptsApprox', 'var'); nptsApprox = 12*365.25*3*24; end
     if ~exist('magPhase_deg', 'var'); magPhase_deg = 0; end
@@ -212,7 +215,7 @@
                     nSW_pcc = 4 / a_AU^2 * ones(1,npts);
                     vSW_kms = 400  * ones(1,npts);
                     [mpBvec, OUTSIDE_MP] = MpauseFld(nSW_pcc, vSW_kms, t_h*3600, xyz_km, ...
-                        Mdip_nT, Odip_km, S3coords, parentName, MPmodel, SPHOUT);
+                        Mdip_nT, Odip_km, parentName, S3coords, MPmodel, coeffPath, SPHOUT);
                     Bvec = Bvec + mpBvec;
                     Bvec(:,OUTSIDE_MP) = 0;
 

comparison/CompareNepModels.m

@@ -177,15 +177,5 @@ function CompareNepModels(LIVE_PLOTS, scName, SEQUENTIAL)
     name{3} = "SPICE in NLS as defined in O8";
 
     legend([scTraj{1} scTraj{2} scTraj{3}], [name{1}, name{2}, name{3}])
-
-    function [x, y, z] = GetDespun(xyz_Rp, despin)
-        r_Rp = sqrt(xyz_Rp(1,:).^2 + xyz_Rp(2,:).^2 + xyz_Rp(3,:).^2);
-        theta = acos(xyz_Rp(3,:) ./ r_Rp);
-        phi = atan2(xyz_Rp(2,:), xyz_Rp(1,:));
-
-        x = r_Rp .* sin(theta) .* cos(phi + despin);
-        y = r_Rp .* sin(theta) .* sin(phi + despin);
-        z = r_Rp .* cos(theta);
-    end
 
-end
+end

comparison/OutputHEALpixTo_n10.m

@@ -61,42 +61,3 @@ function OutputHEALpixTo_n10(hpDir, hpFname, outDir, outFbase)
     end
 
 end
-
-
-function WritePureHarmFile(fpath, header, colNames, outFileBase, n, m, gnm, hnm, r_RP, ...
-    theta_rad, phi_rad, outFlag)
-% Print a formatted file containing the specified pure harmonic.
-%
-% Parameters
-% ----------
-% fpath : char, 1xC
-%   Full path to output file.
-% header : char, 1xD
-%   Description of the output file contents.
-% colNames : char, 1xE
-%   Formatted string describing output columns.
-% outFileBase : char, 1xF
-%   Pattern for output file names to which n and m are appended.
-% n, m : int
-%   Degree and order of output functions in Schmidt normalization.
-% gnm, hnm : double
-%   Value to use for Schmidt semi-normalized gauss coefficients in gauss.
-% r_RP : double, 1xN
-%   Radius in units of :math:`R_P` (planetary radii) at which to evaluate functions.
-% theta_rad : double, 1xN
-%   Colatitude in radians at which to evaluate functions.
-% phi_rad : double, 1xN
-%   East longitude in radians at which to evaluate functions.
-% outFlag : char, 1xG
-%   Strings to append to filenames to distinguish pure g from pure h for each n,m pair.
-
-    coeffsLine = sprintf('%12d,%12d,%12.8f,%12.8f', n, m, gnm, hnm);
-    dlmwrite(fpath, header, 'delimiter', '');
-    dlmwrite(fpath, colNames, 'delimiter','', '-append');
-    dlmwrite(fpath, coeffsLine, 'delimiter','', 'precision',8, '-append');
-
-    [Bvec, ~, ~] = MagFldParent('PureHarmonic', r_RP, theta_rad, phi_rad, 'SphericalHarmonic', ...
-                                'None', 0, 0);
-    save([outFileBase sprintf('n%02dm%02d%s.mat', n, m, outFlag)], 'r_RP', 'theta_rad', ...
-         'phi_rad', 'Bvec');
-end

comparison/WritePureHarmFile.m

@@ -0,0 +1,43 @@
+function WritePureHarmFile(fpath, header, colNames, outFileBase, n, m, gnm, hnm, r_RP, ...
+    theta_rad, phi_rad, outFlag)
+% Print a formatted file containing the specified pure harmonic.
+%
+% Parameters
+% ----------
+% fpath : char, 1xC
+%   Full path to output file.
+% header : char, 1xD
+%   Description of the output file contents.
+% colNames : char, 1xE
+%   Formatted string describing output columns.
+% outFileBase : char, 1xF
+%   Pattern for output file names to which n and m are appended.
+% n, m : int
+%   Degree and order of output functions in Schmidt normalization.
+% gnm, hnm : double
+%   Value to use for Schmidt semi-normalized gauss coefficients in gauss.
+% r_RP : double, 1xN
+%   Radius in units of :math:`R_P` (planetary radii) at which to evaluate functions.
+% theta_rad : double, 1xN
+%   Colatitude in radians at which to evaluate functions.
+% phi_rad : double, 1xN
+%   East longitude in radians at which to evaluate functions.
+% outFlag : char, 1xG
+%   Strings to append to filenames to distinguish pure g from pure h for each n,m pair.
+
+% Part of the PlanetMag framework for evaluation and study of planetary magnetic fields.
+% Created by Corey J. Cochrane and Marshall J. Styczinski
+% Maintained by Marshall J. Styczinski
+% Contact: corey.j.cochrane@jpl.nasa.gov
+%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+
+    coeffsLine = sprintf('%12d,%12d,%12.8f,%12.8f', n, m, gnm, hnm);
+    dlmwrite(fpath, header, 'delimiter', '');
+    dlmwrite(fpath, colNames, 'delimiter','', '-append');
+    dlmwrite(fpath, coeffsLine, 'delimiter','', 'precision',8, '-append');
+
+    [Bvec, ~, ~] = MagFldParent('PureHarmonic', r_RP, theta_rad, phi_rad, 'SphericalHarmonic', ...
+                                'None', 0, 0);
+    save([outFileBase sprintf('n%02dm%02d%s.mat', n, m, outFlag)], 'r_RP', 'theta_rad', ...
+         'phi_rad', 'Bvec');
+end

docs/functions.rst

@@ -6,8 +6,8 @@
 Matlab supporting functions
 ===========================
 
-Calculation functions
----------------------
+Main calculation functions
+--------------------------
 
 .. mat:module:: functions
    :noindex:
@@ -48,9 +48,56 @@ Calculation functions
 .. mat:autofunction:: dLegendreS
 
 
+Magnetopause functions
+----------------------
+
+.. mat:module:: functions.magnetopauseModels
+   :noindex:
+
+:mat:func:`GetMPsurfAB2005`
++++++++++++++++++++++++++++
+
+.. mat:autofunction:: GetMPsurfAB2005
+
+:mat:func:`GetMPsurfS1997`
+++++++++++++++++++++++++++
+
+.. mat:autofunction:: GetMPsurfS1997
+
+:mat:func:`GetMPsurfSM1996`
++++++++++++++++++++++++++++
+
+.. mat:autofunction:: GetMPsurfSM1996
+
+:mat:func:`MPboxHarmonic`
++++++++++++++++++++++++++
+
+.. mat:autofunction:: MPboxHarmonic
+
+:mat:func:`MPsphericalHarmonic`
++++++++++++++++++++++++++++++++
+
+.. mat:autofunction:: MPsphericalHarmonic
+
+
+Khurana and Schwarzl (2005) Jupiter model functions
+---------------------------------------------------
+
+.. mat:module:: functions.KS2005functions
+   :noindex:
+
+:mat:func:`VIP4noDipole`
+++++++++++++++++++++++++
+
+.. mat:autofunction:: VIP4noDipole
+
+
 Additional supporting functions
 -------------------------------
 
+.. mat:module:: functions
+   :noindex:
+
 :mat:func:`Bsph2Bxyz`
 +++++++++++++++++++++
 
@@ -91,6 +138,11 @@ Additional supporting functions
 
 .. mat:autofunction:: GetBplotAndLsqNeptune
 
+:mat:func:`GetDespun`
++++++++++++++++++++++
+
+.. mat:autofunction:: GetDespun
+
 :mat:func:`GetExcitations`
 ++++++++++++++++++++++++++
 
@@ -141,6 +193,11 @@ Additional supporting functions
 
 .. mat:autofunction:: RotateVecSpice
 
+:mat:func:`vlines`
+++++++++++++++++++
+
+.. mat:autofunction:: vlines
+
 Model coefficient functions
 ---------------------------
 

functions/ExcitationSpectrum.m

@@ -1,5 +1,5 @@
 function [Tpeak_h, B0vec, B1vec] = ExcitationSpectrum(moonName, nOsc, rate, Tinterest_h, ...
-    outData, fPatternFT, fPatternTseries, magPhase, DO_EPO, SPHOUT, DO_MPAUSE)
+    outData, coeffPath, fPatternFT, fPatternTseries, magPhase, DO_EPO, SPHOUT, DO_MPAUSE)
 % Calculate a Fast Fourier Transform (FFT) spectrum of magnetic oscillations for a given moon.
 %
 % Note
@@ -26,6 +26,8 @@
 %   synodic period, as evaluated from values returned by GetBodyParams.
 % outData : char, 1xF, default='out'
 %   Directory to use for output of complex spectrum amplitudes.
+% coeffPath : char, 1xE, default='modelCoeffs'
+%   Directory containing model coefficients files.
 % fPatternFT : char, 1xG, default='FTdata'
 %   Pattern for file names of FFT spectrum data saved to disk.
 % fPatternTseries : char, 1xH, default='TseriesData'
@@ -59,6 +61,7 @@
     if ~exist('nOsc', 'var'); nOsc = 5000; end
     if ~exist('rate', 'var'); rate = 100; end
     if ~exist('outData', 'var'); outData = 'out'; end
+    if ~exist('coeffPath', 'var'); coeffPath = 'modelCoeffs'; end
     if ~exist('fPatternFT', 'var'); fPatternFT = 'FTdata'; end
     if ~exist('fPatternTseries', 'var'); fPatternTseries = 'TseriesData'; end
     if ~exist('magPhase', 'var'); magPhase = 0; end
@@ -118,8 +121,8 @@
 
         nSW_pcc = 0.14 * ones(1,npts);
         vSW_kms = 400  * ones(1,npts);
-        [mpBvec, OUTSIDE_MP] = MpauseFld(nSW_pcc, vSW_kms, t_h*3600, xyz_km, ...
-            Mdip_nT, Odip_km, parentName, MPmodel, SPHOUT);
+        [mpBvec, OUTSIDE_MP] = MpauseFld(nSW_pcc, vSW_kms, t_h*3600, xyz_km, Mdip_nT, Odip_km, ...
+            parentName, S3coords, MPmodel, coeffPath, MPmodel, SPHOUT);
         Bvec = Bvec + mpBvec;
         Bvec(:,OUTSIDE_MP) = 0;
 

functions/GetBplotAndLsq.m

@@ -1,5 +1,5 @@
 function GetBplotAndLsq(ets, t_h, r_km, theta, phi, xyz_km, BrSC, BthSC, BphiSC, scName, ...
-    parentName, S3coords, orbStr, opt, MPopt, SEQUENTIAL, jt_h)
+    parentName, S3coords, orbStr, opt, MPopt, SEQUENTIAL, coeffPath, jt_h)
 % Plots and calculates comparisons between modeled and measured magnetic fields.
 %
 % Generates time series data of a specified combination of magnetic field models implemented in
@@ -45,6 +45,8 @@ function GetBplotAndLsq(ets, t_h, r_km, theta, phi, xyz_km, BrSC, BthSC, BphiSC,
 %   available options.
 % SEQUENTIAL : bool, default=0
 %   Whether to plot points by index or hours relative to a reference time.
+% coeffPath : char, 1xE, default='modelCoeffs'
+%   Directory containing model coefficients files.
 % jt_h : double, 1xM, default=[]
 %   Ephemeris times of Juno measurements in TDB hours relative to J2000 to use in data comparisons.
 
@@ -55,6 +57,7 @@ function GetBplotAndLsq(ets, t_h, r_km, theta, phi, xyz_km, BrSC, BthSC, BphiSC,
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
     if ~exist('SEQUENTIAL', 'var'); SEQUENTIAL = 0; end
+    if ~exist('coeffpath', 'var'); coeffPath = 'modelCoeffs'; end
     if ~exist('jt_h', 'var'); jt_h = []; end
 
     npts = length(t_h);
@@ -79,8 +82,8 @@ function GetBplotAndLsq(ets, t_h, r_km, theta, phi, xyz_km, BrSC, BthSC, BphiSC,
 
         nSW_pcc = 0.14 * ones(1,npts);
         vSW_kms = 400  * ones(1,npts);
-        [mpBvec, OUTSIDE_MP] = MpauseFld(nSW_pcc, vSW_kms, t_h*3600, xyz_km, ...
-            Mdip_nT, Odip_km, S3coords, parentName, MPmodel, 1);
+        [mpBvec, OUTSIDE_MP] = MpauseFld(nSW_pcc, vSW_kms, t_h*3600, xyz_km, Mdip_nT, Odip_km, ...
+            parentName, S3coords, MPmodel, coeffPath, 1);
         Bvec = Bvec + mpBvec;
         Bvec(:,OUTSIDE_MP) = 0;
 
@@ -152,4 +155,4 @@ function GetBplotAndLsq(ets, t_h, r_km, theta, phi, xyz_km, BrSC, BthSC, BphiSC,
     BphiLsq = BphiD.^2;
     chi2 = sum([BrLsq, BthLsq, BphiLsq], 'all') / 3 / npts;
     disp(['Overall, for this range and model chi^2 = ' sprintf('%.2f', chi2) '.'])
-end
+end

functions/GetBplotAndLsqMoon.m

@@ -1,5 +1,5 @@
 function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC, BzSC, scName, ...
-    parentName, S3coords, moonName, era, fbStr, opt, MPopt, SEQUENTIAL, dataDir, jt_h)
+    parentName, S3coords, moonName, era, fbStr, opt, MPopt, SEQUENTIAL, dataDir, coeffPath, jt_h)
 % Plots and calculates comparisons between modeled and measured magnetic fields in the vicinity of
 % a target moon.
 %
@@ -59,6 +59,8 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
 %   Whether to plot points by index or hours relative to a reference time (closest approach).
 % dataDir : char, 1xH, default='out'
 %   Name of directory where excitation moments are printed to disk.
+% coeffPath : char, 1xE, default='modelCoeffs'
+%   Directory containing model coefficients files.
 % jt_h : double, 1xM, default=[]
 %   Ephemeris times of Juno measurements in TDB hours relative to J2000 to use in data comparisons.
 
@@ -69,6 +71,7 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
 %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 
     if ~exist('dataDir', 'var'); dataDir = 'out'; end
+    if ~exist('coeffpath', 'var'); coeffPath = 'modelCoeffs'; end
     if ~exist('jt_h', 'var'); jt_h = []; end
     if isempty(jt_h); JUNOTOO=0; else; JUNOTOO=1; end
 
@@ -87,8 +90,8 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
 
         nSW_pcc = 0.14 * ones(1,npts);
         vSW_kms = 400  * ones(1,npts);
-        [mpBvec, OUTSIDE_MP] = MpauseFld(nSW_pcc, vSW_kms, [t_h, jt_h]*3600, xyz_km, ...
-            Mdip_nT, Odip_km, S3coords, parentName, MPmodel, 1);
+        [mpBvec, OUTSIDE_MP] = MpauseFld(nSW_pcc, vSW_kms, [t_h, jt_h]*3600, xyz_km, Mdip_nT, ...
+            Odip_km, parentName, S3coords, MPmodel, coeffPath, 1);
         Bvec = Bvec + mpBvec;
         Bvec(:,OUTSIDE_MP) = 0;
 
@@ -197,7 +200,7 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
     plot(xx, Bx);
     plot(xx, BxSC);
     plot(xx, BxExc);
-    if JUNOTOO; jvlines(xJuno, b000ff); end
+    if JUNOTOO; vlines(xJuno, b000ff); end
     xlabel(xDescrip);
     ylabel('Vector component (nT)');
     legend(magModelDescrip, strcat(scName, ' MAG'), excStr);
@@ -206,7 +209,7 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
     plot(xx, By);
     plot(xx, BySC);
     plot(xx, ByExc);
-    if JUNOTOO; jvlines(xJuno, b000ff); end
+    if JUNOTOO; vlines(xJuno, b000ff); end
     xlabel(xDescrip);
     ylabel('Vector component (nT)');
     legend(magModelDescrip, strcat(scName, ' MAG'), excStr);
@@ -215,7 +218,7 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
     plot(xx, Bz);
     plot(xx, BzSC);
     plot(xx, BzExc);
-    if JUNOTOO; jvlines(xJuno, b000ff); end
+    if JUNOTOO; vlines(xJuno, b000ff); end
     xlabel(xDescrip);
     ylabel('Vector component (nT)');
     legend(magModelDescrip, strcat(scName, ' MAG'), excStr);
@@ -235,7 +238,7 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
     plot(xx, BxDexc);
     plot(xx, ByDexc);
     plot(xx, BzDexc);
-    if JUNOTOO; jvlines(xJuno, b000ff); end
+    if JUNOTOO; vlines(xJuno, b000ff); end
     xlabel(xDescrip);
     ylabel('Component difference (nT)');
     legend('\Delta B_x', '\Delta B_y', '\Delta B_z', '\Delta B^e_x', '\Delta B^e_y', ...
@@ -257,10 +260,3 @@ function GetBplotAndLsqMoon(ets, t_h, r_km, theta, phi, xyz_km, r_RM, BxSC, BySC
     disp(['Excitation moments only, beyond ' num2str(RMmin) ...
         ' R_' moonName(1) ': chi^2 = ' sprintf('%.2f', chi2exc) '.'])
 end
-
-function jvlines(x, color)
-    yLim = get(gca,'ylim');
-    for i=1:length(x)
-        plot([x(i), x(i)], yLim, 'Color', color);
-    end
-end