1 | function AO = measbeta(varargin) |
---|
2 | %MEASBETA - Measure betatron functions |
---|
3 | % |
---|
4 | % INPUTS |
---|
5 | % 1. Quadrupole family name {Default : All} |
---|
6 | % Optional |
---|
7 | % 'Archive', 'Display' |
---|
8 | % Optional override of the mode: |
---|
9 | % 'Online' - Set/Get data online |
---|
10 | % 'Model' - Get the model chromaticity directly from AT (uses modelchro, DeltaRF is ignored) |
---|
11 | % 'Simulator' - Set/Get data on the simulated accelerator using AT (ie, same commands as 'Online') |
---|
12 | % |
---|
13 | % OUPUTS |
---|
14 | % 1. betax - Horizontal beta functions |
---|
15 | % 2. betaz - Vertical beta functions |
---|
16 | % |
---|
17 | % ALGORITHM |
---|
18 | % betax = 4*pi*Dtunex/D(KL) |
---|
19 | % betaz = -4*pi*Dtunez/D(KL) |
---|
20 | % |
---|
21 | % See Also measbeta |
---|
22 | |
---|
23 | % |
---|
24 | % Written by Laurent S. Nadolski |
---|
25 | |
---|
26 | DisplayFlag = 1; |
---|
27 | ArchiveFlag = 1; |
---|
28 | FileName = ''; |
---|
29 | ModeFlag = ''; % model, online, manual, or '' for default mode |
---|
30 | waittime = 10; %seconds taken into account for simulator and Online |
---|
31 | OutputFlag = 1; |
---|
32 | |
---|
33 | for i = length(varargin):-1:1 |
---|
34 | if isstruct(varargin{i}) |
---|
35 | % Ignore structures |
---|
36 | elseif iscell(varargin{i}) |
---|
37 | % Ignore cells |
---|
38 | elseif strcmpi(varargin{i},'Display') |
---|
39 | DisplayFlag = 1; |
---|
40 | varargin(i) = []; |
---|
41 | elseif strcmpi(varargin{i},'NoDisplay') |
---|
42 | DisplayFlag = O; |
---|
43 | varargin(i) = []; |
---|
44 | elseif strcmpi(varargin{i},'NoArchive') |
---|
45 | ArchiveFlag = 0; |
---|
46 | varargin(i) = []; |
---|
47 | elseif strcmpi(varargin{i},'Archive') |
---|
48 | ArchiveFlag = 1; |
---|
49 | varargin(i) = []; |
---|
50 | elseif strcmpi(varargin{i},'Simulator') || strcmpi(varargin{i},'Model') ... |
---|
51 | || strcmpi(varargin{i},'Online') || strcmpi(varargin{i},'Manual') |
---|
52 | ModeFlag = varargin{i}; |
---|
53 | varargin(i) = []; |
---|
54 | end |
---|
55 | end |
---|
56 | |
---|
57 | if strcmpi(ModeFlag,'Model') |
---|
58 | waittime = -1; |
---|
59 | OutputFlag = 0; |
---|
60 | end |
---|
61 | |
---|
62 | |
---|
63 | % Input parsing |
---|
64 | if isempty(varargin) |
---|
65 | QuadFam = findmemberof('QUAD'); |
---|
66 | elseif ischar(varargin{1}) |
---|
67 | QuadFam = {varargin{:}}; |
---|
68 | else |
---|
69 | QuadFam = varargin{:} |
---|
70 | end |
---|
71 | |
---|
72 | if ArchiveFlag |
---|
73 | if isempty(FileName) |
---|
74 | FileName = appendtimestamp(getfamilydata('Default', 'QUADArchiveFile')); |
---|
75 | DirectoryName = getfamilydata('Directory', 'QUAD'); |
---|
76 | if isempty(DirectoryName) |
---|
77 | DirectoryName = fullfile(getfamilydata('Directory','DataRoot'), 'Response', 'BPM'); |
---|
78 | else |
---|
79 | % Make sure default directory exists |
---|
80 | DirStart = pwd; |
---|
81 | [DirectoryName, ErrorFlag] = gotodirectory(DirectoryName); |
---|
82 | cd(DirStart); |
---|
83 | end |
---|
84 | [FileName, DirectoryName] = uiputfile('*.mat', 'Select a Quad File ("Save" starts measurement)', [DirectoryName FileName]); |
---|
85 | if FileName == 0 |
---|
86 | ArchiveFlag = 0; |
---|
87 | disp(' Quadrupole betatron measurement canceled.'); |
---|
88 | return |
---|
89 | end |
---|
90 | FileName = fullfile(DirectoryName, FileName); |
---|
91 | elseif FileName == -1 |
---|
92 | FileName = appendtimestamp(getfamilydata('Default', 'QUADArchiveFile')); |
---|
93 | DirectoryName = getfamilydata('Directory', 'QUAD'); |
---|
94 | FileName = fullfile(DirectoryName, FileName); |
---|
95 | end |
---|
96 | end |
---|
97 | |
---|
98 | % Starting time |
---|
99 | t0 = clock; |
---|
100 | |
---|
101 | nu_start = gettune(ModeFlag); |
---|
102 | |
---|
103 | for k1 = 1:length(QuadFam), |
---|
104 | |
---|
105 | if ~isfamily(QuadFam{k1}) |
---|
106 | error('%s is not a valid Family \n', QuadFam{k1}); |
---|
107 | return; |
---|
108 | end |
---|
109 | |
---|
110 | DeviceList = family2dev(QuadFam{k1}); |
---|
111 | |
---|
112 | % initialize data to zeros |
---|
113 | beta = zeros(length(DeviceList),2); |
---|
114 | beta_vrai = beta; |
---|
115 | tune0 = beta; |
---|
116 | tune1 = beta; |
---|
117 | tune2 = beta; |
---|
118 | dtune = beta; |
---|
119 | |
---|
120 | k3 = 0; |
---|
121 | |
---|
122 | for k2 = 1:length(DeviceList), |
---|
123 | Ic = getam(QuadFam{k1}, DeviceList(k2,:), ModeFlag); |
---|
124 | K = hw2physics(QuadFam{k1}, 'Setpoint', Ic, DeviceList(k2,:)); |
---|
125 | |
---|
126 | if OutputFlag |
---|
127 | fprintf('Measuring Family %s [%d %d] actual current %f A : ... \n', ... |
---|
128 | QuadFam{k1}, DeviceList(k2,:),Ic) |
---|
129 | end |
---|
130 | |
---|
131 | k3 = k3 + 1; |
---|
132 | tune0(k3,:) = gettune(ModeFlag); % Starting time |
---|
133 | |
---|
134 | DeltaI = getfamilydata(QuadFam{k1},'Setpoint','DeltaKBeta')*1.; % Amp |
---|
135 | |
---|
136 | if OutputFlag |
---|
137 | fprintf('Current increment of %d A\n', DeltaI) |
---|
138 | end |
---|
139 | |
---|
140 | stepsp(QuadFam{k1}, DeltaI, DeviceList(k2,:), ModeFlag); % Step value |
---|
141 | sleep(waittime) % wait for quad reaching new setpoint value |
---|
142 | |
---|
143 | tune1(k3,:) = gettune(ModeFlag); % get new tunes |
---|
144 | |
---|
145 | if OutputFlag |
---|
146 | tune1 |
---|
147 | fprintf('Current increment of %d A\n', -2*DeltaI) |
---|
148 | end |
---|
149 | |
---|
150 | stepsp(QuadFam{k1}, -2*DeltaI, DeviceList(k2,:), ModeFlag); % go back to initial values |
---|
151 | sleep(waittime) % wait for quad reaching new setpoint value |
---|
152 | |
---|
153 | tune2(k3,:) = gettune(ModeFlag); % get new tunes |
---|
154 | |
---|
155 | if OutputFlag |
---|
156 | tune2 |
---|
157 | end |
---|
158 | |
---|
159 | if OutputFlag |
---|
160 | fprintf('Current increment of %d A\n', DeltaI) |
---|
161 | end |
---|
162 | |
---|
163 | stepsp(QuadFam{k1}, DeltaI, DeviceList(k2,:), ModeFlag); % go back to initial values |
---|
164 | sleep(waittime) % wait for quad reaching new setpoint value |
---|
165 | |
---|
166 | %% computation part |
---|
167 | |
---|
168 | dtune(k3,:) = tune1(k3,:) - tune2(k3,:); |
---|
169 | |
---|
170 | Leff = getleff(QuadFam{k1}, DeviceList(k2,:)); % Get effective length |
---|
171 | %KL = hw2physics(QuadFam{k1}, 'Setpoint', DeltaK, DeviceList(k2,:))*Leff; |
---|
172 | DeltaKL = 2*DeltaI/Ic*K*Leff; |
---|
173 | |
---|
174 | K1 = hw2physics(QuadFam{k1}, 'Setpoint', Ic+DeltaI, DeviceList(k2,:)); |
---|
175 | K2 = hw2physics(QuadFam{k1}, 'Setpoint', Ic-DeltaI, DeviceList(k2,:)); |
---|
176 | DeltaKL_vrai = (K1-K2)*Leff; |
---|
177 | |
---|
178 | beta(k3,:) = 4*pi*dtune(k3,:)./DeltaKL.*[1 -1]; |
---|
179 | beta_vrai(k3,:) = 4*pi*dtune(k3,:)./DeltaKL_vrai.*[1 -1]; |
---|
180 | |
---|
181 | if OutputFlag |
---|
182 | dtune |
---|
183 | beta |
---|
184 | beta_vrai |
---|
185 | end |
---|
186 | end |
---|
187 | |
---|
188 | % structure to be saved |
---|
189 | AO.FamilyName.(QuadFam{k1}).beta = beta; |
---|
190 | AO.FamilyName.(QuadFam{k1}).beta_vrai = beta_vrai; |
---|
191 | AO.FamilyName.(QuadFam{k1}).dtune = dtune; |
---|
192 | AO.FamilyName.(QuadFam{k1}).tune0 = tune0; |
---|
193 | AO.FamilyName.(QuadFam{k1}).tune1 = tune1; |
---|
194 | AO.FamilyName.(QuadFam{k1}).tune2 = tune2; |
---|
195 | AO.FamilyName.(QuadFam{k1}).deltaI = DeltaI; |
---|
196 | AO.FamilyName.(QuadFam{k1}).DeviceList = DeviceList; |
---|
197 | %AO.FamilyName.(QuadFam{k1}).Position = getspos(QuadFam{k1},DeviceList); |
---|
198 | end |
---|
199 | |
---|
200 | AO.CreatedBy = 'measbeta'; |
---|
201 | AO.GeV = getenergy; |
---|
202 | AO.t = t0; |
---|
203 | AO.tout = etime(clock,t0); |
---|
204 | AO.TimeStamp = datestr(clock); |
---|
205 | |
---|
206 | if ArchiveFlag |
---|
207 | save(FileName,'AO'); |
---|
208 | fprintf('Data save in filename %s \n', FileName); |
---|
209 | end |
---|
210 | |
---|
211 | %% tune variation during measurement |
---|
212 | nu_end = gettune(ModeFlag); |
---|
213 | |
---|
214 | fprintf('Tunes before mesurement nux = %4.4f nuz = %4.4f \n', nu_start); |
---|
215 | fprintf('Tunes after mesurement nux = %4.4f nuz = %4.4f \n', nu_end); |
---|
216 | |
---|
217 | %% raw statistics on beta measurement |
---|
218 | dbxobx = (max(beta_vrai(:,1)-min(beta_vrai(:,1))))./min(beta_vrai(:,1))*100; |
---|
219 | dbzobz = (max(beta_vrai(:,2)-min(beta_vrai(:,2))))./min(beta_vrai(:,2))*100; |
---|
220 | fprintf('maximum betabeat dbxobx = %4.1f %% dbzobz = %4.1f %% \n', dbxobx, dbzobz); |
---|
221 | |
---|
222 | rmsbx = std((beta_vrai(:,1)-mean(beta_vrai(:,1)))./mean(beta_vrai(:,1)))*100; |
---|
223 | rmsbz = std((beta_vrai(:,2)-mean(beta_vrai(:,2)))./mean(beta_vrai(:,2)))*100; |
---|
224 | fprintf('rms betabeat bx = %4.1f %% rms bz = %4.1f %% rms \n', rmsbx, rmsbz); |
---|
225 | |
---|
226 | if DisplayFlag |
---|
227 | plotmeasbeta(AO); |
---|
228 | end |
---|