function [ids_struct_out,ids_struct_out_description] =  tcv_ids_coil(shot, ids_structures, gdat_params, varargin)
%
% [ids_struct_out,ids_struct_out_desccription] =  tcv_get_ids_coil(shot, ids_structures, gdat_params, varargin)
%
%  pf_active coils
%
% error_bar: from gdat_params.error_bar
%            'delta' (default): error_bar to be added inserted in "upper" only as mentioned in description
%            'delta_with_lower' : error_bar (abs) inserted in both lower and upper
%            'added': value already added to data: upper/lower = data +/- error_bar
%

% TODO add logic to get the current in the coils only once.

error_bar = 'delta';
if exist('gdat_params') && isfield(gdat_params,'error_bar') && ~isempty(gdat_params.error_bar)
  error_bar = gdat_params.error_bar;
end

% Coils that can be characterized by R, Z and a
% crosssectional area are described as distinct coils with a single element, this also corresponds to coils with distinct TCV names).
% Elements can be used in the future to refine the spatial
% description of each coil.
% The circular conductor of G-coils approximated by a square of equal
% crosssectional area.
% Each coil has a positive turnsign. The return current in T_003 is
% dealt with in the circuit description and in the machine mapping.
% In practice the dissipated energy in the pf '
% coils is not a relevant limit.
mdsopen(shot)

%% Get power supply/coils names for each circuit.
[pf_def] =  tcv_ids_pf_active_definition();
coil_names2ids = pf_def.coil_names;
mds_paths2ids = pf_def.mds_paths;

ncircuits2ids = numel(coil_names2ids);

%% Get geometrical data from the static tree
% Information from the static tree - Z-zeroed with respect to coils
r_c     = mdsvalue('_r_c=static(''r_c'')'); % R position
z_c     = mdsvalue('_z_c=static(''z_c'')'); % Z position
w_c     = mdsvalue('static(''w_c'')'); % width rectangular description
h_c     = mdsvalue('static(''h_c'')'); % height rectangular description
a_c     = mdsvalue('static(''a_c'')'); % tilt angle for parallelogram representation ( deprecated representation in IDS)
nt_c    = mdsvalue('static(''nt_c'')'); % number of turns
xsect_c = mdsvalue('static(''xsect_c'')');
res_c   = mdsvalue('static(''res_c'')'); % resistence of the coil
sizepfc = numel(r_c);
namepfc = mdsvalue('dim_of(_r_c)');

% Set effective turns in T1 and T2 (see J.-M. Moret, et al., RSI 69 (1998) 2333)
iT = strmatch('T_001',namepfc);
nt_c(iT)=26/68;
iT = strmatch('T_002',namepfc);
nt_c(iT)=42/68;
iT = strmatch('T_003',namepfc);      % Return current
nt_c(iT)=1;

% Approximate circular G-coil conductor with a square of the same crosssection
iG=strmatch('G_00',namepfc);
w_c(iG)=sqrt(xsect_c(iG));
h_c(iG)=sqrt(xsect_c(iG));

%% Put data to ids structure
ids_struct_out(1:sizepfc) = ids_structures;
ind_coil_ids = 0;
for ii=1:ncircuits2ids
    ncoil2ids = numel(coil_names2ids{ii}); % number of coils for a given circuit

    for jj = 1:ncoil2ids
      ind_coil_ids = ind_coil_ids +1;
      ids_struct_out{ind_coil_ids}.name = coil_names2ids{ii}{jj};
      ids_struct_out_description{ind_coil_ids}.name = ['through aa=tcv_ids_pf_active_definition from aa.coil_names  '];
      
      % time-varying data
      if shot == -1
        tmpdata.dim{1} = [];
        tmpdata.data = [];
        warning('no time data loaded for shot %d',shot);
      else
        tmpdata = tdi(mds_paths2ids{ii});
      end
      ids_struct_out{ind_coil_ids}.current.data = tmpdata.data;
      ids_struct_out_description{ind_coil_ids}.current.data = ['from ' mds_paths2ids{ii}];
      ids_struct_out{ind_coil_ids}.current.time = tmpdata.dim{1};

      
      % Find index on static tree
      tmpind = find(strcmp(ids_struct_out{ind_coil_ids}.name, namepfc));
      ids_struct_out{ind_coil_ids}.element{1}.geometry.rectangle.r = r_c(tmpind);
      ids_struct_out_description{ind_coil_ids}.element{1}.geometry.rectangle_r = ['from static(''r_c'')'];
      ids_struct_out{ind_coil_ids}.element{1}.geometry.rectangle.z = z_c(tmpind);
      ids_struct_out_description{ind_coil_ids}.element{1}.geometry.rectangle_z = ['from static(''z_c'')'];
      ids_struct_out{ind_coil_ids}.element{1}.geometry.rectangle.width = w_c(tmpind);
      ids_struct_out_description{ind_coil_ids}.element{1}.geometry.rectangle_width = ['from static(''w_c'')'];
      ids_struct_out{ind_coil_ids}.element{1}.geometry.rectangle.height = h_c(tmpind);
      ids_struct_out_description{ind_coil_ids}.element{1}.geometry.rectangle_height = ['from static(''h_c'')'];
      ids_struct_out{ind_coil_ids}.element{1}.turns_with_sign = nt_c(tmpind);
      ids_struct_out_description{ind_coil_ids}.element{1}.turns_with_sign = ['from static(''nt_c'')'];
      ids_struct_out{ind_coil_ids}.element{1}.geometry.geometry_type = 2; % 1 outline, 2 rectangle, 4 arc of circle
      ids_struct_out_description{ind_coil_ids}.element{1}.geometry.geometry_type = 'rectangle';
      
      ids_struct_out{ind_coil_ids}.resistance = res_c(tmpind);
      ids_struct_out_description{ind_coil_ids}.resistance = ['from static(''res_c'')'];
    end
end

% fixed_error = 200.;
ind_coil_ids = 0;
switch error_bar
 case 'delta'
  for ii=1:ncircuits2ids
    ncoil2ids = numel(coil_names2ids{ii}); % number of coils for a given circuit
    for jj = 1:ncoil2ids
      fixed_error = pf_def.coil_current_error{ii}(jj);
      ind_coil_ids = ind_coil_ids +1;
      ids_struct_out{ind_coil_ids}.current.data_error_upper = fixed_error.*ones(size(ids_struct_out{ind_coil_ids}.current.data));
      ids_struct_out_description{ind_coil_ids}.current_data_error_upper = fixed_error.*ones(size(ids_struct_out{ind_coil_ids}.current.data));
      ids_struct_out_description{ind_coil_ids}.current_data_error_lower = ['not provided since symmetric'];
    end
  end
 case 'delta_with_lower'
  for ii=1:ncircuits2ids
    ncoil2ids = numel(coil_names2ids{ii}); % number of coils for a given circuit
    for jj = 1:ncoil2ids
      fixed_error = pf_def.coil_current_error{ii}(jj);
      ind_coil_ids = ind_coil_ids +1;
      ids_struct_out{ind_coil_ids}.current.data_error_upper = fixed_error.*ones(size(ids_struct_out{ind_coil_ids}.current.data));
      ids_struct_out{ind_coil_ids}.current.data_error_lower = ids_struct_out{ind_coil_ids}.current.data_error_upper;
      ids_struct_out_description{ind_coil_ids}.current_data_error_upper = ['from fixed error value in case ' error_bar];
      ids_struct_out_description{ind_coil_ids}.current_data_error_lower = ['from fixed error value in case ' error_bar];
    end
  end
 case 'added'
  for ii=1:ncircuits2ids
    ncoil2ids = numel(coil_names2ids{ii}); % number of coils for a given circuit
    for jj = 1:ncoil2ids
      fixed_error = pf_def.coil_current_error{ii}(jj);
      ind_coil_ids = ind_coil_ids +1;
      ids_struct_out{ind_coil_ids}.current.data_error_upper = ids_struct_out{ind_coil_ids}.current.data ...
          + fixed_error.*ones(size(ids_struct_out{ind_coil_ids}.current.data));
      ids_struct_out{ind_coil_ids}.current.data_error_lower = ids_struct_out{ind_coil_ids}.current.data ...
          - fixed_error.*ones(size(ids_struct_out{ind_coil_ids}.current.data));
      ids_struct_out_description{ind_coil_ids}.current_data_error_upper = ['from data + fixed error value in case ' error_bar];
      ids_struct_out_description{ind_coil_ids}.current_data_error_lower = ['from data - fixed error value in case ' error_bar];
    end
  end
 otherwise
  error(['tcv_ids_bpol_loop: error_bar option not known: ' error_bar])
end
