From a349ae17f342ea159446aec91f9e615f9c9b3655 Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Thu, 2 Mar 2023 17:11:47 +0100 Subject: [PATCH 01/10] legend on multiple columns in `CompareOptics.m` --- MADX-optics/CompareOptics.m | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/MADX-optics/CompareOptics.m b/MADX-optics/CompareOptics.m index 1b55e32..c37aa29 100644 --- a/MADX-optics/CompareOptics.m +++ b/MADX-optics/CompareOptics.m @@ -73,12 +73,12 @@ function CompareOptics(optics,labels,geometry,whats,myTitle,emig,sigdpp) % - H plane axs{2,1}=subplot(3,nCols,2); PlotOptics(optics(:,:),usrWhats(1),emigUsr,sigdppUsr); - legend(labels,"Location","best"); + legend(labels,"Location","best","NumColumns",ceil(nOpts/10.)); grid on; % - V plane axs{3,1}=subplot(3,nCols,3); PlotOptics(optics(:,:),usrWhats(2),emigUsr,sigdppUsr); - legend(labels,"Location","best"); + legend(labels,"Location","best","NumColumns",ceil(nOpts/10.)); grid on; % - miscellanea linkaxes([axs{:,1}],'x'); From 151bb585ab667236d80626f32a98ae31fe256d9c Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Thu, 2 Mar 2023 17:12:29 +0100 Subject: [PATCH 02/10] better palette in `PlotOptics.m` --- MADX-optics/PlotOptics.m | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/MADX-optics/PlotOptics.m b/MADX-optics/PlotOptics.m index ecc3ecc..359f109 100644 --- a/MADX-optics/PlotOptics.m +++ b/MADX-optics/PlotOptics.m @@ -6,7 +6,7 @@ function PlotOptics(tfsTables,what,emig,sigdpp,avedpp,refTfsTable,whatRef) nTables=size(tfsTables,1); if ( nTables>1 ) - colormap(parula(nTables)); + colormap(jet(nTables)); end for jj=1:nTables if ( exist('refTfsTable','var') ) From 5cae76925237ead9e2981c1d4720386541c3df9e Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Fri, 28 Apr 2023 12:04:45 +0200 Subject: [PATCH 03/10] minor bug fix --- general/Plot2DHistograms.m | 12 ++++++++---- 1 file changed, 8 insertions(+), 4 deletions(-) diff --git a/general/Plot2DHistograms.m b/general/Plot2DHistograms.m index f33d735..3b01836 100644 --- a/general/Plot2DHistograms.m +++ b/general/Plot2DHistograms.m @@ -11,8 +11,10 @@ function Plot2DHistograms(showMe,showMe1DX,showMe1DY,xShow,yShow,xShowLabel,ySho plot(showMe(:,1),showMe(:,2),"k."); end if ( exist("contours","var") ) - for ii=1:size(contours,3) - hold on; plot(contours(:,1,ii),contours(:,2,ii),".-"); + if ( ~ismissing(contours) ) + for ii=1:size(contours,3) + hold on; plot(contours(:,1,ii),contours(:,2,ii),".-"); + end end end % additionals @@ -24,8 +26,10 @@ function Plot2DHistograms(showMe,showMe1DX,showMe1DY,xShow,yShow,xShowLabel,ySho if ( lSquared ) myMin=min(min(xShow),min(yShow)); myMax=max(max(xShow),max(yShow)); if ( exist("contours","var") ) - myMin=min(myMin,min(contours,[],"all")); - myMax=max(myMax,max(contours,[],"all")); + if ( ~ismissing(contours) ) + myMin=min(myMin,min(contours,[],"all")); + myMax=max(myMax,max(contours,[],"all")); + end end xlim([myMin myMax]); ylim([myMin myMax]); end From c03c783a6364e8f5591200e0efdac277b1507bac Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Fri, 28 Apr 2023 12:05:09 +0200 Subject: [PATCH 04/10] implementing parsing of summary files (optics/geometry/rMatrix) at specific element --- MADX-optics/GetColumnsAndMappingTFS.m | 182 ++++++++++++++++++-------- MADX-optics/ParseTfsTable.m | 26 ++-- 2 files changed, 147 insertions(+), 61 deletions(-) diff --git a/MADX-optics/GetColumnsAndMappingTFS.m b/MADX-optics/GetColumnsAndMappingTFS.m index b2b03a3..84a30d2 100644 --- a/MADX-optics/GetColumnsAndMappingTFS.m +++ b/MADX-optics/GetColumnsAndMappingTFS.m @@ -1,5 +1,5 @@ function [ colNames, colUnits, colFacts, mapping, readFormat ] = ... - GetColumnsAndMappingTFS(whichTFS) + GetColumnsAndMappingTFS(whichTFS,genBy) % GetColumnsAndMappingTFS return column names, units and column % mapping of twiss tables % @@ -8,9 +8,12 @@ % % input arguments: % whichData: format of the table: -% 'optics': optics table by MAD-X; -% 'geometry': lattice geometry by MAD-X; -% 'rmatrix': (selected elements of) response matrix by MAD-X; +% 'OPTICS': optics table by MAD-X; +% 'GEOMETRY': lattice geometry by MAD-X; +% 'RMATRIX': (selected elements of) response matrix by MAD-X; +% genBy: format of the table: +% 'TWISS': a TFS table as generated by a MADX TWISS command; +% 'SCAN': a TFS table as generated by a current scan; % % output arguments: % colNames: array with name of columns/variables; @@ -19,64 +22,141 @@ % mapping: on which column a given variable is found; % readFormat: format string necessary to correctly parse the file; % -% This function states that the optics and geometry TFS tables have -% specific formats: -% - TFS table for optics: +% This function states that the optics, geometry and Rmatrix TFS tables +% have specific formats: +% - generated by a MADX TWISS command (genBy=="TWISS") +% . optics: % NAME, KEYWORD, L, S, BETX, ALFX, BETY, ALFY, X, PX, Y, PY, % DX, DPX, DY, DPY, MUX, MUY; -% - TFS table for geometry: +% . geometry: % NAME, KEYWORD, L, S, KICK, HKICK, VKICK, ANGLE, K0L, K1L, K2L, % APERTYPE, APER_1, APER_2, APER_3, APER_4, APOFF_1, APOFF_2; -% - TFS table for response matrix: +% . response matrix: % NAME, KEYWORD, L, S, RE11, RE12, RE21, RE22, RE16, RE26, RE33, RE34, % RE43, RE44, RE36, RE46, RE51, RE52, RE55, RE56, RE66; % +% - generated by a scan (genBy=="SCAN") +% . optics: +% Brho[Tm], BP[mm], myID[], BETX, ALFX, BETY, ALFY, X, PX, Y, PY, +% DX, DPX, DY, DPY, MUX, MUY; +% . geometry: +% Brho[Tm], BP[mm], myID[], KICK, HKICK, VKICK, ANGLE, K0L, K1L, K2L, +% APERTYPE, APER_1, APER_2, APER_3, APER_4, APOFF_1, APOFF_2; +% . response matrix: +% Brho[Tm], BP[mm], myID[], RE11, RE12, RE21, RE22, RE16, RE26, RE33, RE34, +% RE43, RE44, RE36, RE46, RE51, RE52, RE55, RE56, RE66; +% % See also ParseTfsTable. + if (~exist("genBy","var") | ismissing(genBy)), genBy="TWISS"; end % default genBy % from label to data column - switch lower(whichTFS) - case {'opt','optics'} - colNames=[ "NAME" "KEYWORD" "L" "S" "BETX" "ALFX" "BETY" "ALFY" ... - "X" "PX" "Y" "PY" "DX" "DPX" "DY" "DPY" ... - "MUX" "MUY" ]; - colUnits=[ "" "" "m" "m" "m" "" "m" "" ... - "m" "" "m" "" "m" "" "m" "" ... - "2\pi" "2\pi"]; - mapping =[ 1 2 3 4 5 6 7 8 ... - 9 10 11 12 13 14 15 16 ... - 17 18 ]; - colFacts = ones(1,length(colNames)); - readFormat = '%s %s %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f'; - case {'geo','geometry'} - colNames=[ "NAME" "KEYWORD" "L" "S" ... - "KICK" "HKICK" "VKICK" "ANGLE" "K0L" "K1L" "K2L" ... - "APERTYPE" "APER_1" "APER_2" "APER_3" "APER_4" "APOFF_1" "APOFF_2" ]; - colUnits=[ "" "" "m" "m" ... - "rad" "rad" "rad" "rad" "rad" "rad m^{-1}" "rad m^{-2}" ... - "" "m" "m" "m" "m" "m" "m" ]; - mapping =[ 1 2 3 4 ... - 5 6 7 8 9 10 11 ... - 12 13 14 15 16 17 18 ]; - colFacts = ones(1,length(colNames)); - readFormat = '%s %s %f %f %f %f %f %f %f %f %f %s %f %f %f %f %f %f %f'; - case {'rm','rmatrix'} - colNames=[ "NAME" "KEYWORD" "L" "S" ... - "RE11" "RE12" "RE21" "RE22" "RE16" "RE26" ... - "RE33" "RE34" "RE43" "RE44" "RE36" "RE46" ... - "RE51" "RE52" "RE55" "RE56" "RE66" ]; - colUnits=[ "" "" "m" "m" ... - "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... - "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... - "s m^{-1}" "s rad^{-1}" "" "s" "m" ]; - mapping =[ 1 2 3 4 ... - 5 6 7 8 9 10 ... - 11 12 13 14 15 16 ... - 17 18 19 20 21 ]; - colFacts = ones(1,length(colNames)); - readFormat = '%s %s %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f'; + switch upper(genBy) + case "TWISS" + switch upper(whichTFS) + case {'OPT','OPTICS'} + colNames=[ "NAME" "KEYWORD" "L" "S" "BETX" "ALFX" "BETY" "ALFY" ... + "X" "PX" "Y" "PY" "DX" "DPX" "DY" "DPY" ... + "MUX" "MUY" ]; + colUnits=[ "" "" "m" "m" "m" "" "m" "" ... + "m" "" "m" "" "m" "" "m" "" ... + "2\pi" "2\pi"]; + mapping =[ 1 2 3 4 5 6 7 8 ... + 9 10 11 12 13 14 15 16 ... + 17 18 ]; + colFacts = ones(1,length(colNames)); + readFormat = '%s %s %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f'; + case {'GEO','GEOMETRY'} + colNames=[ "NAME" "KEYWORD" "L" "S" ... + "KICK" "HKICK" "VKICK" "ANGLE" "K0L" "K1L" "K2L" ... + "APERTYPE" "APER_1" "APER_2" "APER_3" "APER_4" "APOFF_1" "APOFF_2" ]; + colUnits=[ "" "" "m" "m" ... + "rad" "rad" "rad" "rad" "rad" "rad m^{-1}" "rad m^{-2}" ... + "" "m" "m" "m" "m" "m" "m" ]; + mapping =[ 1 2 3 4 ... + 5 6 7 8 9 10 11 ... + 12 13 14 15 16 17 18 ]; + colFacts = ones(1,length(colNames)); + readFormat = '%s %s %f %f %f %f %f %f %f %f %f %s %f %f %f %f %f %f %f'; + case {'RM','RMATRIX'} + colNames=[ "NAME" "KEYWORD" "L" "S" ... + "RE11" "RE12" "RE21" "RE22" "RE16" "RE26" ... + "RE33" "RE34" "RE43" "RE44" "RE36" "RE46" ... + "RE51" "RE52" "RE55" "RE56" "RE66" ]; + colUnits=[ "" "" "m" "m" ... + "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... + "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... + "s m^{-1}" "s rad^{-1}" "" "s" "m" ]; + mapping =[ 1 2 3 4 ... + 5 6 7 8 9 10 ... + 11 12 13 14 15 16 ... + 17 18 19 20 21 ]; + colFacts = ones(1,length(colNames)); + readFormat = '%s %s %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f'; + otherwise + error('which column mapping for TFS table?'); + end + case "SCAN" + mySep=','; + switch upper(whichTFS) + case {'OPT','OPTICS'} + colNames=[ "BRHO" "BP" "ID" ... + "BETX" "ALFX" "BETY" "ALFY" ... + "X" "PX" "Y" "PY" ... + "DX" "DPX" "DY" "DPY" ... + "MUX" "MUY" ]; + colUnits=[ "Tm" "mm" "" ... + "m" "" "m" "" ... + "m" "" "m" "" ... + "m" "" "m" "" ... + "2\pi" "2\pi"]; + mapping =[ 1 2 3 ... + 4 5 6 7 ... + 8 9 10 11 ... + 12 13 14 15 ... + 16 17 ]; + colFacts = ones(1,length(colNames)); + readFormat = strcat('%f',mySep,'%f',mySep,'%d',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... + '%f',mySep,'%f' ); + case {'GEO','GEOMETRY'} + colNames=[ "BRHO" "BP" "ID" ... + "KICK" "HKICK" "VKICK" "ANGLE" "K0L" "K1L" "K2L" ... + "APERTYPE" "APER_1" "APER_2" "APER_3" "APER_4" "APOFF_1" "APOFF_2" ]; + colUnits=[ "Tm" "mm" "" ... + "rad" "rad" "rad" "rad" "rad" "rad m^{-1}" "rad m^{-2}" ... + "" "m" "m" "m" "m" "m" "m" ]; + mapping =[ 1 2 3 ... + 4 5 6 7 8 9 10 ... + 11 12 13 14 15 16 17 ]; + colFacts = ones(1,length(colNames)); + readFormat = strcat('%f',mySep,'%f',mySep,'%d',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f' ); + case {'RM','RMATRIX'} + colNames=[ "BRHO" "BP" "ID" ... + "RE11" "RE12" "RE21" "RE22" "RE16" "RE26" ... + "RE33" "RE34" "RE43" "RE44" "RE36" "RE46" ... + "RE51" "RE52" "RE55" "RE56" "RE66" ]; + colUnits=[ "Tm" "mm" "" ... + "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... + "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... + "s m^{-1}" "s rad^{-1}" "" "s" "m" ]; + mapping =[ 1 2 3 ... + 4 5 6 7 8 9 ... + 10 11 12 13 14 15 ... + 16 17 18 19 20 ]; + colFacts = ones(1,length(colNames)); + readFormat = strcat('%f',mySep,'%f',mySep,'%d',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f' ); + otherwise + error('which column mapping for TFS table?'); + end otherwise - error('which column mapping for TFS table?'); - return + error("How was the .tfs table generated (TWISS/SCAN)? %s NOT recognised!",genBy); end end \ No newline at end of file diff --git a/MADX-optics/ParseTfsTable.m b/MADX-optics/ParseTfsTable.m index 35809bf..67d6578 100644 --- a/MADX-optics/ParseTfsTable.m +++ b/MADX-optics/ParseTfsTable.m @@ -1,4 +1,4 @@ -function tfsTable = ParseTfsTable(fileNames,whichTFS,nHeader) +function tfsTable = ParseTfsTable(fileNames,whichTFS,genBy) % ParseTfsTable read TFS table from file. % % The function can read multiple files at the same time. @@ -20,31 +20,37 @@ % and the second one the column. % % optional input arguments: -% nHeader: number of header lines; -% if not specified, the default value (i.e. 48) is used; +% genBy: format of the table: +% 'TWISS': a TFS table as generated by a MADX TWISS command; +% 'SCAN': a TFS table as generated by a current scan; % Nota Bene: the column names and formats are considered as part of the % header; % % See also GetColumnsAndMappingTFS, ParseOpticsFileHeader, GetAperture, PlotLattice. - nHeaderUsr=48; - if ( exist('nHeader','var') ) - nHeaderUsr=nHeader; + if (~exist("genBy","var") | ismissing(genBy)), genBy="TWISS"; end % default genBy + switch upper(genBy) + case "TWISS" + nHeader=48; + case "SCAN" + nHeader=1; + otherwise + error("How was the .tfs table generated (TWISS/SCAN)? %s NOT recognised!",genBy); end % get format of the file [ colNames, colUnits, colFacts, mapping, readFormat ] = ... - GetColumnsAndMappingTFS(whichTFS); + GetColumnsAndMappingTFS(whichTFS,genBy); if ( length(fileNames)==1 ) - fprintf('parsing %s file %s ...\n',whichTFS,fileNames); + fprintf('parsing %s file %s generated by %s ...\n',whichTFS,fileNames,genBy); fileID = fopen(fileNames,'r'); - tfsTable = textscan(fileID,readFormat,'HeaderLines',nHeaderUsr); + tfsTable = textscan(fileID,readFormat,'HeaderLines',nHeader); fclose(fileID); else tfsTable=[]; for fileName=fileNames - tmpTfsTable = ParseTfsTable(fileName,whichTFS,nHeaderUsr); + tmpTfsTable = ParseTfsTable(fileName,whichTFS,genBy); tfsTable=[ tfsTable ; tmpTfsTable ]; end fprintf('...acquired %i %s files.\n',length(fileNames),whichTFS); From d60f2ee71dce39f010a48249e09115e18e57eff6 Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Fri, 28 Apr 2023 12:05:51 +0200 Subject: [PATCH 05/10] adding a function to parse a table generated by MADX with magnets currents --- MADX-optics/ParseTfsTableCurrent.m | 31 ++++++++++++++++++++++++++++++ 1 file changed, 31 insertions(+) create mode 100644 MADX-optics/ParseTfsTableCurrent.m diff --git a/MADX-optics/ParseTfsTableCurrent.m b/MADX-optics/ParseTfsTableCurrent.m new file mode 100644 index 0000000..9478a16 --- /dev/null +++ b/MADX-optics/ParseTfsTableCurrent.m @@ -0,0 +1,31 @@ +function [myTable,myColNames,MagNames,iMagNames]=ParseTfsTableCurrent(fileName) + myDelimiter=","; + %% get header from first line + fid=fopen(fileName,'r'); + tmpLine=fgetl(fid); + fclose(fid); + % - crunch columns + myColNames=string(split(tmpLine,myDelimiter)); + MagNames=[]; iMagNames=[]; + for ii=1:length(myColNames) + [myName,myUnit]=DecodeColumn(myColNames(ii)); + if (startsWith(upper(myName),"I_")) + MagNames=[ MagNames extractAfter(upper(myName),"I_") ]; + iMagNames=[ iMagNames ii ]; + end + end + %% get data + myTable=readmatrix(fileName,'HeaderLines',1,'Delimiter',myDelimiter,'FileType','text'); +end + +function [myName,myUnit]=DecodeColumn(myColHeadIN) + myColHead=myColHeadIN; + if (startsWith(myColHead,"#")) + myColHead=split(myColHead,"#"); + myColHead=myColHead(2); + end + tmp=split(myColHead,"["); + myName=tmp(1); + tmp=split(tmp(2),"]"); + myUnit=tmp(1); +end From 4ccb3e0ac58a7123a09316f81543eb02f798db3c Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Thu, 11 May 2023 15:33:10 +0200 Subject: [PATCH 06/10] bug fix in parsing .tfs files from scan ...and a minor bit of code refactoring. --- MADX-optics/GetColumnsAndMappingTFS.m | 40 +++++---------------------- 1 file changed, 7 insertions(+), 33 deletions(-) diff --git a/MADX-optics/GetColumnsAndMappingTFS.m b/MADX-optics/GetColumnsAndMappingTFS.m index 84a30d2..003da87 100644 --- a/MADX-optics/GetColumnsAndMappingTFS.m +++ b/MADX-optics/GetColumnsAndMappingTFS.m @@ -60,10 +60,6 @@ colUnits=[ "" "" "m" "m" "m" "" "m" "" ... "m" "" "m" "" "m" "" "m" "" ... "2\pi" "2\pi"]; - mapping =[ 1 2 3 4 5 6 7 8 ... - 9 10 11 12 13 14 15 16 ... - 17 18 ]; - colFacts = ones(1,length(colNames)); readFormat = '%s %s %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f'; case {'GEO','GEOMETRY'} colNames=[ "NAME" "KEYWORD" "L" "S" ... @@ -72,10 +68,6 @@ colUnits=[ "" "" "m" "m" ... "rad" "rad" "rad" "rad" "rad" "rad m^{-1}" "rad m^{-2}" ... "" "m" "m" "m" "m" "m" "m" ]; - mapping =[ 1 2 3 4 ... - 5 6 7 8 9 10 11 ... - 12 13 14 15 16 17 18 ]; - colFacts = ones(1,length(colNames)); readFormat = '%s %s %f %f %f %f %f %f %f %f %f %s %f %f %f %f %f %f %f'; case {'RM','RMATRIX'} colNames=[ "NAME" "KEYWORD" "L" "S" ... @@ -86,11 +78,6 @@ "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... "s m^{-1}" "s rad^{-1}" "" "s" "m" ]; - mapping =[ 1 2 3 4 ... - 5 6 7 8 9 10 ... - 11 12 13 14 15 16 ... - 17 18 19 20 21 ]; - colFacts = ones(1,length(colNames)); readFormat = '%s %s %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f %f'; otherwise error('which column mapping for TFS table?'); @@ -109,13 +96,7 @@ "m" "" "m" "" ... "m" "" "m" "" ... "2\pi" "2\pi"]; - mapping =[ 1 2 3 ... - 4 5 6 7 ... - 8 9 10 11 ... - 12 13 14 15 ... - 16 17 ]; - colFacts = ones(1,length(colNames)); - readFormat = strcat('%f',mySep,'%f',mySep,'%d',mySep, ... + readFormat = strcat('%f',mySep,'%f',mySep,'%f',mySep, ... '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... @@ -127,11 +108,7 @@ colUnits=[ "Tm" "mm" "" ... "rad" "rad" "rad" "rad" "rad" "rad m^{-1}" "rad m^{-2}" ... "" "m" "m" "m" "m" "m" "m" ]; - mapping =[ 1 2 3 ... - 4 5 6 7 8 9 10 ... - 11 12 13 14 15 16 17 ]; - colFacts = ones(1,length(colNames)); - readFormat = strcat('%f',mySep,'%f',mySep,'%d',mySep, ... + readFormat = strcat('%f',mySep,'%f',mySep,'%f',mySep, ... '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f' ); case {'RM','RMATRIX'} @@ -143,20 +120,17 @@ "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... "" "m rad^{-1}" "rad m^{-1}" "" "m" "rad" ... "s m^{-1}" "s rad^{-1}" "" "s" "m" ]; - mapping =[ 1 2 3 ... - 4 5 6 7 8 9 ... - 10 11 12 13 14 15 ... - 16 17 18 19 20 ]; - colFacts = ones(1,length(colNames)); - readFormat = strcat('%f',mySep,'%f',mySep,'%d',mySep, ... + readFormat = strcat('%f',mySep,'%f',mySep,'%f',mySep, ... '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep, ... - '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f' ); + '%f',mySep,'%f',mySep,'%f',mySep,'%f',mySep,'%f' ); otherwise error('which column mapping for TFS table?'); end otherwise error("How was the .tfs table generated (TWISS/SCAN)? %s NOT recognised!",genBy); end + colFacts = ones(1,length(colNames)); + mapping = 1:length(colNames); -end \ No newline at end of file +end From cb55d2e0d664f540d6fe54ae6ba21343c495953e Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Fri, 12 May 2023 15:52:19 +0200 Subject: [PATCH 07/10] possibility to omit geometric emittance when converti physical units to normalised units and back In this case, an emittance of 1 m rad is considered --- optics/Norm2Phys.m | 1 + optics/Phys2Norm.m | 1 + 2 files changed, 2 insertions(+) diff --git a/optics/Norm2Phys.m b/optics/Norm2Phys.m index aba1c61..fad8619 100644 --- a/optics/Norm2Phys.m +++ b/optics/Norm2Phys.m @@ -9,6 +9,7 @@ % output: % - zOut(float(nPoints,2)): beam population (physical units [m,rad]); + if (~exist("emiG","var")), emiG=1; end % one plane at a time! TT=[1 0; -alpha/beta 1/beta]*sqrt(beta*emiG); zOut=TT*zIn'; diff --git a/optics/Phys2Norm.m b/optics/Phys2Norm.m index 6790070..1b27a04 100644 --- a/optics/Phys2Norm.m +++ b/optics/Phys2Norm.m @@ -9,6 +9,7 @@ % output: % - zOut(float(nPoints,2)): beam population (normalised units [,]); + if (~exist("emiG","var")), emiG=1; end % one plane at a time! TT=[1 0; alpha beta]/sqrt(beta*emiG); zOut=TT*zIn'; From c498dd9368bb59fc5baa1a0d746f42b78f7b828e Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Wed, 17 May 2023 12:31:09 +0200 Subject: [PATCH 08/10] adding plot of theta*L to ShowMCS.m --- educational/ShowMCS.m | 13 +++++++++++++ 1 file changed, 13 insertions(+) diff --git a/educational/ShowMCS.m b/educational/ShowMCS.m index 670e993..401c88d 100644 --- a/educational/ShowMCS.m +++ b/educational/ShowMCS.m @@ -48,6 +48,19 @@ ShowMe(theta0_C*1000,mmEquiv,"\theta [mrad]","H_2O Range [mm]",sprintf("MCS for CARBON ions of %g MeV/u",Ek_C)); set(gca, 'YScale', 'log'); end +%% show effect of scattering +L=0.15; % distance travelled by scattered bea, [m] +if (length(mmEquiv)==1 && (length(Ek_p)>1 || length(Ek_C)>1)) + % as a function of energy + ShowMe(theta0_p*L*1000,Ek_p,"\theta\timesL [mm]","E_k [MeV]",sprintf("MCS for PROTONs after traversing %g mm of water equivalent",mmEquiv)); set(gca, 'YScale', 'log'); + ShowMe(theta0_C*L*1000,Ek_C,"\theta\timesL [mm]","E_k [MeV/u]",sprintf("MCS for CARBON ions after traversing %g mm of water equivalent",mmEquiv)); set(gca, 'YScale', 'log'); +else + % as a function of material thickness + ShowMe(theta0_p*1000,mmEquiv,"\theta\timesL [mm]","H_2O Range [mm]",sprintf("MCS for PROTONs of %g MeV",Ek_p)); set(gca, 'YScale', 'log'); + ShowMe(theta0_C*1000,mmEquiv,"\theta\timesL [mm]","H_2O Range [mm]",sprintf("MCS for CARBON ions of %g MeV/u",Ek_C)); set(gca, 'YScale', 'log'); +end + + %% local functions function theta0=ComputeTheta0(beta,pp,z,x,X0) From 12bf812319cb0193fa77307b3820ef8402be8113 Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Wed, 17 May 2023 12:31:50 +0200 Subject: [PATCH 09/10] general/Plot2DHistograms.m returns axis objects At the same time, more robust checking of missing ellypses contouring distributions --- general/Plot2DHistograms.m | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/general/Plot2DHistograms.m b/general/Plot2DHistograms.m index 3b01836..013f9fd 100644 --- a/general/Plot2DHistograms.m +++ b/general/Plot2DHistograms.m @@ -1,4 +1,4 @@ -function Plot2DHistograms(showMe,showMe1DX,showMe1DY,xShow,yShow,xShowLabel,yShowLabel,contours,lHist,lSquared) +function [axM,axX,axY]=Plot2DHistograms(showMe,showMe1DX,showMe1DY,xShow,yShow,xShowLabel,yShowLabel,contours,lHist,lSquared) if ( ~exist("lHist","var") ), lHist=true; end if ( ~exist("lSquared","var") ), lSquared=false; end @@ -11,7 +11,7 @@ function Plot2DHistograms(showMe,showMe1DX,showMe1DY,xShow,yShow,xShowLabel,ySho plot(showMe(:,1),showMe(:,2),"k."); end if ( exist("contours","var") ) - if ( ~ismissing(contours) ) + if ( ~all(ismissing(contours),"all") ) for ii=1:size(contours,3) hold on; plot(contours(:,1,ii),contours(:,2,ii),".-"); end @@ -26,7 +26,7 @@ function Plot2DHistograms(showMe,showMe1DX,showMe1DY,xShow,yShow,xShowLabel,ySho if ( lSquared ) myMin=min(min(xShow),min(yShow)); myMax=max(max(xShow),max(yShow)); if ( exist("contours","var") ) - if ( ~ismissing(contours) ) + if ( ~all(ismissing(contours),"all") ) myMin=min(myMin,min(contours,[],"all")); myMax=max(myMax,max(contours,[],"all")); end From 8466a95abb9c083852d329c23c33d0c6aa944dca Mon Sep 17 00:00:00 2001 From: Mereghetti Alessio Date: Wed, 31 May 2023 13:58:20 +0200 Subject: [PATCH 10/10] DisplayBeamProfiles.m: Figs are saved only if myFigName var is defined --- DisplayBeamProfiles.m | 45 +++++++++++++++++++++++++++---------------- 1 file changed, 28 insertions(+), 17 deletions(-) diff --git a/DisplayBeamProfiles.m b/DisplayBeamProfiles.m index fb035b0..3668be2 100644 --- a/DisplayBeamProfiles.m +++ b/DisplayBeamProfiles.m @@ -33,21 +33,25 @@ % USER's input data % ------------------------------------------------------------------------- kPath="P:\Accelerating-System\Accelerator-data"; - myTit="Check fibre - Protoni, 320mm"; - monTypes=[ "QPP" "QPP" "SFP" "SFP" ]; % CAM, DDS, GIM, QPP/SFH/SFM/SFP - QBM/PMM/PIB to come - MonPaths=[... - strcat(kPath,"\scambio\Alessio\2023-02-05_check_cablaggi_SFP\HE-007A-CEB\HE_010B_QPP\HOR\PRC-544-230205-0646\") - strcat(kPath,"\scambio\Alessio\2023-02-05_check_cablaggi_SFP\HE-007A-CEB\HE_010B_QPP\HOR\PRC-544-230205-0729\") - strcat(kPath,"\scambio\Alessio\2023-02-05_check_cablaggi_SFP\HE-007A-CEB\HE_012B_SFP\HOR\PRC-544-230205-0742\") - strcat(kPath,"\scambio\Alessio\2023-02-05_check_cablaggi_SFP\HE-007A-CEB\HE_012B_SFP\HOR\PRC-544-230205-0747\") - ]; + monTypes=[ "GIM" ]; % CAM, DDS, GIM, QPP/SFH/SFM/SFP - QBM/PMM/PIB to come myLabels=[... - "HE-010B-QPP - prima di invertire i cavi" - "HE-010B-QPP - dopo aver invertito i cavi" - "HE-012B-SFP - prima di invertire i cavi" - "HE-012B-SFP - dopo aver invertito i cavi" + "H2-009B-GIM" ]; lSkip=false; % DDS summary file: skip first 2 lines (in addition to header line) + myFigPath="."; + % part-dependent stuff + % - protoni +% myFigName="summary_protoni_GIM_2023-05-09.10"; +% myTit="summary 2023-05-09.10 - Protoni"; +% MonPaths=[... +% strcat(kPath,"\Area dati MD\00Summary\Protoni\2023\Maggio\2023.05.09-10\Steering ridotti\GIM\PRC-544-230511-0147_H2-009B-GIM_AllTrig\") +% ]; + % - carbonio + myFigName="summary_carbonio_GIM_2023-05-09.10"; + myTit="summary 2023-05-09.10 - Carbonio"; + MonPaths=[... + strcat(kPath,"\Area dati MD\00Summary\Carbonio\2023\Maggio\2023.05.09-10\Steering ridotti\GIM\PRC-544-230511-0028_H2-009B-GIM_AllTrig\") + ]; end %% check of user input data @@ -148,9 +152,17 @@ end end % - 3D plot of profiles -ShowSpectra(profiles,sprintf("%s - 3D profiles",myTit),addIndex,addLabel,myLabels,strcat(myFigPath,"\3Dprofiles_",myFigName,".fig")); +if (exist("myFigName","var") & ~ismissing(myFigName)) + ShowSpectra(profiles,sprintf("%s - 3D profiles",myTit),addIndex,addLabel,myLabels,strcat(myFigPath,"\3Dprofiles_",myFigName,".fig")); +else + ShowSpectra(profiles,sprintf("%s - 3D profiles",myTit),addIndex,addLabel,myLabels); +end % - statistics on profiles -ShowBeamProfilesSummaryData(BARsProf,FWHMsProf,INTsProf,missing(),addIndex,addLabel,myLabels,missing(),myTit,strcat(myFigPath,"\Stats_",myFigName,".fig")); +if (exist("myFigName","var") & ~ismissing(myFigName)) + ShowBeamProfilesSummaryData(BARsProf,FWHMsProf,INTsProf,missing(),addIndex,addLabel,myLabels,missing(),myTit,strcat(myFigPath,"\Stats_",myFigName,".fig")); +else + ShowBeamProfilesSummaryData(BARsProf,FWHMsProf,INTsProf,missing(),addIndex,addLabel,myLabels,missing(),myTit); +end % - statistics on profiles vs summary files for iDataAcq=1:nDataSets switch upper(monTypes(iDataAcq)) @@ -159,9 +171,8 @@ CompBars=BARsSumm(:,:,iDataAcq); CompBars(:,:,2)=BARsProf(:,:,iDataAcq); CompFwhms=FWHMsSumm(:,:,iDataAcq); CompFwhms(:,:,2)=FWHMsProf(:,:,iDataAcq); CompInts=INTsSumm(:,:,iDataAcq); CompInts(:,:,2)=INTsProf(:,:,iDataAcq); - % CompXs=mmsSumm(:,iDataAcq); - CompXs=(1:size(BARsSumm,1))'; - CompXs(:,2)=addIndex(:,iDataAcq); +% CompXs=mmsSumm(:,iDataAcq); CompXs(:,2)=mmsProf(:,iDataAcq); + CompXs=(1:size(BARsSumm,1))'; CompXs(:,2)=addIndex(:,iDataAcq); ShowBeamProfilesSummaryData(CompBars,CompFwhms,CompInts,missing(),CompXs,addLabel,... ["summary data" "stat on profiles"],missing(),sprintf("%s - %s - summary vs profile stats",myTit,myLabels(iDataAcq))); end