From 2adc1d312260414c1421767152d01fa4aa75dc09 Mon Sep 17 00:00:00 2001 From: yumouwei <46117079+yumouwei@users.noreply.github.com> Date: Thu, 18 Dec 2025 11:53:39 -0500 Subject: [PATCH 1/6] Revamp d3d.get_n1_bradial_parameters method --- disruption_py/machine/d3d/physics.py | 66 ++++++++++++++++++++++++++++ 1 file changed, 66 insertions(+) diff --git a/disruption_py/machine/d3d/physics.py b/disruption_py/machine/d3d/physics.py index 2149ee69e..2e1f6540c 100644 --- a/disruption_py/machine/d3d/physics.py +++ b/disruption_py/machine/d3d/physics.py @@ -898,6 +898,72 @@ def get_z_parameters(params: PhysicsMethodParams): z_cur_norm = z_cur / nominal_flattop_radius return {"zcur": z_cur, "zcur_normalized": z_cur_norm} + @staticmethod + @physics_method( + columns=["n_equal_1_normalized", "n_equal_1_mode"], + tokamak=Tokamak.D3D, + ) + def get_n1_bradial_parameters(params: PhysicsMethodParams): + ''' + TODO: finish docstring + Docstring for get_n1_bradial_parameters + + :param params: Description + :type params: PhysicsMethodParams + + References + https://github.com/MIT-PSFC/disruption-py/issues/261 + https://github.com/MIT-PSFC/disruption-py/blob/matlab/DIII-D/get_n1_bradial_d3d.m + + This method was previously removed in v0.8 + ''' + # The following shots are missing bradial calculations in MDSplus and + # must be loaded from a separate datafile + # TODO: re-check these shot range + if 156199 <= params.shot_id <= 172430: + # TODO: load from tree files on disc + raise NotImplementedError + if 176030 <= params.shot_id <= 176912: + # TODO: load from netcdf file on disc + raise NotImplementedError + + try: + # Get data from the ONFR system + n_equal_1_mode, t_n1 = params.mds_conn.get_data_with_dims( + f"ptdata('onsbradial', {params.shot_id})", + ) # [G], [ms] # TODO: check unit + except mdsExceptions.MdsException: + try: + # Fallback: get data from the legacy DUD system + n_equal_1_mode, t_n1 = params.mds_conn.get_data_with_dims( + f"ptdata('dusbradial', {params.shot_id})", + ) # [G], [ms] # TODO: check unit + params.logger.verbose('n_equal_1_mode: Failed to get ONSBRADIAL signal. Use DUSBRADIAL instead.') + except mdsExceptions.MdsException: + params.logger.debug("n_equal_1_mode: Failed to get either ONSBRADIAL or DUSBRADIAL signal") + return {"n_equal_1_normalized": [np.nan], "n_equal_1_mode": [np.nan]} + t_n1 /= 1e3 # [ms] -> [s] + n_equal_1_mode = interp1(t_n1, n_equal_1_mode, params.times) + + # Calculate n_equal_1_normalized + # TODO: move b_tor calculation to an individual method? + # TODO: shouldn't we interpolate b_tor to n1_mode timebase, compute _norm, then interpolate to output timebase? + try: + b_tor, t_b_tor = params.mds_conn.get_data_with_dims( + f"ptdata('bt', {params.shot_id})" + ) # [ms], [T?] # TODO: check unit -- this should be identical to n1rms_norm computation + t_b_tor /= 1e3 # [ms] -> [s] + b_tor = interp1(t_b_tor, b_tor, params.times) # [T] + n_equal_1_normalized = n_equal_1_mode / np.abs(b_tor) + except mdsExceptions.MdsException as e: + params.logger.warning( + "Failed to get b_tor signal to compute n_equal_1_normalized" + ) + params.logger.opt(exception=True).debug(e) + n_equal_1_normalized = [np.nan] + # TODO: return btor? Yes in CMOD and EAST + return {"n_equal_1_mode": n_equal_1_mode, "n_equal_1_normalized": n_equal_1_normalized} + @staticmethod @physics_method(columns=["n1rms", "n1rms_normalized"], tokamak=Tokamak.D3D) def get_n1rms_parameters(params: PhysicsMethodParams): From c25e4fc87b13368a6235ed88dcbb4d8be71352bb Mon Sep 17 00:00:00 2001 From: yumouwei <46117079+yumouwei@users.noreply.github.com> Date: Thu, 18 Dec 2025 14:05:18 -0800 Subject: [PATCH 2/6] Test method on d3d server, then made improvements --- disruption_py/machine/d3d/physics.py | 39 ++++++++++++++++------------ 1 file changed, 23 insertions(+), 16 deletions(-) diff --git a/disruption_py/machine/d3d/physics.py b/disruption_py/machine/d3d/physics.py index 2e1f6540c..458645185 100644 --- a/disruption_py/machine/d3d/physics.py +++ b/disruption_py/machine/d3d/physics.py @@ -900,23 +900,23 @@ def get_z_parameters(params: PhysicsMethodParams): @staticmethod @physics_method( - columns=["n_equal_1_normalized", "n_equal_1_mode"], + columns=["n_equal_1_normalized", "n_equal_1_mode", "btor"], tokamak=Tokamak.D3D, ) def get_n1_bradial_parameters(params: PhysicsMethodParams): - ''' + """ TODO: finish docstring Docstring for get_n1_bradial_parameters - + :param params: Description :type params: PhysicsMethodParams - + References https://github.com/MIT-PSFC/disruption-py/issues/261 https://github.com/MIT-PSFC/disruption-py/blob/matlab/DIII-D/get_n1_bradial_d3d.m - + This method was previously removed in v0.8 - ''' + """ # The following shots are missing bradial calculations in MDSplus and # must be loaded from a separate datafile # TODO: re-check these shot range @@ -926,32 +926,36 @@ def get_n1_bradial_parameters(params: PhysicsMethodParams): if 176030 <= params.shot_id <= 176912: # TODO: load from netcdf file on disc raise NotImplementedError - + try: # Get data from the ONFR system n_equal_1_mode, t_n1 = params.mds_conn.get_data_with_dims( f"ptdata('onsbradial', {params.shot_id})", - ) # [G], [ms] # TODO: check unit + ) # [G], [ms] except mdsExceptions.MdsException: try: # Fallback: get data from the legacy DUD system n_equal_1_mode, t_n1 = params.mds_conn.get_data_with_dims( f"ptdata('dusbradial', {params.shot_id})", - ) # [G], [ms] # TODO: check unit - params.logger.verbose('n_equal_1_mode: Failed to get ONSBRADIAL signal. Use DUSBRADIAL instead.') + ) # [G], [ms] + params.logger.verbose( + "n_equal_1_mode: Failed to get ONSBRADIAL signal. Use DUSBRADIAL instead." + ) except mdsExceptions.MdsException: - params.logger.debug("n_equal_1_mode: Failed to get either ONSBRADIAL or DUSBRADIAL signal") + params.logger.warning( + "n_equal_1_mode: Failed to get either ONSBRADIAL or DUSBRADIAL signal" + ) return {"n_equal_1_normalized": [np.nan], "n_equal_1_mode": [np.nan]} - t_n1 /= 1e3 # [ms] -> [s] + t_n1 /= 1e3 # [ms] -> [s] n_equal_1_mode = interp1(t_n1, n_equal_1_mode, params.times) - + # Calculate n_equal_1_normalized # TODO: move b_tor calculation to an individual method? # TODO: shouldn't we interpolate b_tor to n1_mode timebase, compute _norm, then interpolate to output timebase? try: b_tor, t_b_tor = params.mds_conn.get_data_with_dims( f"ptdata('bt', {params.shot_id})" - ) # [ms], [T?] # TODO: check unit -- this should be identical to n1rms_norm computation + ) # [ms], [T] t_b_tor /= 1e3 # [ms] -> [s] b_tor = interp1(t_b_tor, b_tor, params.times) # [T] n_equal_1_normalized = n_equal_1_mode / np.abs(b_tor) @@ -961,8 +965,11 @@ def get_n1_bradial_parameters(params: PhysicsMethodParams): ) params.logger.opt(exception=True).debug(e) n_equal_1_normalized = [np.nan] - # TODO: return btor? Yes in CMOD and EAST - return {"n_equal_1_mode": n_equal_1_mode, "n_equal_1_normalized": n_equal_1_normalized} + return { + "n_equal_1_mode": n_equal_1_mode, + "n_equal_1_normalized": n_equal_1_normalized, + "btor": b_tor, + } @staticmethod @physics_method(columns=["n1rms", "n1rms_normalized"], tokamak=Tokamak.D3D) From ba7b1e797afceaffb3b2dab924a149740e9d3d81 Mon Sep 17 00:00:00 2001 From: yumouwei <46117079+yumouwei@users.noreply.github.com> Date: Thu, 18 Dec 2025 14:45:54 -0800 Subject: [PATCH 3/6] Move get_btor to its own physics method --- disruption_py/machine/d3d/physics.py | 66 ++++++++++++++++++---------- 1 file changed, 42 insertions(+), 24 deletions(-) diff --git a/disruption_py/machine/d3d/physics.py b/disruption_py/machine/d3d/physics.py index 458645185..9f8fb839a 100644 --- a/disruption_py/machine/d3d/physics.py +++ b/disruption_py/machine/d3d/physics.py @@ -898,9 +898,42 @@ def get_z_parameters(params: PhysicsMethodParams): z_cur_norm = z_cur / nominal_flattop_radius return {"zcur": z_cur, "zcur_normalized": z_cur_norm} + @staticmethod + @cache_method + @physics_method(columns=["btor"], tokamak=Tokamak.D3D) + def get_btor(params: PhysicsMethodParams): + """ + Calculate the toroidal magnetic field. + + Parameters + ---------- + params : PhysicsMethodParams + Parameters containing MDS connection and shot information + + Returns + ------- + dict + A dictionary containing the toroidal magnetic field (`btor`). + + References + ------- + # TODO: finish references + """ + try: + b_tor, t_b_tor = params.mds_conn.get_data_with_dims( + f"ptdata('bt', {params.shot_id})" + ) # [ms], [T] + except mdsExceptions.MdsException as e: + params.logger.warning("Failed to get b_tor signal") + params.logger.opt(exception=True).debug(e) + return {"btor": [np.nan]} + t_b_tor /= 1e3 # [ms] -> [s] + b_tor = interp1(t_b_tor, b_tor, params.times) + return {"btor": b_tor} + @staticmethod @physics_method( - columns=["n_equal_1_normalized", "n_equal_1_mode", "btor"], + columns=["n_equal_1_normalized", "n_equal_1_mode"], tokamak=Tokamak.D3D, ) def get_n1_bradial_parameters(params: PhysicsMethodParams): @@ -950,25 +983,16 @@ def get_n1_bradial_parameters(params: PhysicsMethodParams): n_equal_1_mode = interp1(t_n1, n_equal_1_mode, params.times) # Calculate n_equal_1_normalized - # TODO: move b_tor calculation to an individual method? - # TODO: shouldn't we interpolate b_tor to n1_mode timebase, compute _norm, then interpolate to output timebase? - try: - b_tor, t_b_tor = params.mds_conn.get_data_with_dims( - f"ptdata('bt', {params.shot_id})" - ) # [ms], [T] - t_b_tor /= 1e3 # [ms] -> [s] - b_tor = interp1(t_b_tor, b_tor, params.times) # [T] - n_equal_1_normalized = n_equal_1_mode / np.abs(b_tor) - except mdsExceptions.MdsException as e: + b_tor = D3DPhysicsMethods.get_btor(params)["btor"] + if np.isnan(b_tor).all(): params.logger.warning( "Failed to get b_tor signal to compute n_equal_1_normalized" ) - params.logger.opt(exception=True).debug(e) n_equal_1_normalized = [np.nan] + n_equal_1_normalized = n_equal_1_mode / np.abs(b_tor) return { "n_equal_1_mode": n_equal_1_mode, "n_equal_1_normalized": n_equal_1_normalized, - "btor": b_tor, } @staticmethod @@ -998,21 +1022,15 @@ def get_n1rms_parameters(params: PhysicsMethodParams): n1rms, t_n1rms = params.mds_conn.get_data_with_dims(r"\n1rms", tree_name="d3d") n1rms *= 1.0e-4 # Gauss -> Tesla t_n1rms /= 1e3 # [ms] -> [s] - n1rms = interp1(t_n1rms, n1rms, params.times) + n1rms = interp1(t_n1rms, n1rms, params.times) # Calculate n1rms_norm - try: - b_tor, t_b_tor = params.mds_conn.get_data_with_dims( - f"ptdata('bt', {params.shot_id})" - ) - t_b_tor /= 1e3 # [ms] -> [s] - b_tor = interp1(t_b_tor, b_tor, params.times) # [T] - n1rms_norm = n1rms / np.abs(b_tor) - except mdsExceptions.MdsException as e: + b_tor = D3DPhysicsMethods.get_btor(params)["btor"] + if np.isnan(b_tor).all(): params.logger.warning( - "Failed to get b_tor signal to compute n1rms_normalized" + "Failed to get b_tor signal to compute n_equal_1_normalized" ) - params.logger.opt(exception=True).debug(e) n1rms_norm = [np.nan] + n1rms_norm = n1rms / np.abs(b_tor) return {"n1rms": n1rms, "n1rms_normalized": n1rms_norm} # TODO: Need to test and unblock recalculating peaking factors From 562e27e76fba518f9a560dc11e1ae5d49805a678 Mon Sep 17 00:00:00 2001 From: yumouwei <46117079+yumouwei@users.noreply.github.com> Date: Thu, 18 Dec 2025 14:46:13 -0800 Subject: [PATCH 4/6] Fix n_equal_1_mode unit --- disruption_py/machine/d3d/physics.py | 1 + 1 file changed, 1 insertion(+) diff --git a/disruption_py/machine/d3d/physics.py b/disruption_py/machine/d3d/physics.py index 9f8fb839a..047def789 100644 --- a/disruption_py/machine/d3d/physics.py +++ b/disruption_py/machine/d3d/physics.py @@ -980,6 +980,7 @@ def get_n1_bradial_parameters(params: PhysicsMethodParams): ) return {"n_equal_1_normalized": [np.nan], "n_equal_1_mode": [np.nan]} t_n1 /= 1e3 # [ms] -> [s] + n_equal_1_mode *= 1.0e-4 # [G] -> [T] n_equal_1_mode = interp1(t_n1, n_equal_1_mode, params.times) # Calculate n_equal_1_normalized From fbd4b24866946f2b7282ca79b9e529742f6d1f0e Mon Sep 17 00:00:00 2001 From: yumouwei <46117079+yumouwei@users.noreply.github.com> Date: Thu, 18 Dec 2025 14:56:43 -0800 Subject: [PATCH 5/6] Add docstrings --- disruption_py/machine/d3d/physics.py | 36 ++++++++++++++++++---------- 1 file changed, 23 insertions(+), 13 deletions(-) diff --git a/disruption_py/machine/d3d/physics.py b/disruption_py/machine/d3d/physics.py index 047def789..dcc36b49f 100644 --- a/disruption_py/machine/d3d/physics.py +++ b/disruption_py/machine/d3d/physics.py @@ -917,7 +917,9 @@ def get_btor(params: PhysicsMethodParams): References ------- - # TODO: finish references + - original sources: [get_n1_bradial_d3d.m](https://github.com/MIT-PSFC/disruption + -py/blob/matlab/DIII-D/get_n1_bradial_d3d.m), [get_n1rms_d3d.m](https://github.com + /MIT-PSFC/disruption-py/blob/matlab/DIII-D/get_n1rms_d3d.m) """ try: b_tor, t_b_tor = params.mds_conn.get_data_with_dims( @@ -938,26 +940,34 @@ def get_btor(params: PhysicsMethodParams): ) def get_n1_bradial_parameters(params: PhysicsMethodParams): """ - TODO: finish docstring - Docstring for get_n1_bradial_parameters + This method obtaines the n=1 Bradial information (`n_equal_1_mode`) from + the ESLD coils (units of Tesla). It also calculates the normalized Bradial, + dividing by the toroidal B-field (`n_equal_1_normalized`). - :param params: Description - :type params: PhysicsMethodParams + Parameters + ---------- + params : PhysicsMethodParams + The parameters containing the MDSplus connection, shot id and more. - References - https://github.com/MIT-PSFC/disruption-py/issues/261 - https://github.com/MIT-PSFC/disruption-py/blob/matlab/DIII-D/get_n1_bradial_d3d.m + Returns + ------- + dict + A dictionary containing `n_equal_1_mode` and `n_equal_1_normalized`. - This method was previously removed in v0.8 + References + ------- + - original source: [get_n1_bradial_d3d.m](https://github.com/MIT-PSFC/disruption-py + /blob/matlab/DIII-D/get_n1_bradial_d3d.m) + - pull requests: + - issues: #[261](https://github.com/MIT-PSFC/disruption-py/issues/261) """ # The following shots are missing bradial calculations in MDSplus and # must be loaded from a separate datafile - # TODO: re-check these shot range if 156199 <= params.shot_id <= 172430: - # TODO: load from tree files on disc + # TODO: load data from custom tree structures raise NotImplementedError if 176030 <= params.shot_id <= 176912: - # TODO: load from netcdf file on disc + # TODO: load data from NetCDF file raise NotImplementedError try: @@ -1023,7 +1033,7 @@ def get_n1rms_parameters(params: PhysicsMethodParams): n1rms, t_n1rms = params.mds_conn.get_data_with_dims(r"\n1rms", tree_name="d3d") n1rms *= 1.0e-4 # Gauss -> Tesla t_n1rms /= 1e3 # [ms] -> [s] - n1rms = interp1(t_n1rms, n1rms, params.times) + n1rms = interp1(t_n1rms, n1rms, params.times) # Calculate n1rms_norm b_tor = D3DPhysicsMethods.get_btor(params)["btor"] if np.isnan(b_tor).all(): From 19f6fd6dd49b341ea1c3ab89d8916b99c38dc552 Mon Sep 17 00:00:00 2001 From: yumouwei <46117079+yumouwei@users.noreply.github.com> Date: Thu, 18 Dec 2025 15:37:47 -0800 Subject: [PATCH 6/6] Add PR reference --- disruption_py/machine/d3d/physics.py | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/disruption_py/machine/d3d/physics.py b/disruption_py/machine/d3d/physics.py index dcc36b49f..1340ce766 100644 --- a/disruption_py/machine/d3d/physics.py +++ b/disruption_py/machine/d3d/physics.py @@ -958,7 +958,7 @@ def get_n1_bradial_parameters(params: PhysicsMethodParams): ------- - original source: [get_n1_bradial_d3d.m](https://github.com/MIT-PSFC/disruption-py /blob/matlab/DIII-D/get_n1_bradial_d3d.m) - - pull requests: + - pull requests: #[500](https://github.com/MIT-PSFC/disruption-py/pull/500) - issues: #[261](https://github.com/MIT-PSFC/disruption-py/issues/261) """ # The following shots are missing bradial calculations in MDSplus and