%**************************************************************************
%This program compares the requested field changes to the actual field
%changes caused by a shim.
%
%The program uses the convention that z is along the insertion device, y is
%in the main magnetic field direction, and x is in the nominally field free
%direction.  xyz make a right handed coordinate system.
%
%The planar ID analysis and planar ID shimming programs must be run before
%this program can be run.
%
%Input -----
%
%The following structures are imported:
%
%param, structure containing required program parameters
%	input_file_dir      %Input file directory
%   input_file_name     %Input file name
%   input_file_type     %Input file type, gives file format
%   corr_required       %Logical variable, apply a correction to the measured fields, true or false
%   corr_bx_mult        %Multiplicative correction factor for Bx
%   corr_by_mult        %Multiplicative correction factor for By
%   corr_bx_add         %Additive constant correction for Bx (T)
%   corr_by_add         %Additive constant correction for By (T)
%   corr_bxy_add_file   %Correction file name including directory, additive for Bx and By, 'none' for no correction
%   corr_bxy_add_file_type  %Correction file type, gives format
%   gamma               %Lorentz factor
%   nom_lambda_u        %Nominal undulator period (m)
%   nom_k_value         %Nominal undulator K parameter value
%   nom_lambda_r        %Nominal radiation wavelength (m)
%   num_poles           %Number of poles in the undulator
%   num_end_poles       %Number of poles used to determine end effects
%   cell_length         %Cell length, undulator length plus drift length at each end
%   nom_phase_cell      %Nominal phase change through the cell (rad)
%   save_as_final_result%Save output as a final result, true or false
%   result_file_dir     %Directory where final results are saved
%   result_file_name    %Name of file for final result
%   fig_close_on_exit   %Close figures when program exits, true or false
%   fig_vis             %Figures visible on screen, true or false
%
%header, cell array containing the header of the data file
%
%field, structure containing magnetic field information
%   zm                  %Z positions at which the field measurements are made (m)
%   z                   %Translated Z, with z=0 at undulator center (m)
%   bx                  %Measured horizontal field (T)
%   by                  %Measured vertical field (T)
%   z_peak              %Z positions of the By peak fields (m)
%   by_peak             %Fitted By peak values (T)
%   z_zero              %Z positions of the By zero points (m)
%   lambda_u            %Measured average undulator period (m)
%   int_bx_dz           %First integral of Bx (Tm)
%   int_by_dz           %First integral of By (Tm)
%   int_int_bx_dz_dz    %Second integral of Bx (Tm^2)
%   int_int_by_dz_dz    %Second integral of By (Tm^2)
%   bx_filt             %Filtered measured horizontal field (T)
%   by_filt             %Filtered measured vertical field (T)
%   int_bx_dz_filt      %Filtered first integral of Bx (Tm)
%   int_by_dz_filt      %Filtered first integral of By (Tm)
%   int_int_bx_dz_dz_filt  %Filtered second integral of Bx (Tm^2)
%   int_int_by_dz_dz_filt  %Filtered second integral of By (Tm^2)
%
%traj, structure containing trajectory information
%   z                   %Translated Z, with z=0 at undulator center (m)
%   x                   %Trajectory x position (m)
%   x_p                 %First derivative dx/dz
%   x_pp                %Second derivative d(dx/dz)/dz (m^-1)
%   y                   %Trajectory y position (m)
%   y_p                 %First derivative dy/dz
%   y_pp                %Second derivative d(dy/dz)/dz (m^-1)
%   x_filt              %Filtered trajectory x position (m)
%   x_p_filt            %Filtered first derivative dx/dz
%   x_pp_filt           %Filtered second derivative d(dx/dz)/dz (m^-1)
%   y_filt              %Filtered trajectory y position (m)
%   y_p_filt            %Filtered first derivative dy/dz
%   y_pp_filt           %Filtered second derivative d(dy/dz)/dz (m^-1)
%   x_slope_cor         %X trajectory minus linear fit
%   x_p_slope_cor       %X trajectory slope minus linear fit slope
%   y_slope_cor         %Y trajectory minus linear fit
%   y_p_slope_cor       %Y trajectory slope minus linear fit slope
%
%rad, structure containing radiation information
%	z                   %Translated Z, from traj structure (m)
%	z_peak              %z at the field peaks (m)
%	z_peak_no_ends      %z at the field peaks excluding the ends (m)
%   lambda_u            %Measured average undulator period (m)
%	slip                %Slippage, slip = 0 at z = 0 (m)
%	slip_peak           %Slippage at the field peaks (m)
%	slip_peak_no_ends	%Slippage at the field peaks excluding the ends (m)
%	slope               %Slope of the slip_peak vs z_peak (no ends) linear fit
%	slip_norm           %Normalized slippage, slip_norm = 0 at z = 0 (m)
%	slip_norm_peak      %Normalized slippage at the field peaks (m)
%	slip_norm_peak_no_ends	%Normalized slippage at the field peaks excluding the ends (m)
%	slope_norm          %Slope of the slip_norm_peak vs z_peak (no ends) linear fit
%	k_eff               %K value for undulator
%	b_eff               %Effective field strength in undulator (T)
%	lambda_r            %Radiation wavelength (m)
%	phase               %Phase, phase = 0 at cell boundary (rad)
%	phase_peak          %Phase at the field peaks (rad)
%	phase_peak_no_ends	%Phase at the field peaks excluding the ends (rad)
%	phase_error_peak_no_ends	%Phase errors, error at first pole after end defines zero error (rad)
%	phase_cell          %Total phase change between the light and the electrons in the cell (rad)
%	phase_nom           %Phase relative to nominal radiation wave, phase = 0 at cell boundary (rad)
%	phase_nom_peak      %Phase relative to nominal radiation wave at the field peaks (rad)
%	phase_nom_peak_no_ends  %Phase relative to nominal radiation wave at the field peaks excluding the ends (rad)
%	phase_nom_error_peak_no_ends  %Phase errors relative to nominal radiation wave, error at first pole after end defines zero error (rad)
%	phase_nom_cell      %Total phase change between the nominal radiation wave and the electrons in the cell (rad)
%
%shim_param, structure containing shim parameters
%	input_file_dir      %Input file directory
%   input_file_name     %Input file name
%   traj_num_seg        %Number of segments in the segmented trajectory
%   traj_shift_max      %Max trajectory shift when a shim is moved
%   traj_shift_accum_max  %Max trajectory shift accumulation before a shim is placed
%   x_traj_shim_par     %X trajectory shim parameters, determined experimentally
%   y_traj_shim_par     %Y trajectory shim parameters, determined experimentally
%   phaserr_num_seg     %Number of segments in the segmented phase error
%   phaserr_pos_delta_max  %Max positive phase error change before a shim is added
%   phaserr_neg_delta_max  %Max negative phase error change before a shim is added
%   phaserr_gap_shim_width %Effective width of a gap shim (m)
%   phaserr_gap_shim_nom_gap %Nominal undulator gap, used in shim calculations (m)
%   phaserr_gap_shim_corr_fact %Multiplicative correction to the gap change accounting for actual shape
%   phaserr_phase_shim_par   %Phase error phase shim parameters, determined experimentally
%   k_shim_cant_angle   %Undulator pole cant angle (rad)
%   phase_match_phase   %Desired phase on the pole after the end poles (rad)
%
%setgap_shim, structure containing shim information to set the gap
%   z_peak_no_ends      %Peak field z positions, no ends
%   dbrel_no_ends       %Relative field change required at the poles
%   dgap_no_ends        %Required gap change at the poles
%   dx_no_ends          %Required tapered shim motions at the poles
%
%seg_traj, structure containing the segmented trajectory information
%   x_seg               %structure containing the x segmented trajectory
%   y_seg               %structure containing the y segmented trajectory
%x_seg, y_seg
%   traj_name           %'X' or 'Y'
%   num_seg             %Number of segments
%   z_start             %Z position of the start of each segment
%   z_end               %Z position of the end of each segment
%   traj_slope_orig     %Trajectory slope of the segment of the original trajectory
%   traj_start_orig     %Trajectory position at the start of the segment of the original trajectory
%   traj_end_orig       %Trajectory position at the end of the segment of the original trajectory
%   shim_slope          %Shim at each segment, effect on slope
%   shim_bl             %Shim at each segment, integrated field strength
%   traj_slope_shimmed  %Trajectory slope of the segment after shimming
%   traj_start_shimmed  %Trajectory position at the start of the segment after shimming
%   traj_end_shimmed    %Trajectory position at the end of the segment after shimming
%
%traj_shim, structure containing the trajectory shim information
%   x_shim              %Structure containing the x trajectory shim information
%   y_shim              %Structure containing the y trajectory shim information
%x_shim, structure containing the x-trajectory shim information
%   num_shim            %Number of shims
%   pole_num            %Pole number to place the shim
%   pole_z_pos          %z position of the pole
%   pole_by_sign        %Sign of the field in the gap at the pole_num location
%   shim_slope          %Shim effect on trajectory slope
%   shim_bl             %Shim integrated field strength
%   shim_type           %Shim type (types given in parameter file)
%   shim_dim            %Shim physical dimension
%   shim_d_dim          %Shim physical dimension change for pre-existing shims
%y_shim, structure containing the y-trajectory shim information (The shim uses two poles.)
%   num_shim            %Number of shims
%   pole_1_num          %Pole number of upstream pole in shim
%   pole_2_num          %Pole number of downstream pole in shim
%   pole_1_z_pos        %Z position of upstream pole in shim
%   pole_2_z_pos        %Z position of downstream pole in shim
%   pole_1_by_sign      %Sign of the field at the upstream pole
%   pole_2_by_sign      %Sign of the field at the downstream pole
%   pole_1_loc          %Location of shim on upstream pole, +1 upper, -1 lower
%   pole_2_loc          %Location of shim on downstream pole, +1 upper, -1 lower
%   shim_slope          %Shim effect on trajectory slope
%   shim_bl             %Shim integrated field strength
%   shim_type           %Shim type
%   shim_dim            %Shim physical dimension
%   shim_d_dim          %Shim physical dimension change for pre-existing shims
%
%seg_phaserr, structure containing the segmented phase error information
%   num_seg             %Number of segments
%   z_start             %Z position of the start of each segment
%   z_end               %Z position of the end of each segment
%   phaserr_slope_orig  %Phase error slope of the segment before shimming
%   phaserr_start_orig  %Phase error at the start of the segment before shimming
%   phaserr_end_orig    %Phase error at the end of the segment before shimming
%   gap_shim_num        %Gap shim number, gap shims can extend over several segments
%   gap_shim_slope      %Gap shim at each segment, effect on slope
%   gap_shim_dbrel      %Gap shim at each segment, delta B / B
%   phase_shim_delta    %Phase shim at each segment, change in phase
%   phase_shim_dbrel    %Phase shim at each segment, delta B / B
%   phaserr_slope_shimmed  %Phase error slope of the segment after shimming
%   phaserr_start_shimmed  %Phase error at the start of the segment after shimming
%   phaserr_end_shimmed    %Phase error at the end of the segment after shimming
%
%phaserr_shim, structure containing the shim information
%   num_gap_shim        %Number of gap shims
%   gap_shim_z_start    %Z position the gap shim starts at
%   gap_shim_z_end      %Z position the gap shim ends at
%   gap_shim_pole_start %Pole number the gap shim starts at
%   gap_shim_pole_end   %Pole number the gap shim ends at
%   gap_shim_slope      %Phase error slope change from the gap shim
%   gap_shim_dbrel      %Fractional field change from the gap shim
%   gap_shim_dim        %Gap size change from the shim
%   num_phase_shim      %Number of phase shims
%   phase_shim_z        %Z position of the phase shim
%   phase_shim_pole_start  %Pole number the phase shim starts at
%   phase_shim_pole_end %Pole number the phase shim ends at
%   phase_shim_delta    %Phase error change from the phase shim
%   phase_shim_dbrel    %Fractional field change from the phase shim
%   phase_shim_dim      %Dimension of the phase shim
%
%k_shim, structure containing shim information to get the proper K value
%   dbrel               %Required dB / B to get the proper K value
%   dgap                %Required gap change (m)
%   dxbeam              %Change in x of the beam required to give the proper K value
%
%phasematch_shim, structure containing the phase matching shim information
%   expect_phase_poles  %expected phases on the poles based on the phase match phase
%   expect_phase_cell   %expected phase of the cell based on the pahse match phase
%   phase_err_enter     %phase errors on the poles entering the undulator (rad)
%   phase_err_exit      %phase errors on the poles leaving the undulator (rad)
%   des_phase_change_enter  %desired phase change from the shim at the entrance (rad)
%   des_phase_change_exit   %desired phase change from the shim at the exit (rad)
%   shim_dim_enter      %size of the shim required to give the phase change at the entrance (m)
%   shim_dim_exit       %size of the shim required to give the phase change at the exit (m)
%
%Zachary Wolf, 04/02/07

function [] = shim_effect_anal()

%Clear the workspace
clear

%Set the parameter values
shim_effect_anal_param;

%Output files
plt_file = [INPUT_FILE_DIR_AFT, 'shim_effect_anal.ps'];

%Initialize the plots
shim_effect_plot_init(plt_file, SHIM_NAME, ANAL_NAME);

%Input the analysis results for the run before the shim
anal_file_bef = [INPUT_FILE_DIR_BEF, 'id_anal.mat'];
load(anal_file_bef, 'param', 'header', 'field', 'traj', 'rad');
shim_file_bef = [INPUT_FILE_DIR_BEF, 'id_shim.mat'];
load(shim_file_bef, 'setgap_shim', 'seg_traj', 'traj_shim', 'seg_phaserr', 'phaserr_shim', 'k_shim', 'phasematch_shim');
param_bef = param;
header_bef = header;
field_bef = field;
traj_bef = traj;
rad_bef = rad;
setgap_shim_bef = setgap_shim;
seg_traj_bef = seg_traj;
traj_shim_bef = traj_shim;
seg_phaserr_bef = seg_phaserr;
phaserr_shim_bef = phaserr_shim;
k_shim_bef = k_shim;
phasematch_shim_bef = phasematch_shim;

%Input the analysis results for the run after the shim
anal_file_aft = [INPUT_FILE_DIR_AFT, 'id_anal.mat'];
load(anal_file_aft, 'param', 'header', 'field', 'traj', 'rad');
shim_file_aft = [INPUT_FILE_DIR_AFT, 'id_shim.mat'];
load(shim_file_aft, 'setgap_shim', 'seg_traj', 'traj_shim', 'seg_phaserr', 'phaserr_shim', 'k_shim', 'phasematch_shim');
param_aft = param;
header_aft = header;
field_aft = field;
traj_aft = traj;
rad_aft = rad;
setgap_shim_aft = setgap_shim;
seg_traj_aft = seg_traj;
traj_shim_aft = traj_shim;
seg_phaserr_aft = seg_phaserr;
phaserr_shim_aft = phaserr_shim;
k_shim_aft = k_shim;
phasematch_shim_aft = phasematch_shim;

%Set gap shim
if (strcmp(SHIM_NAME, 'SET_GAP') == 1)
    dbrel_no_ends_requested = setgap_shim_bef.dbrel_no_ends;
    dbrel_actual = (field_aft.by_peak - field_bef.by_peak) ./ field_bef.by_peak;
    dbrel_no_ends_actual = dbrel_actual((param_bef.num_end_poles + 1) : (end - param_bef.num_end_poles));
    pole_num = 1 : param_bef.num_poles;
    pole_num_no_ends = pole_num((param_bef.num_end_poles + 1) : (end - param_bef.num_end_poles));
	shim_effect_plot_setgap(plt_file, ANAL_NAME, pole_num_no_ends, dbrel_no_ends_requested, dbrel_no_ends_actual);
end

%X trajectory shim
if (strcmp(SHIM_NAME, 'X_TRAJ') == 1)
    z = field_bef.z;
    Iy_bef = field_bef.int_by_dz_filt;
    Iy_aft = field_aft.int_by_dz_filt;
    shim_bl = zeros(1, length(z));
    for ii = 1 : traj_shim_bef.x_shim.num_shim
        [min_val, index] = min(abs(z - traj_shim_bef.x_shim.pole_z_pos(ii)));
        shim_bl(index) = traj_shim_bef.x_shim.shim_bl(ii);
    end
    Iy_shim = cumsum(shim_bl);
    shim_effect_plot_xtraj(plt_file, ANAL_NAME, z, Iy_bef, Iy_aft, Iy_shim)
end

%Y trajectory shim
if (strcmp(SHIM_NAME, 'Y_TRAJ') == 1)
    z = field_bef.z;
    Ix_bef = field_bef.int_bx_dz_filt;
    Ix_aft = field_aft.int_bx_dz_filt;
    shim_bl = zeros(1, length(z));
    for ii = 1 : traj_shim_bef.y_shim.num_shim
        [min_val, index] = min(abs(z - (traj_shim_bef.y_shim.pole_1_z_pos(ii) + traj_shim_bef.y_shim.pole_2_z_pos(ii)) / 2.));
        shim_bl(index) = traj_shim_bef.y_shim.shim_bl(ii);
    end
    Ix_shim = cumsum(shim_bl);
    shim_effect_plot_ytraj(plt_file, ANAL_NAME, z, Ix_bef, Ix_aft, Ix_shim)
end

%Phase shim
if (strcmp(SHIM_NAME, 'PHASE') == 1)
    pole_num = 1 : param_bef.num_poles;
    dbrel_actual = (field_aft.by_peak - field_bef.by_peak) ./ field_bef.by_peak;
    dbrel_requested = zeros(1, length(pole_num));
    for ii = 1 : phaserr_shim_bef.num_gap_shim
        index_array = find((pole_num >= phaserr_shim_bef.gap_shim_pole_start(ii)) & (pole_num <= phaserr_shim_bef.gap_shim_pole_end(ii)));
        dbrel_requested(index_array) = phaserr_shim_bef.gap_shim_dbrel(ii);
    end
    for ii = 1 : phaserr_shim_bef.num_phase_shim
        index_array = find((pole_num >= phaserr_shim_bef.phase_shim_pole_start(ii)) & (pole_num <= phaserr_shim_bef.phase_shim_pole_end(ii)));
        dbrel_requested(index_array) = phaserr_shim_bef.phase_shim_dbrel(ii);
    end
	shim_effect_plot_phaserr(plt_file, ANAL_NAME, pole_num, dbrel_requested, dbrel_actual);
end

%Set K
if (strcmp(SHIM_NAME, 'SET_K') == 1)
    pole_num = 1 : param_bef.num_poles;
    pole_num_no_ends = pole_num((param_bef.num_end_poles + 1) : (end - param_bef.num_end_poles));
    dbrel_no_ends_requested = k_shim_bef.dbrel * ones(1, length(pole_num_no_ends));
    dbrel_actual = (field_aft.by_peak - field_bef.by_peak) ./ field_bef.by_peak;
    dbrel_no_ends_actual = dbrel_actual((param_bef.num_end_poles + 1) : (end - param_bef.num_end_poles));
	shim_effect_plot_setk(plt_file, ANAL_NAME, pole_num_no_ends, dbrel_no_ends_requested, dbrel_no_ends_actual);
end

%Phasematch shim
if (strcmp(SHIM_NAME, 'PHASEMATCH') == 1)
    dbrel_actual = (field_aft.by_peak - field_bef.by_peak) ./ field_bef.by_peak;
    pole_num = 1 : param_bef.num_poles;
	shim_effect_plot_phasematch(plt_file, ANAL_NAME, pole_num, dbrel_actual);
end

%Done
