%**************************************************************************
%This function performs a phase analysis of a device.
%The slippage, the distance between the light and the electrons, is calculated.
%The phase is analyzed, including the total phase through the device.
%
%Input:
%param, parameter structure
%dat_file, output data file name
%plt_file, output plot file name
%field, field structure
%traj, trajectory structure
%
%Output:
%phas, structure containing phase 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)
%	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)
%	nom_lambda_r        %Wavelength of nominal radiation wave (m)
%	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_cell      %Total phase change between the nominal radiation wave and the electrons in the cell (rad)
%   phase_integral      %Int_-inf^z [(Int_-inf^z'' Bx(z') dz')^2 + (Int_-inf^z'' By(z') dz')^2] dz'' (T^2m^3)
%   phase_integral_end  %Int_-inf^+inf [(Int_-inf^z'' Bx(z') dz')^2 + (Int_-inf^z'' By(z') dz')^2] dz'' (T^2m^3)
%   phase_integral_cell %Int_z_enter^z_exit [(Int_-inf^z'' Bx(z') dz')^2 + (Int_-inf^z'' By(z') dz')^2] dz'' (T^2m^3)
%
%Zachary Wolf
%5/20/11

function [phas] = phase_anal_phase(param, dat_file, plt_file, field, traj)

%Extract required information
gamma = param.gamma;
nom_lambda_r = param.nom_lambda_r;
cell_length = param.cell_length;
z = field.z;
z_peak = field.z_peak;
z_peak_no_ends = field.z_peak_no_ends;
bx = field.bx;
by = field.by;
int_bx_dz = field.int_bx_dz;
int_by_dz = field.int_by_dz;
vz = traj.vz;

%Slippage -----

%Calculate the slippage as a function of z
slip_vz = slip_calc_slip_vz(z, vz);
% slip_slope_o1 = slip_calc_slip_slope_o1(gamma, z, x_p, y_p);  %trajectory slope, order 1
% slip_slope = slip_calc_slip_slope(gamma, z, x_p, y_p);
slip = slip_vz;

%Debug
% slip_diff = slip_vz - slip_slope_o1;
% slip_plot_slip_diff(plt_file, z, slip_diff);
% slip_diff = slip_vz - slip_slope;
% slip_plot_slip_diff(plt_file, z, slip_diff);
% slip_diff = slip_slope - slip_slope_o1;
% slip_plot_slip_diff(plt_file, z, slip_diff);

%Find the slippage at the field peaks
slip_peak = slip_calc_slip_at_peaks(z, slip, z_peak);
slip_peak_no_ends = slip_calc_slip_at_peaks(z, slip, z_peak_no_ends);

%Plot the slippage
slip_plot_slip(plt_file, z, slip);
%slip_plot_slip_pole(plt_file, z, slip, z_peak);

%Phase To Nominal Wave -----

%Calculate the phase relative to the nominal radiation wave as a function of z, phase(z(1)) = 0 ref
phase_nom = phase_calc_phase(nom_lambda_r, slip);

%Calculate the phase relative to the nominal radiation wave at the field peaks
phase_nom_peak = phase_calc_phase(nom_lambda_r, slip_peak);

%Calculate the phase relative to the nominal radiation wave at the field peaks, neglect the ends
phase_nom_peak_no_ends = phase_calc_phase(nom_lambda_r, slip_peak_no_ends);

%Calculate the phase relative to the nominal radiation wave at the cell entrance
phase_nom_cell_entrance = phase_calc_phase_cell_entrance(nom_lambda_r, gamma, cell_length, z, phase_nom);

%Calculate the phase relative to the nominal radiation wave at the cell exit
phase_nom_cell_exit = phase_calc_phase_cell_exit(nom_lambda_r, gamma, cell_length, z, phase_nom);

%Make the cell entrance phase_nom = 0 reference
phase_nom = phase_nom - phase_nom_cell_entrance;
phase_nom_peak = phase_nom_peak - phase_nom_cell_entrance;
phase_nom_peak_no_ends = phase_nom_peak_no_ends - phase_nom_cell_entrance;
phase_nom_cell_exit = phase_nom_cell_exit  - phase_nom_cell_entrance;
phase_nom_cell_entrance = 0;  %Keep at end

%Calculate the phase change through the cell
phase_nom_cell = phase_nom_cell_exit;

%Plot the phase
phase_plot_phase_nom(plt_file, z, phase_nom);

%Write the cell phase results to the data file
phase_dat_phase_nom_cell(dat_file, nom_lambda_r, cell_length, phase_nom_cell);

%Phase Integral -----

%Calculate the phase integral
phase_integral = int_calc_int_y_dx(z, int_bx_dz.^2 + int_by_dz.^2);
phase_integral_end = phase_integral(end);
phase_integral_cell_entrance = spline(z, phase_integral, - cell_length / 2);
phase_integral_cell_exit = spline(z, phase_integral, cell_length / 2);
phase_integral_cell = phase_integral_cell_exit - phase_integral_cell_entrance;

%Plot the phase integral
phase_plot_phase_integral(plt_file, z, phase_integral);

%Write the phase integral to the data file
phase_dat_phase_integral(dat_file, phase_integral_end, phase_integral_cell);

%Results -----

%Save values in the phas structure
phas.z = z;
phas.z_peak = z_peak;
phas.z_peak_no_ends = z_peak_no_ends;
phas.slip = slip;
phas.slip_peak = slip_peak;
phas.slip_peak_no_ends = slip_peak_no_ends;
phas.nom_lambda_r = nom_lambda_r;
phas.phase_nom = phase_nom;
phas.phase_nom_peak = phase_nom_peak;
phas.phase_nom_peak_no_ends = phase_nom_peak_no_ends;
phas.phase_nom_cell = phase_nom_cell;
phas.phase_integral = phase_integral;
phas.phase_integral_end = phase_integral_end;
phas.phase_integral_cell = phase_integral_cell;

%Done
