function [R,Z] = gyrotronRadiusTDM(theta,d,m)
%[R,z] = gyrotronRadiusTDM(theta,d,m)
% Calculates the radius of multi-section cylindrical shell 
% and gives the corresponding axial z-values for a user-
% specified number of interpolation points in each section.
% The inputs are 'Theta' (angle), 'Distance', and 'Mesh'
% hence the 'TDM' in the function name.
%
% Inputs: (must be the same size)
%   'd': [distance] the corresponding axial length 
%   vector defining the length of each section.
%   'theta': [radians] a vector of theta angles 
%   that describe the increase in radius with increasing 
%   axial distance.
%   'm': [] the mesh size, that is, the number of points
%   of interpolation for each section.
% 
% Outputs:
%   'R': [distance] the radius vector at the axial 
%   points given by the vector 'z'.
%   'z': [distance] the axial distance vector.
%
% Example:
%   d     = [2 1 1.5  1  1  1  1]; 
%   theta = [1 2 0.7 -1 -2 -1  0]*pi/180; 
%   m     = [1 3 1    2  3  1  1]*20;
%   [R,z] = gyrotronRadiusTDM(theta,d,m);
%   plot(z,R,'b.',z,R,'g');
%
% Last edit by Eunmi Choi and Colin Joye 7/15/03

% 'g' is the radius offset for each new section to make it continuous
g = cumsum( d.*tan(theta) );
g = [0 g(1:(end-1))];

% must extend distance array for first point at origin.
d = [0 d]; 

% initialize output vectors.
R = [];
Z = [];

% loop to generate outputs:
%   'inc' is the step size based on the input mesh number
%   'z' is the local values to generate the incremental radius
for k = 1:length(g),
  inc = d(k+1)/m(k);
  z   = [0 : inc : (d(k+1)-inc)];
  R   = [R z*tan(theta(k))+g(k)];
  Z   = [Z z+sum(d(1:k))];
end
% ------------ End gyrotronRadiusTDM.m -------------