%**************************************************************************
%pm_block_sim
%This program simulates the magnetic field from a permanent magnet
%block.  Magnet block parameters are entered.  Equivalent charges are
%calculated to replace the magnet blocks.  The field from each charge is
%calculated and added to the field from the other charges.
%
%The following structures are used.
%
%param, structure containing program parameters
%   output_file_dir %Directory for output files
%   num_block_xdiv  %Number of divisions in x for charge calculation
%   num_block_ydiv  %Number of divisions in y for charge calculation
%   num_block_zdiv  %Number of divisions in z for charge calculation
%   br              %Permanent magnet remnant field (T)
%   block_x_width   %Magnet block width in x direction (m)
%   block_y_width   %Magnet block width in y direction (m)
%   block_z_width   %Magnet block width in z direction (m)
%   x_field_pts     %X positions where field is calculated (m)
%   y_field_pts     %Y positions where field is calculated (m)
%   z_field_pts     %Z positions where field is calculated (m)
%   field_pts_x     %Vector of positions where field is calculated, x coordinate (m)
%   field_pts_y     %Vector of positions where field is calculated, y coordinate (m)
%   field_pts_z     %Vector of positions where field is calculated, z coordinate (m)
%
%block, structure containing the permanent magnet block information
%   num_block       %Number of blocks
%   mx              %X component of the magnetic moment density (A/m)
%   my              %Y component of the magnetic moment density (A/m)
%   mz              %Z component of the magnetic moment density (A/m)
%   x_center        %X position of the block center (m)
%   y_center        %Y position of the block center (m)
%   z_center        %Z position of the block center (m)
%   x_width         %Block width in the x direction (m)
%   y_width         %Block width in the y direction (m)
%   z_width         %Block width in the z direction (m)
%
%charge, structure containing the charge information
%   num_charge      %Number of charges
%   q               %Magnetic charge, M.n dA (Am)
%   x               %X position of charge (m)
%   y               %Y position of charge (m)
%   z               %Z position of charge (m)
%
%field_pos, structure giving the positions at which to calculate the field
%   x               %X position of the field point (m)
%   y               %Y position of the field point (m)
%   z               %Z position of the field point (m)
%
%mag_field, structure containing the calculated magnetic field
%   x               %X position of the field point (m)
%   y               %Y position of the field point (m)
%   z               %Z position of the field point (m)
%   bx              %Bx magnetic field component at (x, y, z) (T)
%   by              %By magnetic field component at (x, y, z) (T)
%   bz              %Bz magnetic field component at (x, y, z) (T)
%
%Zachary Wolf, 4/11/08

%Clear all variables
clear;

%Set the working directory
cd('c:\Simulation\PM Block Sim');

%Global
global mu_0;
mu_0 = 4 * pi * 10^-7;  %(Tm/A)

%Program parameters
[param] = pm_block_sim_get_param();

%Initialize output files
[dat_file, plt_file] = pm_block_sim_file_init(param);

%Initialize the permanent magnet blocks
[block] = pm_block_sim_block_init(param);

%Calculate charges on blocks
[charge] = pm_calc_charge(param, block);

%Add radiation damage
%for ii = 1 : length(charge.q)
%    if charge.x(ii)^2 + charge.y(ii)^2 + charge.z(ii)^2 < 0.005^2
%        charge.q(ii) = 0;
%    end
%end

%Calculate fields from the charges
[mag_field] = pm_calc_mag_field(param, charge);

%Clean up the plots
%index_array = find((mag_field.z > -.02 - .001) & (mag_field.z < -.02 + .001) & abs(mag_field.y) <= .05);
%mag_field.bx(index_array) = 0;
%mag_field.by(index_array) = 0;
%mag_field.bz(index_array) = 0;

%Plot the field
pm_plot_bz_vs_xy_on_z_slice(plt_file, param, block, charge, mag_field, 0.008);
pm_plot_b_on_x_slice(plt_file, block, charge, mag_field, 0);
pm_plot_b_on_z_slice(plt_file, block, charge, mag_field, 0.008);
pm_plot_bz_vs_x_at_yz(plt_file, block, charge, mag_field, 0, 0.008);
pm_plot_bz_vs_x_at_yz(plt_file, block, charge, mag_field, 0.022, 0.008);

%Write the results to a data file
%ppm_dat_byx_vs_z(dat_file, mag_field, 0, 0);

%Output the results to the output file directory
mat_file = strcat(lower(param.output_file_dir), 'pm_block_sim.mat');
save(mat_file, 'param', 'block', 'charge', 'mag_field');

%Done
