Home > pubtools > calc_TouschekPM.m

calc_TouschekPM

PURPOSE ^

tauT = calc_TouschekPM(TD,dppPM,Trf,Ib,U0,coupling, sigE, emit_x)

SYNOPSIS ^

function tauT = calc_TouschekPM(TD,dppPM,Trf,Ib,U0,coupling, sigE, emit_x)

DESCRIPTION ^

tauT = calc_TouschekPM(TD,dppPM,Trf,Ib,U0,coupling, sigE, emit_x)
tauT = calc_TouschekPM(TD,dppPM,alpha,Ib,U0,coupling, sigE, emit_x)
        Ib, mA, single bunch current
        U0, MeV, one-turn energy loss
       emit_x, nm-rad
       coupling, average emit_y/emit_x
       TD, lattice function structure (same as twissring output) 
       dppPM, Nx2, positive/negative momentum aperture
       Trf, a structure or a scalar, when a structure, it consists of [Frequency, HarmNumber, Energy, alpha,
       Voltage], when a scalar it is alpha of the lattice (other
       parameters are assumed default of SPEAR3).
ex:    
        load('dppAP_LAT_DW_withIDs_Mar20_07_normrf.mat')
        [td, tune,chrom] = twissring(THERING,0,indextab, 'chrom', 1e-5);
        dppPM = [deltap' deltam'];
        tauT = calc_TouschekPM(td,dppPM,a1,100/280,1.04,0.064e-2, 0.001, 18)/3600; %hrs
        disp(tauT)

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function tauT = calc_TouschekPM(TD,dppPM,Trf,Ib,U0,coupling, sigE, emit_x)
0002 %tauT = calc_TouschekPM(TD,dppPM,Trf,Ib,U0,coupling, sigE, emit_x)
0003 %tauT = calc_TouschekPM(TD,dppPM,alpha,Ib,U0,coupling, sigE, emit_x)
0004 %        Ib, mA, single bunch current
0005 %        U0, MeV, one-turn energy loss
0006 %       emit_x, nm-rad
0007 %       coupling, average emit_y/emit_x
0008 %       TD, lattice function structure (same as twissring output)
0009 %       dppPM, Nx2, positive/negative momentum aperture
0010 %       Trf, a structure or a scalar, when a structure, it consists of [Frequency, HarmNumber, Energy, alpha,
0011 %       Voltage], when a scalar it is alpha of the lattice (other
0012 %       parameters are assumed default of SPEAR3).
0013 %ex:
0014 %        load('dppAP_LAT_DW_withIDs_Mar20_07_normrf.mat')
0015 %        [td, tune,chrom] = twissring(THERING,0,indextab, 'chrom', 1e-5);
0016 %        dppPM = [deltap' deltam'];
0017 %        tauT = calc_TouschekPM(td,dppPM,a1,100/280,1.04,0.064e-2, 0.001, 18)/3600; %hrs
0018 %        disp(tauT)
0019 
0020 e0 = 1.6e-19; %Coulomb
0021 cspeed = 299792458; 
0022 r0 = 2.82e-15; %m
0023 
0024 U0=U0*1e6;
0025 emit_x = emit_x*1.0e-9; %convert nm-rad to m-rad
0026 
0027 %cavity related parameters
0028 if isstruct(Trf)
0029     freq = Trf.Frequency;
0030     harm = Trf.HarmNumber;
0031     E0 = Trf.Energy;
0032     alpha = Trf.alpha;
0033     Vrf = Trf.Voltage;
0034     circ = Trf.circum;
0035 else
0036     freq = 352.202e6; %Hz
0037     harm = 992;
0038     E0 = 6.04e9; %eV
0039     alpha = Trf;
0040     Vrf = 9e6;
0041     circ = 844.39;
0042     
0043 end
0044 gamma = E0/0.511e6;
0045 N0 = 0.001/(freq/harm)/e0; %Number of particle per 1mA bunch.
0046 
0047 
0048 %bunch length
0049 phi_s = asin(U0/Vrf);
0050 nus = sqrt(harm*Vrf*alpha*cos(phi_s)/2/pi/E0);
0051 
0052 
0053 sigZ = sigE/nus*harm*alpha/2/pi/freq*cspeed;
0054 
0055 %rf bucket height
0056 delta_max_rf = sqrt(2*U0/pi/alpha/harm/E0)*sqrt( sqrt((Vrf/U0).^2-1) - acos(U0./Vrf));
0057 %---------------------------------
0058 
0059 %beam size around the ring
0060 %[td, tune,chrom] = twissring(THERING,0,1:length(THERING)+1, 'chrom', 1e-5);
0061 td = TD; 
0062 Dx = cat(2, td.Dispersion)';
0063 betxy = cat(1, td.beta);
0064 alfxy = cat(1, td.alpha);
0065 
0066 spos = cat(1,td.SPos);
0067 
0068 sigX = sqrt(betxy(:,1)*emit_x+Dx(:,1).^2*sigE^2);
0069 
0070 sigY = sqrt(betxy(:,2)*emit_x*coupling);
0071 sigXp = sqrt(emit_x*(1+alfxy(:,1).^2)./betxy(:,1)+Dx(:,2).^2*sigE^2);
0072 %--------------------------------
0073 
0074 curH = (Dx(:,1).^2 + (betxy(:,1).*Dx(:,2)+alfxy(:,1).*Dx(:,1)).^2)./betxy(:,1);
0075 
0076 %delta_max_perp = hori_acceptance./sqrt(curH);
0077 deltap = dppPM(:,1);
0078 deltam = dppPM(:,2);
0079 
0080 delta_maxp = min([deltap, ones(size(curH))*delta_max_rf]')';
0081 delta_maxm = min([-deltam, ones(size(curH))*delta_max_rf]')';
0082 
0083 xip = (delta_maxp/gamma.*betxy(:,1)./sigX).^2;
0084 xim = (delta_maxm/gamma.*betxy(:,1)./sigX).^2;
0085 Dvalp = funcD(xip);
0086 Dvalm = funcD(xim);
0087 
0088 ds = diff(spos);
0089 n=1:length(ds);
0090 
0091 avgfacp = sum(Dvalp(n)./sigX(n)./sigY(n)/sigZ./delta_maxp(n).^3.*ds)/circ;
0092 avgfacm = sum(Dvalm(n)./sigX(n)./sigY(n)/sigZ./delta_maxm(n).^3.*ds)/circ;
0093 
0094 lossrate = Ib*N0*r0^2*cspeed/8/gamma^2/pi*(avgfacp+avgfacm)/2.;
0095 tauT = 1/lossrate;
0096 
0097 if 0
0098    figure
0099    h = plot(spos, delta_maxp, spos, -delta_maxm); %delta_max_rf*ones(size(spos)));
0100    set(gca,'fontsize', 16,'xlim',[0,240])
0101    xlabel('s (m)')
0102    ylabel('\delta_{max}')
0103    grid
0104    set(gca,'ylim',[-0.03,0.03]);
0105    
0106 end
0107 
0108 function D=funcD(xi)
0109 %a look-up table
0110 DfunTable = [
0111 %xi                Dfunc
0112 0.000500    0.123802    
0113 0.001000    0.153464    
0114 0.001500    0.172578    
0115 0.002000    0.186757    
0116 0.002500    0.198008    
0117 0.003000    0.207298    
0118 0.003500    0.215179    
0119 0.004000    0.221992    
0120 0.004500    0.227968    
0121 0.005000    0.233269    
0122 0.005500    0.238015    
0123 0.006000    0.242294    
0124 0.006500    0.246176    
0125 0.007000    0.249717    
0126 0.007500    0.252961    
0127 0.008000    0.255944    
0128 0.008500    0.258697    
0129 0.009000    0.261244    
0130 0.009500    0.263607    
0131 0.010000    0.265805    
0132 0.010500    0.267852    
0133 0.011000    0.269763    
0134 0.011500    0.271549    
0135 0.012000    0.273221    
0136 0.012500    0.274788    
0137 0.013000    0.276259    
0138 0.013500    0.277640    
0139 0.014000    0.278938    
0140 0.014500    0.280159    
0141 0.015000    0.281308    
0142 0.015500    0.282391    
0143 0.016000    0.283411    
0144 0.016500    0.284372    
0145 0.017000    0.285278    
0146 0.017500    0.286132    
0147 0.018000    0.286938    
0148 0.018500    0.287698    
0149 0.019000    0.288415    
0150 0.019500    0.289090    
0151 0.020000    0.289727    
0152 0.020500    0.290327    
0153 0.021000    0.290893    
0154 0.021500    0.291425    
0155 0.022000    0.291926    
0156 0.022500    0.292397    
0157 0.023000    0.292840    
0158 0.023500    0.293256    
0159 0.024000    0.293646    
0160 0.024500    0.294011    
0161 0.025000    0.294352    ];
0162 ximin = DfunTable(1,1);
0163 ximax = DfunTable(end,1);
0164 xi(find(xi<ximin)) = ximin;
0165 xi(find(xi>ximax)) = ximax;
0166 D = interp1(DfunTable(:,1), DfunTable(:,2), xi,'linear');
0167

Generated on Thu 24-Aug-2017 18:47:33 by m2html © 2005