1 | function meascmhysteresis(BPMList, CMFamily, CMList, MaxChange, NSteps) |
---|
2 | %MEASCMHYSTERESIS - Measure corrector magnet hysteresis |
---|
3 | % |
---|
4 | % meascmhysteresis(BPMList, MaxChange, NSteps) |
---|
5 | % |
---|
6 | % INPUTS |
---|
7 | % 1. BPMList - BPM list |
---|
8 | % 2. CMFamily - Corrector family |
---|
9 | % 3. CMList - Corrector device list |
---|
10 | % 4. MaxChange - Maximum vertical change from the starting point |
---|
11 | % 5. NSteps - Number of step from starting point to the maximum |
---|
12 | % |
---|
13 | % NOTES |
---|
14 | % 1. This function looks for a bump coefficient file corresponding to BPM list. |
---|
15 | % For instance, VerticalBumpCoef-7_8-6_1 is the vertical file for [7 8;6 1]. |
---|
16 | % If you don't want to use this file, delete it and a new new will be created. |
---|
17 | % All files are stored in <DataRoot><aperturescan> |
---|
18 | |
---|
19 | % |
---|
20 | % Written by Gregory J. Portmann |
---|
21 | % Adapted by Laurent S. Nadolski |
---|
22 | |
---|
23 | BPMxFamily = gethbpmfamily; |
---|
24 | BPMyFamily = getvbpmfamily; |
---|
25 | HCMFamily = gethcmfamily; |
---|
26 | VCMFamily = getvcmfamily; |
---|
27 | |
---|
28 | |
---|
29 | %%%%%%%%%%%%%%%%%% |
---|
30 | % Input checking % |
---|
31 | %%%%%%%%%%%%%%%%%% |
---|
32 | if nargin < 1 |
---|
33 | BPMList = [7 6; 8 1]; |
---|
34 | end |
---|
35 | if nargin < 2 |
---|
36 | CMFamily = VCMFamily; |
---|
37 | end |
---|
38 | if nargin < 3 |
---|
39 | if size(BPMList,1) == 1 |
---|
40 | CMList = [4 4]; |
---|
41 | else |
---|
42 | CMList = [-4 -3 -2 -1 1 2 3 4]; |
---|
43 | end |
---|
44 | end |
---|
45 | if nargin < 4 |
---|
46 | MaxChange = 1.5; |
---|
47 | end |
---|
48 | if nargin < 5 |
---|
49 | NSteps = 5; % Steps from 0 to MaxChange |
---|
50 | end |
---|
51 | |
---|
52 | |
---|
53 | % Corrector starting points |
---|
54 | HCM0 = getsp(HCMFamily); |
---|
55 | VCM0 = getsp(VCMFamily); |
---|
56 | |
---|
57 | |
---|
58 | % Offset |
---|
59 | Xoffset = getoffset(BPMxFamily, BPMList); |
---|
60 | Yoffset = getoffset(BPMyFamily, BPMList); |
---|
61 | |
---|
62 | % Gain |
---|
63 | Xgain = getphysdata(BPMxFamily, 'Gain', BPMList); |
---|
64 | Ygain = getphysdata(BPMyFamily, 'Gain', BPMList); |
---|
65 | |
---|
66 | |
---|
67 | BPMxIndex = findrowindex(BPMList, getlist(BPMxFamily)); |
---|
68 | BPMyIndex = findrowindex(BPMList, getlist(BPMyFamily)); |
---|
69 | |
---|
70 | BPMxList = getlist(BPMxFamily); |
---|
71 | BPMyList = getlist(BPMyFamily); |
---|
72 | |
---|
73 | FileNumberString = sprintf('-%d_%d', BPMList'); |
---|
74 | DirectoryName = [fullfile(getfamilydata('Directory','DataRoot'),'aperturescan'),filesep]; |
---|
75 | |
---|
76 | |
---|
77 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
78 | % Compute corrector coefficients % |
---|
79 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
80 | BPMWeight = 40; |
---|
81 | SVDIndex = 1e-3; |
---|
82 | |
---|
83 | if strcmpi(CMFamily, HCMFamily) |
---|
84 | HorizontalBumpFile = sprintf('%sHorizontalBumpCoef%s', DirectoryName, FileNumberString); |
---|
85 | if size(BPMList,1) == 1 && isempty(CMList) |
---|
86 | % For single BPM just use one corrector |
---|
87 | % Pick the corrector based on the response matrix |
---|
88 | R = getbpmresp('Struct'); |
---|
89 | [i, iNotFound] = findrowindex(BPMList, R(1,1).Monitor.DeviceList); |
---|
90 | m = R(1,1).Data(i,:); |
---|
91 | [MaxValue, j] = max(abs(m)); |
---|
92 | HCMDeviceList = R(1,1).Actuator.DeviceList(j,:); |
---|
93 | %DeltaHCM = (1/m(j)) ./ Xgain; % Amps/mm (at the BPM BPMList) |
---|
94 | DeltaHCM = (1/max(abs(m))) ./ Xgain; % Amps/mm (maximum at many BPM) |
---|
95 | |
---|
96 | elseif size(CMList,1) == 1 && size(CMList,2) <= 2 |
---|
97 | % For single corrector |
---|
98 | R = getbpmresp('Struct'); |
---|
99 | [i, iNotFound] = findrowindex(BPMList, R(1,1).Monitor.DeviceList); |
---|
100 | [j, jNotFound] = findrowindex(CMList, R(1,1).Actuator.DeviceList); |
---|
101 | m = R(1,1).Data(i,j); |
---|
102 | DeltaHCM = (1/m) ./ Xgain; % Amps/mm |
---|
103 | |
---|
104 | elseif exist([HorizontalBumpFile,'.mat'], 'file') |
---|
105 | fprintf(' Loading horizontal corrector magnet bump coefficients from %s\n', ... |
---|
106 | HorizontalBumpFile); |
---|
107 | load(HorizontalBumpFile); |
---|
108 | |
---|
109 | else |
---|
110 | % Set the hysteresis |
---|
111 | fprintf(' Finding horizontal corrector magnet bump coefficients\n'); |
---|
112 | DeltaX = [.2;.2] ./ Xgain; |
---|
113 | x0 = getam(BPMxFamily, BPMList); |
---|
114 | %[HOCS, RF, HOCS0] = setorbitbump(BPMxFamily, BPMList, DeltaX, HCMFamily, CMList, 1, SVDIndex, BPMWeight, 'Inc', 'NoDisplay'); |
---|
115 | HOCS0 = setorbitbump(BPMxFamily, BPMList, DeltaX, HCMFamily, CMList, 1, ... |
---|
116 | BPMWeight, 'Inc', 'NoDisplay'); |
---|
117 | x1 = getam(BPMxFamily, BPMList); |
---|
118 | |
---|
119 | % Get a clean bump |
---|
120 | figure(1); |
---|
121 | clf reset |
---|
122 | x2 = getam(BPMxFamily, BPMList); |
---|
123 | %[HOCS, RF, HOCS0] = setorbitbump(BPMxFamily, BPMList, DeltaX, HCMFamily, CMList, 1, SVDIndex, BPMWeight, 'Inc', 'NoDisplay'); |
---|
124 | HOCS = setorbitbump(BPMxFamily, BPMList, DeltaX, HCMFamily, CMList, 1, ... |
---|
125 | BPMWeight, 'Inc', 'NoDisplay'); |
---|
126 | x3 = getam(BPMxFamily, BPMList); |
---|
127 | drawnow; |
---|
128 | DeltaHCM = (HOCS.CM.Data - HOCS0.CM.Data)/.2; % amps/mm (real units) |
---|
129 | BPMDeviceList = HOCS.BPM.DeviceList; |
---|
130 | HCMDeviceList = HOCS.CM.DeviceList; |
---|
131 | setsp(HCMFamily, HCM0); |
---|
132 | save(HorizontalBumpFile, 'DeltaHCM', 'HCMDeviceList', 'BPMDeviceList', 'HOCS', 'HOCS0'); |
---|
133 | end |
---|
134 | |
---|
135 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
136 | % Linearity and hystersis check % |
---|
137 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
138 | |
---|
139 | % Starting corrector values |
---|
140 | hcm0 = getsp(HCMFamily, HCMDeviceList); |
---|
141 | |
---|
142 | % Set hysteresis in an upward direction |
---|
143 | setsp(HCMFamily, hcm0 - MaxChange * DeltaHCM, HCMDeviceList, -2); |
---|
144 | setsp(HCMFamily, hcm0, HCMDeviceList); |
---|
145 | |
---|
146 | % Build loop |
---|
147 | Loop = [linspace(0,MaxChange,NSteps+1) linspace(MaxChange-MaxChange/NSteps,-MaxChange,2*NSteps) linspace(-MaxChange+MaxChange/NSteps,0,NSteps)]; |
---|
148 | LoopTotal = [Loop Loop(2:end)]; |
---|
149 | |
---|
150 | i = 0; |
---|
151 | for y = LoopTotal |
---|
152 | i = i + 1; |
---|
153 | setsp(HCMFamily, hcm0 + y * DeltaHCM, HCMDeviceList, -2); |
---|
154 | DCCT(i) = getdcct; |
---|
155 | HCMsp(:,i) = getsp(HCMFamily, HCMDeviceList); |
---|
156 | HCMam(:,i) = getam(HCMFamily, HCMDeviceList); |
---|
157 | %BPMx(:,i) = raw2real(BPMxFamily, getx); |
---|
158 | %BPMy(:,i) = raw2real(BPMyFamily, getz); |
---|
159 | BPMx(:,i) = getx; |
---|
160 | BPMy(:,i) = getz; |
---|
161 | end |
---|
162 | |
---|
163 | % Reset correctors |
---|
164 | setsp(HCMFamily, HCM0); |
---|
165 | |
---|
166 | elseif strcmpi(CMFamily, VCMFamily) |
---|
167 | VerticalBumpFile = sprintf('VerticalBumpCoef%s', FileNumberString); |
---|
168 | if size(BPMList,1) == 1 && isempty(CMList) |
---|
169 | % For single BPM just use one corrector |
---|
170 | % Pick the corrector based on the response matrix |
---|
171 | R = getbpmresp('Struct'); |
---|
172 | [i, iNotFound] = findrowindex(BPMList, R(2,2).Monitor.DeviceList); |
---|
173 | m = R(2,2).Data(i,:); |
---|
174 | [MaxValue, j] = max(abs(m)); |
---|
175 | VCMDeviceList = R(2,2).Actuator.DeviceList(j,:); |
---|
176 | %DeltaVCM = (1/m(j)) ./ Ygain; % Amps/mm (at the BPM BPMList) |
---|
177 | DeltaVCM = (1/max(abs(m))) ./ Ygain; % Amps/mm (maximum at many BPM) |
---|
178 | |
---|
179 | elseif size(CMList,1) == 1 && size(CMList,2) <= 2 |
---|
180 | % For single corrector |
---|
181 | R = getbpmresp('Struct'); |
---|
182 | [i, iNotFound] = findrowindex(BPMList, R(2,2).Monitor.DeviceList); |
---|
183 | [j, jNotFound] = findrowindex(CMList, R(2,2).Actuator.DeviceList); |
---|
184 | m = R(2,2).Data(i,j); |
---|
185 | DeltaVCM = (1/m) ./ Ygain; % Amps/mm |
---|
186 | VCMDeviceList = CMList; |
---|
187 | |
---|
188 | elseif exist([VerticalBumpFile,'.mat'], 'file') |
---|
189 | fprintf(' Loading vertical corrector magnet bump coefficients from %s\n', VerticalBumpFile); |
---|
190 | load(VerticalBumpFile); |
---|
191 | else |
---|
192 | % Set the hysteresis |
---|
193 | fprintf(' Finding vertical corrector magnet bump coefficients\n'); |
---|
194 | DeltaY = [.2;.2] ./ Ygain; |
---|
195 | y0 = getam(BPMyFamily, BPMList); |
---|
196 | %[VOCS, RF, VOCS0] = setorbitbump(BPMyFamily, BPMList, DeltaY, VCMFamily, CMList, 1, SVDIndex, BPMWeight, 'Inc', 'NoDisplay'); |
---|
197 | % in one iteration |
---|
198 | VOCS0 = setorbitbump(BPMyFamily, BPMList, DeltaY, VCMFamily, CMList, ... |
---|
199 | 1, BPMWeight, 'Inc', 'NoDisplay'); |
---|
200 | y1 = getam(BPMyFamily, BPMList); |
---|
201 | |
---|
202 | % Get a clean bump |
---|
203 | %vcm0= getsp(VCMFamily, VOCS.CM.DeviceList); |
---|
204 | figure(1); |
---|
205 | clf reset |
---|
206 | y2 = getam(BPMyFamily, BPMList); |
---|
207 | VOCS = setorbitbump(BPMyFamily, BPMList, DeltaY, VCMFamily, CMList, ... |
---|
208 | 5, BPMWeight, 'Inc', 'NoDisplay'); |
---|
209 | % [VOCS, RF, VOCS0] = setorbitbump(BPMyFamily, BPMList, DeltaY, VCMFamily, CMList, 5, SVDIndex, BPMWeight, 'Inc', 'NoDisplay'); |
---|
210 | y3 = getam(BPMyFamily, BPMList); |
---|
211 | drawnow; |
---|
212 | %vcm1= getsp(VCMFamily, VOCS.CM.DeviceList); |
---|
213 | DeltaVCM = (VOCS.CM.Data - VOCS0.CM.Data) / .2; % amps/mm (real units) |
---|
214 | BPMDeviceList = VOCS.BPM.DeviceList; |
---|
215 | VCMDeviceList = VOCS.CM.DeviceList; |
---|
216 | setsp(VCMFamily, VCM0); |
---|
217 | save(VerticalBumpFile, 'DeltaVCM', 'VCMDeviceList', 'BPMDeviceList', 'VOCS', 'VOCS0'); |
---|
218 | end |
---|
219 | |
---|
220 | |
---|
221 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
222 | % Linearity and hystersis check % |
---|
223 | %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% |
---|
224 | |
---|
225 | % Starting corrector values |
---|
226 | vcm0 = getsp(VCMFamily, VCMDeviceList); |
---|
227 | |
---|
228 | % Set hysteresis in an upward direction |
---|
229 | %setsp(VCMFamily, vcm0 - MaxChange * DeltaVCM, VCMDeviceList, -2); |
---|
230 | %setsp(VCMFamily, vcm0, VCMDeviceList); |
---|
231 | |
---|
232 | % Build loop |
---|
233 | Loop = [linspace(0,MaxChange,NSteps+1) ... |
---|
234 | linspace(MaxChange-MaxChange/NSteps,-MaxChange,2*NSteps) ... |
---|
235 | linspace(-MaxChange+MaxChange/NSteps,0,NSteps)]; |
---|
236 | LoopTotal = [Loop Loop(2:end)]; |
---|
237 | |
---|
238 | i = 0; |
---|
239 | for y = LoopTotal |
---|
240 | i = i + 1; |
---|
241 | setsp(VCMFamily, vcm0 + y * DeltaVCM, VCMDeviceList, -2); |
---|
242 | DCCT(i) = getdcct; |
---|
243 | VCMsp(:,i) = getsp(VCMFamily, VCMDeviceList); |
---|
244 | VCMam(:,i) = getam(VCMFamily, VCMDeviceList); |
---|
245 | %BPMx(:,i) = raw2real(BPMxFamily, getx); |
---|
246 | %BPMy(:,i) = raw2real(BPMyFamily, getz); |
---|
247 | BPMx(:,i) = getx; |
---|
248 | BPMy(:,i) = getz; |
---|
249 | |
---|
250 | end |
---|
251 | |
---|
252 | % Reset correctors |
---|
253 | setsp(VCMFamily, VCM0); |
---|
254 | |
---|
255 | else |
---|
256 | error('Corrector magnet family unknown'); |
---|
257 | end |
---|
258 | |
---|
259 | |
---|
260 | % Save data |
---|
261 | HysteresisFilename = sprintf('Hysteresis%s', FileNumberString); |
---|
262 | HysteresisFilename = appendtimestamp(HysteresisFilename); |
---|
263 | save(HysteresisFilename); |
---|
264 | fprintf(' Hysteresis data saved to %s.mat\n', HysteresisFilename); |
---|
265 | |
---|
266 | |
---|
267 | % Plot results |
---|
268 | clf reset |
---|
269 | LineType = '.-b'; |
---|
270 | if size(BPMList,1) == 1 |
---|
271 | |
---|
272 | subplot(2,1,1); |
---|
273 | plot(LoopTotal, BPMy(BPMyIndex(1),:),LineType); |
---|
274 | xlabel(sprintf('Goal Orbit at %s(%d,%d) [mm]', BPMyFamily, BPMList(1,:))); |
---|
275 | ylabel(sprintf('%s(%d,%d) [mm]', BPMyFamily, BPMList(1,:))); |
---|
276 | |
---|
277 | subplot(2,1,2); |
---|
278 | plot(LoopTotal, BPMy(BPMyIndex(1),:)-LoopTotal,LineType); |
---|
279 | xlabel(sprintf('Goal Orbit at %s(%d,%d) [mm]', BPMyFamily, BPMList(1,:))); |
---|
280 | ylabel(sprintf('%s(%d,%d) Error [mm]', BPMyFamily, BPMList(1,:))); |
---|
281 | title('Hysteresis Loops'); |
---|
282 | |
---|
283 | else |
---|
284 | |
---|
285 | subplot(2,2,1); |
---|
286 | plot(LoopTotal, BPMy(BPMyIndex(1),:),LineType); |
---|
287 | xlabel(sprintf('Goal Orbit at %s(%d,%d) [mm]', BPMyFamily, BPMList(1,:))); |
---|
288 | ylabel(sprintf('%s(%d,%d) [mm]', BPMyFamily, BPMList(1,:))); |
---|
289 | grid on; |
---|
290 | |
---|
291 | subplot(2,2,3); |
---|
292 | plot(LoopTotal, BPMy(BPMyIndex(1),:)-LoopTotal,LineType); |
---|
293 | xlabel(sprintf('Goal Orbit at %s(%d,%d) [mm]', BPMyFamily, BPMList(1,:))); |
---|
294 | ylabel(sprintf('%s(%d,%d) Error [mm]', BPMyFamily, BPMList(1,:))); |
---|
295 | grid on; |
---|
296 | |
---|
297 | subplot(2,2,2); |
---|
298 | plot(LoopTotal, BPMy(BPMyIndex(2),:),LineType); |
---|
299 | xlabel(sprintf('Goal Orbit at %s(%d,%d) [mm]', BPMyFamily, BPMList(2,:))); |
---|
300 | ylabel(sprintf('%s(%d,%d) [mm]', BPMyFamily, BPMList(2,:))); |
---|
301 | grid on; |
---|
302 | |
---|
303 | subplot(2,2,4); |
---|
304 | plot(LoopTotal, BPMy(BPMyIndex(2),:)-LoopTotal,LineType); |
---|
305 | xlabel(sprintf('Goal Orbit at %s(%d,%d) [mm]', BPMyFamily, BPMList(2,:))); |
---|
306 | ylabel(sprintf('%s(%d,%d) Error [mm]', BPMyFamily, BPMList(2,:))); |
---|
307 | grid on; |
---|
308 | |
---|
309 | h = addlabel(.5, 1, 'Hysteresis Loops'); |
---|
310 | set(h, 'Fontsize', 12); |
---|
311 | orient tall |
---|
312 | |
---|
313 | end |
---|
314 | |
---|