diff --git a/matlab/TCV_IMAS/tcv_get_ids_bolometer.m b/matlab/TCV_IMAS/tcv_get_ids_bolometer.m new file mode 100644 index 0000000000000000000000000000000000000000..2d34f4b5e53925dbf645a541f4facad0c5d07525 --- /dev/null +++ b/matlab/TCV_IMAS/tcv_get_ids_bolometer.m @@ -0,0 +1,188 @@ +function [ids_bolometer,ids_bolometer_description] = ... + tcv_get_ids_bolometer(shot, ids_bolometer_empty, gdat_params,varargin) +% +% [ids_bolometer,ids_bolometer_description] = ... +% tcv_get_ids_bolometer(shot, ids_bolometer_empty,varargin); +% +% Get the bolometer diagnostics data +% +% ids_bolometer_empty should at least be the empty bolometer ids structure +% in input +% +% gdat_params: gdat_data.gdat_params to get all params passed from original +% call, in particular error_bar options +% +% Date: 20/09/24 +% Author: L. Simons + +% Initialise structure with default ids_properties +% shot=70777; gdat_params=[]; varargin={}; +% aa=gdat([],'ids','source','bolometer'); ids_bolometer_empty=aa.bolometer; +if exist('gdat_params') + [ids_bolometer, ~] = tcv_ids_headpart(shot, ids_bolometer_empty, ... + 'bolometer','homogeneous_time',1,'gdat_params',gdat_params, ... + varargin{:}); +else + [ids_bolometer, ~] = tcv_ids_headpart(shot,ids_bolometer_empty,'bolometer','homogeneous_time',0,varargin{:}); +end + +ids_bolometer_description = struct(); + +% Load the bolometer geometry +bolo_geom=bolou_load_geometry(); + + +params_eff.data_request='\tcv_shot::top.results.bolo_u.intensity'; +params_eff.trialindx=1; +params_eff.machine='tcv'; +params_eff.doplot=0; +bolo_u_intensity = gdat(shot,params_eff); +params_eff.data_request='\tcv_shot::top.results.bolo_u.prad'; +bolo_u_prad = gdat(shot,params_eff); +params_eff.data_request='\tcv_shot::top.results.bolo_u.prad_core'; +bolo_u_prad_core = gdat(shot,params_eff); +params_eff.data_request='\tcv_shot::top.results.bolo_u.confidence'; +bolo_u_confidence = gdat(shot,params_eff); + +status = ~ischar(bolo_u_intensity.data) & ~ischar(bolo_u_prad_core.data); +imas_version_number=getenv('IMAS_VERSION'); +if status + nchannel = numel(bolo_u_intensity.x); + ids_bolometer.channel(1:nchannel) = ids_bolometer.channel(1); + for ii = 1:nchannel + %% Fill geometry information + % Fill information from bolo_geom=bolou_load_geometry(); + ids_bolometer.channel{ii}.name = bolo_geom.channel{ii}.name; + if strcmp(imas_version_number(1),'3') + ids_bolometer.channel{ii}.identifier = bolo_geom.channel{ii}.identifier; + elseif strcmp(imas_version_number(1),'4') + ids_bolometer.channel{ii}.description = bolo_geom.channel{ii}.identifier; + end + ids_bolometer.channel{ii}.detector.geometry_type= bolo_geom.channel{ii}.detector.geometry_type; + ids_bolometer.channel{ii}.detector.centre.phi= bolo_geom.channel{ii}.detector.centre.phi; + ids_bolometer.channel{ii}.detector.centre.r= bolo_geom.channel{ii}.detector.centre.r; + ids_bolometer.channel{ii}.detector.centre.z= bolo_geom.channel{ii}.detector.centre.z; + ids_bolometer.channel{ii}.detector.x1_unit_vector.x= bolo_geom.channel{ii}.detector.x1_unit_vector.x; + ids_bolometer.channel{ii}.detector.x1_unit_vector.y= bolo_geom.channel{ii}.detector.x1_unit_vector.y; + ids_bolometer.channel{ii}.detector.x1_unit_vector.z= bolo_geom.channel{ii}.detector.x1_unit_vector.z; + ids_bolometer.channel{ii}.detector.x2_unit_vector.x= bolo_geom.channel{ii}.detector.x2_unit_vector.x; + ids_bolometer.channel{ii}.detector.x2_unit_vector.y= bolo_geom.channel{ii}.detector.x2_unit_vector.y; + ids_bolometer.channel{ii}.detector.x2_unit_vector.z= bolo_geom.channel{ii}.detector.x2_unit_vector.z; + ids_bolometer.channel{ii}.detector.x3_unit_vector.x= bolo_geom.channel{ii}.detector.x3_unit_vector.x; + ids_bolometer.channel{ii}.detector.x3_unit_vector.y= bolo_geom.channel{ii}.detector.x3_unit_vector.y; + ids_bolometer.channel{ii}.detector.x3_unit_vector.z= bolo_geom.channel{ii}.detector.x3_unit_vector.z; + ids_bolometer.channel{ii}.detector.x1_width= bolo_geom.channel{ii}.detector.x1_width; + ids_bolometer.channel{ii}.detector.x2_width= bolo_geom.channel{ii}.detector.x2_width; + ids_bolometer.channel{ii}.detector.surface= bolo_geom.channel{ii}.detector.surface; + + ids_bolometer.channel{ii}.aperture{1}.geometry_type= bolo_geom.channel{ii}.aperture{1}.geometry_type; + ids_bolometer.channel{ii}.aperture{1}.centre.r= bolo_geom.channel{ii}.aperture{1}.centre.r; + ids_bolometer.channel{ii}.aperture{1}.centre.z= bolo_geom.channel{ii}.aperture{1}.centre.z; + ids_bolometer.channel{ii}.aperture{1}.x1_unit_vector.x= bolo_geom.channel{ii}.aperture{1}.x1_unit_vector.x; + ids_bolometer.channel{ii}.aperture{1}.x1_unit_vector.y= bolo_geom.channel{ii}.aperture{1}.x1_unit_vector.y; + ids_bolometer.channel{ii}.aperture{1}.x1_unit_vector.z= bolo_geom.channel{ii}.aperture{1}.x1_unit_vector.z; + ids_bolometer.channel{ii}.aperture{1}.x2_unit_vector.x= bolo_geom.channel{ii}.aperture{1}.x2_unit_vector.x; + ids_bolometer.channel{ii}.aperture{1}.x2_unit_vector.y= bolo_geom.channel{ii}.aperture{1}.x2_unit_vector.y; + ids_bolometer.channel{ii}.aperture{1}.x2_unit_vector.z= bolo_geom.channel{ii}.aperture{1}.x2_unit_vector.z; + ids_bolometer.channel{ii}.aperture{1}.x3_unit_vector.x= bolo_geom.channel{ii}.aperture{1}.x3_unit_vector.x; + ids_bolometer.channel{ii}.aperture{1}.x3_unit_vector.y= bolo_geom.channel{ii}.aperture{1}.x3_unit_vector.y; + ids_bolometer.channel{ii}.aperture{1}.x3_unit_vector.z= bolo_geom.channel{ii}.aperture{1}.x3_unit_vector.z; + ids_bolometer.channel{ii}.aperture{1}.x1_width= bolo_geom.channel{ii}.aperture{1}.x1_width; + ids_bolometer.channel{ii}.aperture{1}.x2_width= bolo_geom.channel{ii}.aperture{1}.x2_width; + ids_bolometer.channel{ii}.aperture{1}.surface= bolo_geom.channel{ii}.aperture{1}.surface; + + ids_bolometer.channel{ii}.etendue = bolo_geom.channel{ii}.etendue; + ids_bolometer.channel{ii}.etendue_method.name = bolo_geom.channel{ii}.etendue_method.name; + ids_bolometer.channel{ii}.etendue_method.description = bolo_geom.channel{ii}.etendue_method.description; + + ids_bolometer.channel{ii}.line_of_sight.first_point.phi = bolo_geom.channel{ii}.line_of_sight.first_point.phi; + ids_bolometer.channel{ii}.line_of_sight.first_point.r = bolo_geom.channel{ii}.line_of_sight.first_point.r; + ids_bolometer.channel{ii}.line_of_sight.first_point.z = bolo_geom.channel{ii}.line_of_sight.first_point.z; + + ids_bolometer.channel{ii}.line_of_sight.second_point.phi = bolo_geom.channel{ii}.line_of_sight.second_point.phi; + ids_bolometer.channel{ii}.line_of_sight.second_point.r = bolo_geom.channel{ii}.line_of_sight.second_point.r; + ids_bolometer.channel{ii}.line_of_sight.second_point.z = bolo_geom.channel{ii}.line_of_sight.second_point.z; + + %% Fill measurement information + ids_bolometer.channel{ii}.power.data = bolo_u_intensity.data(ii,:); + ids_bolometer_description.channel{ii}.power.data = ... + ['From results.bolo_u.intensity data, Radiance measured by ' ... + 'each bolometer foil in units of [W/m^{2}sr]']; + ids_bolometer.channel{ii}.power.t = bolo_u_intensity.t; + ids_bolometer_description.channel{ii}.power.data = ... + ['From results.bolo_u.intensity dim{1}, time basis of radiance']; + ids_bolometer.channel{ii}.validity = 0; % FIXME: Always valid data + ids_bolometer.channel{ii}.power_radiated_total = bolo_u_prad.data; + ids_bolometer_description.channel{ii}.power_radiated_total = ... + ['From results.bolo_u.prad data, Total radiated power measured by ' ... + 'the bolometer system in units of [kW]']; +% ids_bolometer.channel{ii}.power_radiated_total_error_upper = ... +% 0.*bolo_u_prad.data; % FIXME: Zero error +% ids_bolometer.channel{ii}.power_radiated_total_error_lower = ... +% 0.*bolo_u_prad.data; % FIXME: Zero error + ids_bolometer.channel{ii}.power_radiated_inside_lcfs = ... + bolo_u_prad_core.data; + ids_bolometer_description.channel{ii}.power_radiated_inside_lcfs.data = ... + ['From results.bolo_u.prad_core data, Total radiated power ' ... + 'radiated from inside the last closed flux surface in units of [kW]']; +% ids_bolometer.channel{ii}.power_radiated_inside_lcfs_error_upper = ... +% 0.*bolo_u_prad_core.data; +% ids_bolometer.channel{ii}.power_radiated_inside_lcfs_error_lower = ... +% 0.*bolo_u_prad_core.data; + ids_bolometer.channel{ii}.power_radiated_validity = 0*bolo_u_prad.data; + + ids_bolometer.power_radiated_inside_lcfs = bolo_u_prad_core.data; + ids_bolometer.power_radiated_total = bolo_u_prad.data; + ids_bolometer.power_radiated_validity = 0*bolo_u_prad.data; + end + ids_bolometer.code.name = 'rc_gti_prep'; + ids_bolometer.code.description = ... + ['rc_gti_prep, RADCAM gitlab repo, calls GTI to generate emissivity profiles with liuqe.'] + ids_bolometer.code.commit=''; + ids_bolometer.code.version=''; + ids_bolometer.code.repository= ... + 'https://gitlab.epfl.ch/spc/tcv/diag/radcam/'; + ids_bolometer.code.parameters=['rc_gti_prep(shot,[0 2.2],''bolo'', 10, 1, 0.04, 30,', ... + 'false, false,false, [], ''automated'', 1)']; + + %% Code legacy for TCV_EQ; FBTE + % Call to gti: gti_get_disc + % Call to TCV_eq: temp=tdi('TCV_EQ("psi")',disc.s.equilsrc); + % https://spcwiki.epfl.ch/wiki/Tcv_eq + % a=gdat(shot,'TCV_EQ("psi")'); + ids_bolometer.code.library{1}.name = 'fbte'; + ids_bolometer.code.library{1}.description = 'Magnetic equilibrium'; + ids_bolometer.code.library{1}.commit = ''; + ids_bolometer.code.library{1}.version = ''; + ids_bolometer.code.library{1}.repository = ... + 'https://gitlab.epfl.ch/spc/tcv/tbx/meq'; + ids_bolometer.code.library{1}.parameters = 'TCV_EQ("psi")'; + + + %% Code legacy for GTI + ids_bolometer.code.library{2}.name = 'GTI'; + ids_bolometer.code.library{2}.description = ... + 'General Tomographic Inversion'; + ids_bolometer.code.library{2}.commit = ''; + [a1,a2]=unix('git rev-parse --verify HEAD'); + [b1,b2]=unix('git rev-parse --abbrev-ref HEAD'); + ids_bolometer.code.library{2}.version = sprintf('%s_branch_%s',strtrim(a2),strtrim(b2)); + ids_bolometer.code.library{2}.repository = ... + 'https://gitlab.epfl.ch/spc/tcv/analysis/gti'; + ids_bolometer.code.library{2}.parameters = ''; + + + % Translate bolo_u_confidence values to ids_bolometer.code.output_flag values + if bolo_u_confidence.data == 0 + ids_bolometer.code.output_flag=1; + elseif bolo_u_confidence.data == 1 + ids_bolometer.code.output_flag=0; + else + ids_bolometer.code.output_flag=-abs(bolo_u_confidence.data); + end + ids_bolometer.time = bolo_u_intensity.t; +else + warning('Failed to load data for shot %d',shot); +end + +end diff --git a/matlab/TCV_IMAS/tcv_get_ids_nbi.m b/matlab/TCV_IMAS/tcv_get_ids_nbi.m index eb6777a57379ad39f8a123d920802c65c8fca11e..9d5cf6a56080a327c02766f4cba32082639886f6 100644 --- a/matlab/TCV_IMAS/tcv_get_ids_nbi.m +++ b/matlab/TCV_IMAS/tcv_get_ids_nbi.m @@ -88,7 +88,7 @@ for iunit=1:nb_units params_eff.data_request = ['\results::' results_subname{iunit} ':powr_tcv']; pow=gdat_tcv(shot,params_eff); if ischar(pow.data) - pow.data=0; + pow.data=[]; end ids_nbi.unit{iunit}.power_launched.data = pow.data*1e6; ids_nbi.unit{iunit}.power_launched.time = pow.t; @@ -96,12 +96,18 @@ for iunit=1:nb_units %% energy params_eff.data_request = ['\results::' results_subname{iunit} ':energy']; en=gdat_tcv(shot,params_eff); + if ischar(en.data) + en.data=[]; + end ids_nbi.unit{iunit}.energy.data = en.data*1e3; ids_nbi.unit{iunit}.energy.time = en.t; ids_nbi_description.unit{iunit}.energy = params_eff.data_request; %% power & current fractions params_eff.data_request = ['\results::' results_subname{iunit} ':fraction']; p_frac=gdat(shot,params_eff); + if ischar(p_frac.data) + p_frac.data=[]; + end ids_nbi.unit{iunit}.beam_power_fraction.time = p_frac.t; ids_nbi_description.unit{iunit}.beam_power_fraction = params_eff.data_request; if ~isempty(p_frac.data) && size(p_frac.data,2)>=3 diff --git a/matlab/TCV_IMAS/tcv_get_ids_pf_passive.m b/matlab/TCV_IMAS/tcv_get_ids_pf_passive.m new file mode 100644 index 0000000000000000000000000000000000000000..c95518dce9741b4ebb630f56503004e5df02eeac --- /dev/null +++ b/matlab/TCV_IMAS/tcv_get_ids_pf_passive.m @@ -0,0 +1,85 @@ +function [ids_pf_passive,ids_pf_passive_description,varargout] = tcv_get_ids_pf_passive(shot, ids_pf_passive_empty, gdat_params, varargin) +% +% gdat_params: gdat_data.gdat_params to get all params passed from original call, in particular error_bar options +% + +% Input pharser +if exist('gdat_params','var') + [ids_pf_passive, params_pf_passive] = tcv_ids_headpart(shot, ids_pf_passive_empty,'pf_passive','homogeneous_time',0,'gdat_params',gdat_params,varargin{:}); +else + [ids_pf_passive, params_pf_passive] = tcv_ids_headpart(shot, ids_pf_passive_empty,'pf_passive','homogeneous_time',0,varargin{:}); + aa=gdat_tcv; + gdat_params = aa.gdat_params; % to get default params +end +params_eff_ref = gdat_params; params_eff_ref.doplot=0; +try params_eff_ref=rmfield(params_eff_ref,'source');catch;end % make sure no source (from ids def) + + +% Get subfield + +G.rv = gdat(shot, 'STATIC("R_V")'); % Vessel radial coord 256x1 +G.wv = gdat(shot, 'STATIC("A_V")'); % Vessel conductor width 256x1 +G.zv = gdat(shot, 'STATIC("Z_V")'); % Vessel vert. coord 256x1 +G.hv = gdat(shot, 'STATIC("B_V")'); % Vessel conductor heigth 256x1 +G.Rv = gdat(shot, 'STATIC("RES_V")'); % Vessel resistivity 1 in 256x1 in ohm +G.res = gdat(shot, 'STATIC("RHO_V")'); % Vessel resistance 1 in ohm*m +G.dimv = gdat(shot, 'STATIC("N_V")'); % Number of turns + +Uf = gdat(shot, 'TCV_IDEALLOOP("VLOOP")'); % Flux loop voltage +Tvs = gdat(shot, 'STATIC("T_V_S")'); % Transfer matrix segments to vessel elements +Rs = gdat(shot, 'STATIC("RES_S")'); % Resistance of each segment + +Is = -Uf.data./Rs.data'; % Measured current in each segment +Iv = Tvs.data*Is'; % Mapped to each vessel element + +%ids_core_profiles.profiles_1d{1}.ion{1}.state = {}; +%ids_core_profiles.profiles_1d{1}.neutral{1}.state = {}; +%ids_core_profiles.profiles_1d(1:length(ids_core_profiles.time)) = ids_core_profiles.profiles_1d(1); + +loop_template = ids_pf_passive.loop{1}; + +for ii=1:G.dimv.data + % Duplicate loop substructure + ids_pf_passive.loop{ii} = loop_template; + + ids_pf_passive.loop{ii}.element{1}.geometry.geometry_type = 3; + ids_pf_passive.loop{ii}.element{1}.geometry.oblique.r = G.rv.data(ii)-G.wv.data(ii)/2; + ids_pf_passive.loop{ii}.element{1}.geometry.oblique.z = G.zv.data(ii)-G.hv.data(ii)/2; + ids_pf_passive.loop{ii}.element{1}.geometry.oblique.length_alpha = G.wv.data(ii); + ids_pf_passive.loop{ii}.element{1}.geometry.oblique.length_beta = G.hv.data(ii); + ids_pf_passive.loop{ii}.element{1}.geometry.oblique.alpha = 0; + ids_pf_passive.loop{ii}.element{1}.geometry.oblique.beta = 0; + ids_pf_passive.loop{ii}.element{1}.turns_with_sign = 1; + + ids_pf_passive.loop{ii}.element{1}.name = num2str(ii); + ids_pf_passive.loop{ii}.name = num2str(ii); % might not be correct + + ids_pf_passive.loop{ii}.resistivity = G.Rv.data(ii)/2/pi/G.rv.data(ii)*G.wv.data(ii)*G.hv.data(ii); + ids_pf_passive.loop{ii}.resistance = G.res.data; + ids_pf_passive.loop{ii}.current = Iv(ii,:); + ids_pf_passive.loop{ii}.time = Uf.t; + +end + +%% Code + +if nargin <= 2 + ids_pf_passive.code.name = ['tcv_get_ids_pf_passive, within gdat, with shot= ' num2str(params_pf_passive.shot) ]; +else + ids_pf_passive.code.name = ['tcv_get_ids_pf_passive, within gdat, with shot= ' num2str(params_pf_passive.shot) '; varargin: ' varargin{:}]; +end +ids_pf_passive_description.code.name = ids_pf_passive.code.name; + +ids_pf_passive.code.output_flag = zeros(size(ids_pf_passive.time)); + + + +% Not sure about this, copied from active but might need changes +% cocos automatic transform +if exist('ids_generic_cocos_nodes_transformation_symbolic','file') + [ids_pf_passive,cocoscoeff]=ids_generic_cocos_nodes_transformation_symbolic(ids_pf_passive,'pf_passive',gdat_params.cocos_in, ... + gdat_params.cocos_out,gdat_params.ipsign_out,gdat_params.b0sign_out,gdat_params.ipsign_in,gdat_params.b0sign_in, ... + gdat_params.error_bar,gdat_params.nverbose); + +end + diff --git a/matlab/TCV_IMAS/tcv_get_ids_summary.m b/matlab/TCV_IMAS/tcv_get_ids_summary.m index bda11cc2eceab12cb590b9a086397b16e8967159..0afcb6cc1df3c0cfe1fbcaf2f93baac3b0663bc0 100644 --- a/matlab/TCV_IMAS/tcv_get_ids_summary.m +++ b/matlab/TCV_IMAS/tcv_get_ids_summary.m @@ -250,13 +250,12 @@ for i=1:numel(global_quantities_fieldnames) if ~any(strcmp(global_quantities_fieldnames{i},special_fields)) if ~isstruct(ids_summary.global_quantities.(global_quantities_fieldnames{i}).value) && ... ~isempty(global_quantities.(global_quantities_fieldnames{i}).data) - % setup mask to get rid of nans - mask = ~isnan(global_quantities.(global_quantities_fieldnames{i}).data); - % interpolate quantity on ids_summary.time + %Use interpos_nan to get rid of nans. Arrays will be empty if filled with nans ids_summary.global_quantities.(global_quantities_fieldnames{i}).value = ... - interpos(21,global_quantities.(global_quantities_fieldnames{i}).t(mask),... - global_quantities.(global_quantities_fieldnames{i}).data(mask), ... + interpos_nan(21,global_quantities.(global_quantities_fieldnames{i}).t,... + global_quantities.(global_quantities_fieldnames{i}).data, ... ids_summary.time); + ids_summary.global_quantities.(global_quantities_fieldnames{i}).source = ['gdat request: ' global_quantities_desc.(global_quantities_fieldnames{i})]; elseif ~isempty(global_quantities.(global_quantities_fieldnames{i}).data) special_fields{end+1} = global_quantities_fieldnames{i};