0001 function tauT = calc_TouschekPM(TD,dppPM,Trf,Ib,U0,coupling, sigE, emit_x)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020 e0 = 1.6e-19;
0021 cspeed = 299792458;
0022 r0 = 2.82e-15;
0023
0024 U0=U0*1e6;
0025 emit_x = emit_x*1.0e-9;
0026
0027
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;
0037 harm = 992;
0038 E0 = 6.04e9;
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;
0046
0047
0048
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
0056 delta_max_rf = sqrt(2*U0/pi/alpha/harm/E0)*sqrt( sqrt((Vrf/U0).^2-1) - acos(U0./Vrf));
0057
0058
0059
0060
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
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);
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
0110 DfunTable = [
0111
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