source: MML/trunk/machine/SOLEIL/StorageRing/bba_mml/quadgetdata.m @ 4

Last change on this file since 4 was 4, checked in by zhangj, 11 years ago

Initial import--MML version from SOLEIL@2013
File size: 13.4 KB
Line 
1function [X, Y, BPMxFamily, QUADxFamily, BPMyFamily, QUADyFamily, DateMatX, DateMatY QUADxDevList QUADyDevList] = quadgetdata(DirName, PlotFlag)
2%QUADGETDATA - Collect the date from a quadrupole center run.
3%              When all the quadrupole center data files are stored in a directory this function
4%              will go through all the files and pull out some of the important information.
5%
6%  [X, Y, BPMxFamily, QUADxFamily, BPMyFamily, QUADyFamily, XfileDate, YfilesDate] = quadgetdata(DirName)
7%
8%  INPUTS
9%  1. DirName - Directory name to look for quadrupole center files 
10%               [] to browse {Default}
11%  2. PlotFlag - 0 to just get data without plotting results, else, plot results
12% 
13%  OUTPUTS
14%  1. X - Horizontal output matrix (format below)
15%  2. Y - Vertical output matrix (format below)
16%
17%             1&2     3       4        5     6       7         8         9
18%  Output = [BPMDev Center CenterSTD BPMpos DCCT BPMATIndex Quadpos QuadATIndex]
19%
20%  3-6. BPMxFamily, QUADxFamily, BPMyFamily, QUADyFamily - Family name for each center measurement
21%  7-8. XfileDate, YfilesDate - Date string for the file name
22%  9-10 QUADxDevList, QUADyDevList - device list of quad
23%
24%  See Also quadplotall
25
26%
27%%  Written by Greg Portmann
28
29
30if nargin == 0
31    DirName = '';
32end
33
34if nargin < 2
35    PlotFlag = 1;
36end
37
38if isempty(DirName)
39    DirName = [getfamilydata('Directory', 'DataRoot'),'QMS\'];
40    DirName = uigetdir(DirName, 'Select directory where the QMS data are located');
41    if ~ischar(DirName)
42        fprintf('Cancel \n');
43        return;
44    end
45end
46
47StartDir = pwd;
48cd(DirName);
49
50Files = dir;
51
52X = [];
53Y = [];
54BPMxFamily = [];
55QUADxFamily = [];
56BPMyFamily = [];
57QUADyFamily = [];
58QUADxDevList = [];
59QUADyDevList = [];
60DateMatX = [];
61DateMatY = [];
62GoodIndexX = [];
63CorrDeltaX = [];
64QuadDeltaX = [];
65GoodIndexY = [];
66CorrDeltaY = [];
67QuadDeltaY = [];
68OutlierFactorX = [];
69OutlierFactorY = [];
70MaxCorrOrbitVariationX = [];
71MaxQuadOrbitVariationX = [];
72MaxCorrOrbitVariationY = [];
73MaxQuadOrbitVariationY = [];
74MaxTuneY = [];
75MaxTuneX = [];
76MaxHysteresisY = [];
77MaxHysteresisX = [];
78
79% do not forget sorting part !!!
80
81
82for i = 1:length(Files)
83    if exist(Files(i).name) == 2
84        clear QMS DelHCM DelVCM
85        try
86            load(Files(i).name)
87
88            %Files(i).name
89
90            % Old middle layer
91            if exist('DelHCM', 'var')
92                clear q
93                [q.Center, q.CenterSTD, ErrorString] = quadhplt(Files(i).name);
94                q.BPMFamily = 'BPMx';
95                q.BPMDev = BPMnum + [0 1];
96                q.DCCT = NaN;
97                [q.QuadFamily, q.QuadDev] = bpm2quad(q.BPMFamily, q.BPMDev);
98
99                X = [X; q.BPMDev q.Center q.CenterSTD getspos(q.BPMFamily, q.BPMDev) min(q.DCCT) family2atindex(q.BPMFamily, q.BPMDev) getspos(q.QuadFamily, q.QuadDev) family2atindex(q.QuadFamily, q.QuadDev)];
100                BPMxFamily = strvcat(BPMxFamily, q.BPMFamily);
101                QUADxFamily = strvcat(QUADxFamily, q.QuadFamily);
102                DateMatX = [DateMatX; Files(i).date];
103               
104            elseif exist('DelVCM', 'var')
105                clear q
106                [q.Center, q.CenterSTD, ErrorString] = quadvplt(Files(i).name);
107                q.BPMFamily = 'BPMy';
108                q.BPMDev = BPMnum + [0 1];
109                q.DCCT = NaN;
110                [q.QuadFamily, q.QuadDev] = bpm2quad(q.BPMFamily, q.BPMDev);
111
112                Y = [Y; q.BPMDev q.Center q.CenterSTD getspos(q.BPMFamily, q.BPMDev) min(q.DCCT) family2atindex(q.BPMFamily, q.BPMDev) getspos(q.QuadFamily, q.QuadDev) family2atindex(q.QuadFamily, q.QuadDev)];
113                BPMyFamily = strvcat(BPMyFamily, q.BPMFamily);
114                QUADyFamily = strvcat(QUADyFamily, q.QuadFamily);
115                DateMatY = [DateMatY; Files(i).date];
116
117            else
118
119                % New middle layer
120                q = QMS;
121
122                %fprintf('   %d.  %s(%d,%d)\n', i, q.QuadFamily, q.QuadDev);
123
124                if q.QuadPlane == 1 % H-plane
125                    X = [X; q.BPMDev q.Center q.CenterSTD getspos(q.BPMFamily, q.BPMDev) min(q.DCCT) family2atindex(q.BPMFamily, q.BPMDev) getspos(q.QuadFamily, q.QuadDev) family2atindex(q.QuadFamily, q.QuadDev)];
126                    BPMxFamily = strvcat(BPMxFamily, q.BPMFamily);
127                    QUADxFamily = strvcat(QUADxFamily, q.QuadFamily);
128                    QUADxDevList = [QUADxDevList; q.QuadDev];
129                    DateMatX = [DateMatX; Files(i).date];
130                    GoodIndexX = [ GoodIndexX; size(q.GoodIndex,1)];
131                    CorrDeltaX = [ CorrDeltaX; q.CorrDelta];
132                    QuadDeltaX = [ QuadDeltaX; q.QuadDelta];
133                    OutlierFactorX = [OutlierFactorX; q.OutlierFactor];
134                    [BPMelem1, iNotFound] = findrowindex(q.BPMDev, q.BPMDevList);
135                    MaxCorrOrbitVariationX = [MaxCorrOrbitVariationX; max((q.x2(BPMelem1,:)+q.x1(BPMelem1,:))/2) - min((q.x2(BPMelem1,:)+q.x1(BPMelem1,:))/2)];
136                    MaxQuadOrbitVariationX = [MaxQuadOrbitVariationX; max(max(q.x2-q.x1, [], 2))-min(max(q.x2-q.x1, [], 2))];
137                    MaxTuneX = [MaxTuneX; max(abs(q.Tune1 - q.Tune2), [], 2)'];
138                   
139                   
140                else % V-plane
141                    Y = [Y; q.BPMDev q.Center q.CenterSTD getspos(q.BPMFamily, q.BPMDev) min(q.DCCT) family2atindex(q.BPMFamily, q.BPMDev) getspos(q.QuadFamily, q.QuadDev) family2atindex(q.QuadFamily, q.QuadDev)];
142                    BPMyFamily = strvcat(BPMyFamily, q.BPMFamily);
143                    QUADyFamily = strvcat(QUADyFamily, q.QuadFamily);
144                    QUADyDevList = [QUADyDevList; q.QuadDev];
145                    DateMatY = [DateMatY; Files(i).date];
146                    GoodIndexY = [ GoodIndexY; size(q.GoodIndex,1)];
147                    CorrDeltaY = [ CorrDeltaY; q.CorrDelta];
148                    QuadDeltaY = [ QuadDeltaY; q.QuadDelta];
149                    OutlierFactorY= [OutlierFactorY; q.OutlierFactor];
150                    [BPMelem1, iNotFound] = findrowindex(q.BPMDev, q.BPMDevList);
151                    MaxCorrOrbitVariationY = [MaxCorrOrbitVariationY; max((q.y2(BPMelem1,:)+q.y1(BPMelem1,:))/2) - min((q.y2(BPMelem1,:)+q.y1(BPMelem1,:))/2)];
152                    MaxQuadOrbitVariationY = [MaxQuadOrbitVariationY; max(max(q.y2-q.y1, [], 2))-min(max(q.y2-q.y1, [], 2))];
153                    MaxTuneY = [MaxTuneY; max(abs(q.Tune1 - q.Tune2), [], 2)'];
154                end
155                if any(q.DCCT<5)
156                    fprintf('   %s(%d,%d) (%s) shows a beam less than 5 mamps during the experiment!\n', q.BPMFamily, q.BPMDev, Files(i).name);
157                end
158            end
159        catch
160            % not a bba file
161        end
162    end
163end
164
165[DateMatX, iX] = sortrows(DateMatX);
166[DateMatY, iY] = sortrows(DateMatY);
167X = X(iX,:);
168Y = Y(iY,:);
169BPMxFamily = BPMxFamily(iX,:);
170BPMyFamily = BPMyFamily(iY,:);
171QUADxFamily = QUADxFamily(iX,:);
172QUADyFamily = QUADyFamily(iY,:);
173QUADxDevList = QUADxDevList(iX,:);
174QUADyDevList = QUADyDevList(iY,:);
175
176DateMatX = DateMatX(iX,:);
177DateMatY = DateMatY(iY,:);
178GoodIndexX = GoodIndexX(iX,:);
179CorrDeltaX = CorrDeltaX(iX,:);
180QuadDeltaX = QuadDeltaX(iX,:);
181GoodIndexY = GoodIndexY(iY,:);
182CorrDeltaY = CorrDeltaY(iY,:);
183QuadDeltaY = QuadDeltaY(iY,:);
184OutlierFactorX = OutlierFactorX(iX,:);
185OutlierFactorY = OutlierFactorY(iY,:);
186MaxCorrOrbitVariationX = MaxCorrOrbitVariationX(iX,:);
187MaxQuadOrbitVariationX = MaxQuadOrbitVariationX(iX,:);
188MaxCorrOrbitVariationY = MaxCorrOrbitVariationY(iY,:);
189MaxQuadOrbitVariationY = MaxQuadOrbitVariationY(iY,:);
190MaxTuneY = MaxTuneY(iY,:);
191MaxTuneX = MaxTuneX(iX,:);
192% MaxHysteresisY = MaxHysteresisY(iY,:);
193% MaxHysteresisX = MaxHysteresisX(iX,:);
194
195
196
197cd(StartDir);
198
199
200valGoodX = 30;
201SelectIdx = find(GoodIndexX<valGoodX);
202
203fprintf('%3d BPMx with less than %d good BPMs\n', length(SelectIdx), valGoodX);
204for k =1:length(SelectIdx),
205    fprintf('%3s(%2d,%2d): %3d (Center = %+.4f +/- %.4f OutlierFactor=%2d, Dquad= %.3f Dorb=%+.3f Dcorr = %+.3f Dorb =%+.3f |Dnux| = %+.4f |Dnuy| = %+.4f)\n', ...
206        QUADxFamily(SelectIdx(k),:), ...
207        QUADxDevList(SelectIdx(k),:), GoodIndexX(SelectIdx(k)), ...
208        X(SelectIdx(k),3), X(SelectIdx(k),4), OutlierFactorX(SelectIdx(k)), ...
209        QuadDeltaX(SelectIdx(k)), MaxQuadOrbitVariationX(SelectIdx(k)), ...,
210        CorrDeltaX(SelectIdx(k)), MaxCorrOrbitVariationX(SelectIdx(k)), ...
211        MaxTuneX(SelectIdx(k),:));
212end
213
214% BIZARRE 2 April 2011
215%valGoodX = 3e-3;
216SelectIdx = find(X(:,4)>valGoodX);
217
218fprintf('%3d BPMx offsets determined within more than %d mm\n', length(SelectIdx), valGoodX);
219for k =1:length(SelectIdx),
220    fprintf('%3s(%2d,%2d): %3d (Center = %+.4f +/- %.4f OutlierFactor=%2d, Dquad= %.3f Dorb=%+.3f Dcorr = %.3f Dorb =%+.3f |Dnux| = %+.4f |Dnuy| = %+.4f)\n', ...
221        QUADxFamily(SelectIdx(k),:), ...
222        QUADxDevList(SelectIdx(k),:), GoodIndexX(SelectIdx(k)), ...
223        X(SelectIdx(k),3), X(SelectIdx(k),4), OutlierFactorX(SelectIdx(k)), ...
224        QuadDeltaX(SelectIdx(k)), MaxQuadOrbitVariationX(SelectIdx(k)), ...,
225        CorrDeltaX(SelectIdx(k)), MaxCorrOrbitVariationX(SelectIdx(k)), ...
226        MaxTuneX(SelectIdx(k),:));
227end
228
229
230
231valGoodY = 30;
232SelectIdx = find(GoodIndexY<valGoodY);
233
234
235fprintf('%3d BPMy with less than %d good BPMs\n', length(SelectIdx), valGoodY);
236for k =1:length(SelectIdx),
237    fprintf('%3s(%2d,%2d): %3d (Center = %+.4f +/- %.4f OutlierFactor=%2d, Dquad= %.3f Dorb=%+.3f Dcorr = %.3f Dorb =%+.3f |Dnux| = %+.4f |Dnuy| = %+.4f)\n', ...
238        QUADyFamily(SelectIdx(k),:), QUADyDevList(SelectIdx(k),:), GoodIndexY(SelectIdx(k)), ...
239        Y(SelectIdx(k),3), Y(SelectIdx(k),4), OutlierFactorY(SelectIdx(k)), ...
240        QuadDeltaY(SelectIdx(k)), MaxQuadOrbitVariationY(SelectIdx(k)), ...
241        CorrDeltaY(SelectIdx(k)), MaxCorrOrbitVariationY(SelectIdx(k)), ...
242        MaxTuneY(SelectIdx(k),:));
243end
244
245valGoodY = 3e-3;
246SelectIdx = find(Y(:,4)>valGoodY);
247
248fprintf('%3d BPMy offset determined within more than %d mm\n', length(SelectIdx), valGoodY);
249
250for k =1:length(SelectIdx),
251    fprintf('%3s(%2d,%2d): %3d (Center = %+.4f +/- %.4f OutlierFactor=%2d, Dquad= %.3f Dorb=%+.3f Dcorr = %.3f Dorb =%+.3f |Dnux| = %+.4f |Dnuy| = %+.4f)\n', ...
252        QUADyFamily(SelectIdx(k),:), QUADyDevList(SelectIdx(k),:), GoodIndexY(SelectIdx(k)), ...
253        Y(SelectIdx(k),3), Y(SelectIdx(k),4), OutlierFactorY(SelectIdx(k)), ...
254        QuadDeltaY(SelectIdx(k)), MaxQuadOrbitVariationY(SelectIdx(k)), ...
255        CorrDeltaY(SelectIdx(k)), MaxCorrOrbitVariationY(SelectIdx(k)), ...
256        MaxTuneY(SelectIdx(k),:));
257end
258
259if PlotFlag
260
261    % Plot data
262    L = getfamilydata('Circumference');
263
264    figure;
265    subplot(2,1,1);
266    if ~isempty(X)
267        errorbar(X(:,5), X(:,3), X(:,4), '.b');
268    end
269    ylabel('Horizontal [mm]');
270    xaxis([0 L]);
271    title('BPM Offset');
272
273    subplot(2,1,2);
274    if ~isempty(Y)
275        errorbar(Y(:,5), Y(:,3), Y(:,4), '.b');
276    end
277    xlabel('BPM Position [meters]');
278    ylabel('Vertical [mm]');
279    xaxis([0 L]);
280
281    figure;
282    subplot(2,1,1);
283    if ~isempty(X)
284        plot(X(:,5), GoodIndexX,'.');
285    end
286    ylabel('Horizontal Good Index');
287    xaxis([0 L]);
288    title('Good index');
289
290    subplot(2,1,2);
291    if ~isempty(Y)
292        plot(Y(:,5), GoodIndexY,'.');
293    end
294    xlabel('BPM Position [meters]');
295    ylabel('Vertical Good Index');
296    xaxis([0 L]);
297
298    figure;
299    subplot(2,1,1);
300    if ~isempty(X)
301        plot(X(:,5), CorrDeltaX,'.');
302    end
303    ylabel('Horizontal corrector current (A)');
304    xaxis([0 L]);
305    title('Corrector strength');
306
307    subplot(2,1,2);
308    if ~isempty(Y)
309        plot(Y(:,5), CorrDeltaY,'.');
310    end
311    xlabel('BPM Position [meters]');
312    ylabel('Vertical corrector current (A)');
313    xaxis([0 L]);
314
315    figure;
316    subplot(2,1,1);
317    if ~isempty(X)
318        plot(X(:,5), QuadDeltaX,'.');
319    end
320    ylabel('Horizontal corrector current (A)');
321    xaxis([0 L]);
322    title('Quadrupole strength');
323
324    subplot(2,1,2);
325    if ~isempty(Y)
326        plot(Y(:,5), QuadDeltaY,'.');
327    end
328    xlabel('BPM Position [meters]');
329    ylabel('Vertical corrector current (A)');
330    xaxis([0 L]);
331
332    figure;
333    subplot(2,1,1);
334    if ~isempty(X)
335        plot(X(:,5), OutlierFactorX,'.');
336    end
337    ylabel('Horizontal OutlierFactor');
338    xaxis([0 L]);
339    title('OutlierFactor');
340
341    subplot(2,1,2);
342    if ~isempty(Y)
343        plot(Y(:,5), OutlierFactorY,'.');
344    end
345    xlabel('BPM Position [meters]');
346    ylabel('Vertical OutlierFactor');
347    xaxis([0 L]);
348   
349    figure;
350    subplot(2,1,1);
351    if ~isempty(X)
352        plot(X(:,5), MaxTuneX,'.');
353    end
354    ylabel('H-BBA: tuneshifts');
355    xaxis([0 L]);
356    title('Tune shifts');
357    legend('Dnux', 'Dnuz')
358   
359    subplot(2,1,2);
360    if ~isempty(Y)
361        plot(Y(:,5), MaxTuneY,'.');
362    end
363    xlabel('BPM Position [meters]');
364    ylabel('V-BBA: tuneshifts');
365    xaxis([0 L]);
366    legend('Dnux', 'Dnuz')
367
368    figure;
369    subplot(2,1,1);
370    if ~isempty(X)
371        plot(X(:,5), MaxQuadOrbitVariationX,'.');
372    end
373    ylabel('H-BBA: X (mm)');
374    xaxis([0 L]);
375    title('Orbit shift due to quads');
376   
377    subplot(2,1,2);
378    if ~isempty(Y)
379        plot(Y(:,5), MaxQuadOrbitVariationY,'.');
380    end
381    xlabel('BPM Position [meters]');
382    ylabel('V-BBA: : Z (mm)');
383    xaxis([0 L]);
384
385    figure;
386    subplot(2,1,1);
387    if ~isempty(X)
388        plot(X(:,5), MaxCorrOrbitVariationX,'.');
389    end
390    ylabel('H-BBA: X (mm)');
391    xaxis([0 L]);
392    title('Orbit shift due to correctors');
393   
394    subplot(2,1,2);
395    if ~isempty(Y)
396        plot(Y(:,5), MaxCorrOrbitVariationY,'.');
397    end
398    xlabel('BPM Position [meters]');
399    ylabel('V-BBA: : Z (mm)');
400    xaxis([0 L]);
401end
Note: See TracBrowser for help on using the repository browser.