0001 function locoplotrms(FileName, IterationNumber, PlotType)
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018
0019
0020
0021
0022
0023
0024
0025
0026
0027
0028
0029 if ~nargin==3
0030 error('Requires 3 inputs (see help locoplotrms).');
0031 end
0032
0033 if isempty(FileName)
0034 return;
0035 end
0036
0037 if iscell(FileName)
0038 BPMData = FileName{1};
0039 CMData = FileName{2};
0040 LocoMeasData = FileName{3};
0041 LocoModel = FileName{4};
0042 FitParameters = FileName{5};
0043 LocoFlags = FileName{6};
0044 RINGData = FileName{7};
0045 elseif isstr(FileName)
0046 try
0047 load(FileName);
0048 catch
0049 fprintf(' LOCOPLOT: File does not exist or is not a *.mat file type.\n'); return;
0050 end
0051 else
0052 error('Input problem');
0053 end
0054
0055 if length(BPMData) > 1
0056 IterationNumber = IterationNumber + 1;
0057 if IterationNumber > length(BPMData)
0058 fprintf(' LOCOPLOTRMS: The data file only has %d iterations. Hence, the input InterationNumber cannot be %d.\n', length(BPMData)-1, IterationNumber-1);
0059 return
0060 end
0061
0062 BPMData = BPMData(IterationNumber);
0063 CMData = CMData(IterationNumber);
0064 LocoModel = LocoModel(IterationNumber);
0065 FitParameters = FitParameters(IterationNumber);
0066 LocoFlags = LocoFlags(IterationNumber);
0067 end
0068
0069
0070 if isstr(BPMData)
0071 FileName = BPMData;
0072 IterationNumber = CMData;
0073 PlotType = LocoMeasData;
0074
0075 if isempty(FileName)
0076 return
0077 end
0078 try
0079 load(FileName);
0080 catch
0081 fprintf(' File does not exist or is not a *.mat file type.\n');
0082 cla;
0083 return
0084 end
0085
0086 IterationNumber = IterationNumber + 1;
0087 BPMData = BPMData(IterationNumber);
0088 CMData = CMData(IterationNumber);
0089 LocoModel = LocoModel(IterationNumber);
0090 LocoFlags = LocoFlags(IterationNumber);
0091 else
0092
0093 BPMData = BPMData(end);
0094 CMData = CMData(end);
0095 LocoModel = LocoModel(end);
0096 LocoFlags = LocoFlags(end);
0097 end
0098
0099 if isempty(LocoModel.M)
0100 fprintf(' No model available.\n');
0101 return
0102 end
0103
0104 Mmodel = LocoModel.M;
0105 Outliers = LocoModel.OutlierIndex;
0106 ChiSquare = LocoModel.ChiSquare;
0107 Mmeas = LocoMeasData.M;
0108
0109
0110 BPMstd = LocoMeasData.BPMSTD([BPMData.HBPMGoodDataIndex length(BPMData.HBPMIndex)+BPMData.VBPMGoodDataIndex]);
0111 Mstd = BPMstd * ones(1,size(Mmodel,2));
0112 Mmeas = Mmeas([BPMData.HBPMGoodDataIndex length(BPMData.HBPMIndex)+BPMData.VBPMGoodDataIndex], [CMData.HCMGoodDataIndex length(CMData.HCMIndex)+CMData.VCMGoodDataIndex]);
0113
0114 NHBPM = length(BPMData.HBPMGoodDataIndex);
0115 NVBPM = length(BPMData.VBPMGoodDataIndex);
0116 NBPM = NHBPM + NVBPM;
0117 NHCM = length(CMData.HCMGoodDataIndex);
0118 NVCM = length(CMData.VCMGoodDataIndex);
0119
0120
0121
0122 if strcmp(lower(LocoFlags.ClosedOrbitType), 'fixedmomentum')
0123
0124 HCMEnergyShift = CMData.HCMEnergyShift(CMData.HCMGoodDataIndex);
0125 VCMEnergyShift = CMData.VCMEnergyShift(CMData.VCMGoodDataIndex);
0126
0127
0128 for j = 1:length(FitParameters.Params)
0129 RINGData = locosetlatticeparam(RINGData, FitParameters.Params{j}, FitParameters.Values(j));
0130 end
0131 AlphaMCF = locomcf(RINGData);
0132 EtaXmcf = -AlphaMCF * LocoMeasData.RF * LocoMeasData.Eta(BPMData.HBPMGoodDataIndex) / LocoMeasData.DeltaRF;
0133 EtaYmcf = -AlphaMCF * LocoMeasData.RF * LocoMeasData.Eta(length(BPMData.HBPMIndex)+BPMData.VBPMGoodDataIndex) / LocoMeasData.DeltaRF;
0134
0135 for i = 1:length(HCMEnergyShift)
0136 Mmodel(:,i) = Mmodel(:,i) + HCMEnergyShift(i) * [EtaXmcf; EtaYmcf];
0137 end
0138
0139 for i = 1:length(VCMEnergyShift)
0140 Mmodel(:,NHCM+i) = Mmodel(:,NHCM+i) + VCMEnergyShift(i) * [EtaXmcf; EtaYmcf];
0141 end
0142 end
0143
0144 M = Mmeas(:)-Mmodel(:);
0145
0146
0147 M = reshape((M ./ Mstd(:)) .^ 2, NHBPM+NVBPM, NHCM+NVCM);
0148
0149 M11 = M(1:NHBPM , 1:NHCM );
0150 M12 = M(1:NHBPM , NHCM+1:NHCM+NVCM);
0151 M21 = M(NHBPM+1:NBPM, 1:NHCM );
0152 M22 = M(NHBPM+1:NBPM, NHCM+1:NHCM+NVCM);
0153
0154 BPMRMSxx = sqrt(sum(M11') / NHCM);
0155 BPMRMSyy = sqrt(sum(M22') / NVCM);
0156
0157 HCMRMS = sqrt(sum(M11) / NHBPM);
0158 VCMRMS = sqrt(sum(M22) / NVBPM);
0159
0160
0161 switch PlotType
0162 case 1
0163 plot(BPMData.HBPMGoodDataIndex, BPMRMSxx,'b');
0164 title('Horizontal Response Matrix RMS by Row');
0165 ylabel('\fontsize{10}\surd\Sigma(Mmeas-Mmodel)^2/\sigma^2(BPMx)/NHCM');
0166 xlabel('Horizontal BPM Number');
0167 axis tight
0168
0169 case 2
0170 plot(CMData.HCMGoodDataIndex, HCMRMS,'b');
0171 title('Horizontal Response Matrix RMS by Column');
0172 ylabel('\fontsize{10}\surd\Sigma(Mmeas-Mmodel)^2/\sigma^2(BPMx)/NHBPM');
0173 xlabel('Horizontal Corrector Number');
0174 axis tight
0175
0176 case 3
0177 plot(BPMData.VBPMGoodDataIndex, BPMRMSyy,'b');
0178 title('Vertical Response Matrix RMS by Row');
0179 ylabel('\fontsize{10}\surd\Sigma(Mmeas-Mmodel)^2/\sigma^2(BPMy)/NVCM');
0180 xlabel('Vertical BPM Number');
0181 axis tight
0182
0183 case 4
0184 plot(CMData.VCMGoodDataIndex, VCMRMS,'b');
0185 title('Vertical Response Matrix RMS by Column');
0186 ylabel('\fontsize{10}\surd\Sigma(Mmeas-Mmodel)^2/\sigma^2(BPMy)/NVBPM');
0187 xlabel('Vertical Corrector Number');
0188 axis tight
0189 end
0190
0191
0192
0193 Outliers = LocoModel.OutlierIndex;
0194 if ~isempty(Outliers)
0195 M = M(:);
0196 M(Outliers) = 0;
0197
0198
0199 M = reshape(M, NHBPM+NVBPM, NHCM+NVCM);
0200
0201 M11 = M(1:NHBPM , 1:NHCM );
0202 M12 = M(1:NHBPM , NHCM+1:NHCM+NVCM);
0203 M21 = M(NHBPM+1:NBPM, 1:NHCM );
0204 M22 = M(NHBPM+1:NBPM, NHCM+1:NHCM+NVCM);
0205
0206 BPMRMSxx = sqrt(sum(M11') / NHCM);
0207 BPMRMSyy = sqrt(sum(M22') / NVCM);
0208
0209 HCMRMS = sqrt(sum(M11) / NHBPM);
0210 VCMRMS = sqrt(sum(M22) / NVBPM);
0211
0212 switch PlotType
0213 case 1
0214 hold on;
0215 plot(BPMData.HBPMGoodDataIndex, BPMRMSxx, '--r');
0216 hold off;
0217 axis tight
0218 legend('Full Data Set','Outliers Removed',0);
0219 case 2
0220 hold on;
0221 plot(CMData.HCMGoodDataIndex, HCMRMS, '--r');
0222 hold off;
0223 axis tight
0224 legend('Full Data Set','Outliers Removed',0);
0225 case 3
0226 hold on;
0227 plot(BPMData.VBPMGoodDataIndex, BPMRMSyy, '--r');
0228 hold off;
0229 axis tight
0230 legend('Full Data Set','Outliers Removed',0);
0231 case 4
0232 hold on;
0233 plot(CMData.VCMGoodDataIndex, VCMRMS, '--r');
0234 hold off;
0235 axis tight
0236 legend('Full Data Set','Outliers Removed',0);
0237 end
0238 end