Millepede-II V04-17-03
Internal silicon tracker test case (-t=BRLF)
Author
C. Kleinwort, DESY, 2022
Remarks
Click (on triangle): for (un)folding details.

Setup

Detector

Test case with a simple silicon strip tracker with ten equidistant, parallel planes in a beam telescope like configuration without magnetic field (implementation → More). The strip sensor have a size of 2 cm in X and 4 cm in Y and a thickness corresponding to 2% of a radiation length. The planes are 10 cm apart and oriented perpendicular to the Z-axis (as average beam direction). Each plane (1..10) has a layer of 50 sensors (0..49, 10 columns(X), 5 rows(Y)) measuring in the X-direction (with 20 µm resolution). Every third plane (1,4,7,10) has an additional ±5° stereo layer. Ten thousand tracks are simulated (with $\log(p)$ uniform in 10..100 GeV/c).

Track model

The broken lines track model is used to describe the multiple scattering in the silicon sensors [ref 4,6]. For the seed trajectory the true track parameters at the first layer (before any scattering) are used. Other models with different handling of multiple scattering are available too (see implementation → More).

Alignment

As geometry reference the average corrections for some columns of modules in the first and last plane are forced to zero (using (four) linear equality constraints) to avoid correlations with the four degrees of freedom provided by the track model.

There are 700 alignment parameters in total:

  • 500 (all planes) in X-direction (global label: 50*#plane+#sensor)
  • 200 (stereo planes) in Y-direction (global label: 50*#plane+#sensor+1000)

Misalignment

All modules (except those appearing in the constraints) are mislaligned with a RMS of 100 µm in X and Y.

Size of alignment problem

The number of global parameters (here 700) and the number of constraints (NCGB, here 4) determine the (maximal) size of the global matrix. The actual (total) number of global parameters (identified by their global labels) appearing in the binary files (NTGB) may be smaller. In addition global parameters may have only few entries (number of appearances) producing possibly inaccurate results. An entries cut is available to fix parameters with too few entries. Only the (NVGB) variable global parameters are used by pede. The number of fit parameters (NFGB) is NVGB-NCGB for elimination of constraints or NVGB+NCGB for Lagrange multipliers.

Preparations

  1. Download the software package from the DESY gitlab server to target directory, e.g. (shallow clone):
     git clone --depth 1 --branch V04-12-01 \
         https://gitlab.desy.de/claus.kleinwort/millepede-ii.git target
    
  2. Create Pede executable (in target directory):
     make pede
    
  3. Optionally cleanup:
     rm mp2???.???
    

Create pede files

     ./pede -t=BRLF > pede.dump
  1. Creates the pede input binary (mp2tst.bin) and text files:

     cksum mp2???.???
     2762434670      448 mp2con.txt
     3457733667     1971 mp2str.txt
     2736700272 14919648 mp2tst.bin
    
    • constraints file mp2con.txt
       Constraint  0.0
            54   1.00000
            64   1.00000
            74   1.00000
            84   1.00000
            94   1.00000
       Constraint  0.0
          1054   1.00000
          1064   1.00000
          1074   1.00000
          1084   1.00000
          1094   1.00000
       Constraint  0.0
           504   1.00000
           514   1.00000
           524   1.00000
           534   1.00000
           544   1.00000
       Constraint  0.0
          1504   1.00000
          1514   1.00000
          1524   1.00000
          1534   1.00000
          1544   1.00000
            
    • steering file mp2str.txt
      *        *** Default test steering file ***
      
      *            Additinal files
      fortranfiles ! following bin files are fortran
      mp2con.txt   ! constraints text file
      mp2tst.bin   ! binary data file
      Cfiles       ! following bin files are Cfiles
      
      *            Selection of variable global parameters
      *entries  10 ! lower   limit on number of entries/parameter
      *entries  25 ! default limit on number of entries/parameter
      *entries  50 ! higher  limit on number of entries/parameter
      
      *            Initial parameter values, presigmas
      *Parameter        ! define parameter attributes (start  of list)
      *205  -0.01   0.  ! start value
      *206   0.01   0.  ! start value
      *215   0.01  -1.  ! fix parameter at value
      *216   0.    -1.  ! fix parameter at value
      
      *            Handling of outliers, tails etc
      *hugecut 50.0            !cut factor in iteration 0
      *chisqcut 30.0 6.0       ! cut factor in iterations 1 and 2
      *outlierdownweighting  2 ! number of internal iterations (> 1)
      *dwfractioncut      0.2  ! 0 < value < 0.5
      *presigma           0.01 ! default value for presigma
      *regularisation 1.0      ! regularisation factor
      *regularisation 1.0 0.01 ! regularisation factor, pre-sigma
      
      *            Solution methods
      method fullMINRES       3 0.01 ! (approximate) MINRES, full storage
      method sparseMINRES     3 0.01 ! (approximate) MINRES, sparse storage
      *mrestol      1.0D-8           ! epsilon for MINRES convergence
      *bandwidth 0                   ! width of MINRES precond. band matrix
      method diagonalization 3 0.001 ! diagonalization (-> millepede.eve)
      method decomposition   3 0.001 ! Cholesky decomposition
      method inversion       3 0.001 ! Gauss matrix inversion
      * last method is applied
      
      *            Additional output. monitoring
      printcounts           ! print number of entries
      *monitorresiduals     ! poor man's DMR (-> millepede.mon)
      *printrecord   1  2   ! debug printout for records
      *printrecord  -1 -1   ! debug printout for bad data records
      
      end ! optional for end-of-data
            

    With the python2 script readMilleBinary.py binary files can be read and records printed in text form: The broken lines track model defines at each (scattering) layer offsets in X (odd local labels) and Y-direction (even local labels) as fit parameeters. The measurement in each layer depends ("local derivatives") on a linear combination of the offsets in that layer depending on the orientation in the XY-plane (0.°, ±5°).

    Detected Fortran binary file
     === NR  1 185
     -g- meas.  1 73 1 1 -0.00794364139438 0.00200000009499
     local   [1]
     local   [1.0]
     global  [73]
     global  [1.0]
     -g- meas.  2 73 2 2 -0.0057476432994 0.00200000009499
     local   [3, 4]
     local   [0.9961847066879272, 0.08726999908685684]
     global  [73, 1073]
     global  [0.9961847066879272, 0.08726999908685684]
     ...
     -l- meas.  15 1 3 0 0.0 2.02223127417e-05
     local   [1, 3, 5]
     local   [2.0, -2.1052632331848145, 0.10526315867900848]
     -l- meas.  16 2 3 0 0.0 2.02223127417e-05
     local   [2, 4, 6]
     local   [2.0, -2.1052632331848145, 0.10526315867900848]
     ...
        

    Kinks constructed from three adjacent layers are describing the multiple scattering in the middle layer and are connecting the layers (in the XZ and YZ planes). The local derivatives depend on the distances between the layers.

     python2 tools/readMilleBinary.py mp2tst.bin -n 3 > mp2binary.txt
    

    The first three tracks are printed into mp2binary.txt.

  2. Runs pede and creates output files (millepede.???):

    • pede.dump: Redirected standard output, large overlap with millepede.log. Information from the preparation and execution of the global fit. E.g.:
      • Pede start: Version, configuration and start time.
         ($Id: cc523e423547ec1e335b86a2ead55e4b2e6b3040 $)
         using OpenMP (TM)
         compiled with gcc 12.2.0
        
           <  Millepede II-P starting ... Tue Nov 22 10:20:58 2022
                                          HOSTNAME
                
           head -8 pede.dump
        
      • End of preparation: Number of global parameters (total, variable, fit), constraints (number, rank, ..), solution method, storage mode and size of global matrix.
             NTGB =       700 = total number of parameters
                              (all parameters, appearing in binary files)
             NVGB =       683 = number of variable parameters
                              (appearing in fit matrix/vectors)
             NAGB =       683 = number of all parameters
                              (including Lagrange multiplier or reduced)
           NTPGRP =       700 = total number of parameter groups
           NVPGRP =       683 = number of variable parameter groups
             NFGB =       679 = number of fit parameters
             NOFF =    232903 = max number of off-diagonal elements
             NCGB =         4 = number of constraints
            NAGBN =        18 = max number of global parameters in an event
            NALCN =        28 = max number of local parameters in an event
            NAEQN =        38 = max number of equations in an event
           NCACHE =  25000000 = number of words for caching
                              cache splitting   56.7 %   3.0 %  40.2 %
        
        
         Solution method and matrix-storage mode:
              METSOL = 1:  matrix inversion
                           with           3  iterations
              MATSTO = 1:  full symmetric matrix, (n*n+n)/2 elements
         Convergence assumed, if expected dF <  0.1000E-02
         Constraints handled by elimination
        
         Rank of product matrix of constraints is           4  for           4  constraint equations
        
         QL decomposition of constraints matrix
            largest  |eigenvalue| of L:   -2.2360679774997898
            smallest |eigenvalue| of L:   -2.2360679774997898
        
         Size of global matrix:           1  MB
        
         _______________________________LOOP2-end__________________________________
                
           grep -B33 "LOOP2-end" pede.dump
        
        For sparse storage information about the sparsity structure is added (use -B48).
      • End of global fit: Quality: final Sum(Chi^2)/Sum(Ndf) (should be close to 1.), number of rejected records/tracks (if any, should be <1%).
                                                                               00000   S
          2  3 0.10259E+06 0.28E-14                0  0           0.0       0  00000 F
        
         Sum(Chi^2)/Sum(Ndf) =   102589.35168164221
                             / (       99992 -         679 )
                             =   1.0329901592101960
        
        
         _____________________________Iteration-end________________________________
                
           grep -A5 -B3 "/Sum(Ndf)" pede.dump
        
      • Pede end: End time and maximum memory allocation.
        
           <  Millepede II-P ending   ... Tue Nov 22 10:20:58 2022
        
              Peak dynamic memory allocation:    0.102192 GB
        
                
           tail -5 pede.dump
        
      • In this special case (creation of input files) the simulated and fitted (mis)alignment parameters (and their differences) are printed close to the end of the file.
           grep -A708 "Misalignment test Si tracker" pede.dump
        
        From this a steering file (fragment) mp2param.txt defining the simulated values as start values ("starting from truth") can be created, e.g.:
           grep -A708 ... | awk 'BEGIN {print "Parameter"} NR>9 {print $2,$3,0.}' > mp2param.txt
        
    • millepede.end: Single line text file with exit code and message. Be careful with results for exit code > 1.
    • millepede.log: Log file, large overlap with pede.dump.
    • millepede.his: Internal histogram file. Script for viewing with ROOT available.
    • millepede.res: Result file with alignment corrections from global fit.
      Contains one line for each global parameter with three input values and for variable parameters up to three output values:
       Parameter   ! first 3 elements per line are significant (if used as input)
              50     0.12387E-01    0.0000       0.12387E-01   0.38268E-03         185
              51    -0.10100E-02    0.0000      -0.10100E-02   0.31170E-03         350
              52     0.35959E-02    0.0000       0.35959E-02   0.28972E-03         366
              53     0.10121E-01    0.0000       0.10121E-01   0.28810E-03         331
       ...
             150      0.0000        0.0000
       ...
            
      • Label.
      • Total value (start value (usually 0.), updated from fit).
      • Pre-sigma (usually zero or negative to fix parameter).
      • Correction value from fit.
      • Error of correction for some solution methods (e.g. inversion).
      • Number of entries (only with option printcounts).

    (Output text files in draft manual.)

Save steering and output files:

   mkdir run0; cp mp2???.txt run0; cp pede.dump run0; cp millepede.??? run0

Exercises

1. Check input

Run pede in special mode to check input only (no solution calculated):

   ./pede -C mp2str.txt > pede.dump

The standard output in pede.dump contains now a lot of information about the (linear equality) constraints:

  • After "LOOP2" start: For each constraint (in order of appearance in steering files) the number of all and of variable global parameters involved.
     _________________________________LOOP2____________________________________
    
     constraint     1 :         5 parameters,        5 variable
     constraint     2 :         5 parameters,        5 variable
     constraint     3 :         5 parameters,        5 variable
     constraint     4 :         5 parameters,        5 variable
     PRPCON:           4  constraints accepted
       
    Empty constraints (without variable global parameters) will cause singularities.
     grep -A6 '_LOOP2_' pede.dump
    
  • Around "PRPCON: constraints split": The constraints are sorted by (global) label range, split into disjoint groups and combined into non overlapping blocks. Definition of sorted constraints (with reference to steering file (index, first line number)), groups (with first/last constraint) and blocks (with first/last group) and label ranges:
      Cons. sorted       index file: index, first line first label  last label
      Cons. group        index first cons.  last cons. first label  last label
      Cons. block        index first group  last group first label  last label
      Cons. sorted           1           1           1          54          94
      Cons. group            1           1           1          54          94
      Cons. block            1           1           1          54          94
      Cons. sorted           2           3          13         504         544
      Cons. group            2           2           2         504         544
      Cons. block            2           2           2         504         544
      Cons. sorted           3           2           7        1054        1094
      Cons. group            3           3           3        1054        1094
      Cons. block            3           3           3        1054        1094
      Cons. sorted           4           4          19        1504        1544
      Cons. group            4           4           4        1504        1544
      Cons. block            4           4           4        1504        1544
    
     PRPCON: constraints split into            4 (disjoint) groups,
             groups  combined  into            4 (non overlapping) blocks
             max group size (cons., par.)            1          41
             max block size (cons., par.)            1          41
             total block matrix sizes                       164                    4
       
     grep -A4 -B16 'PRPCON: constraints' pede.dump
    
  • Before "Rank of product": For each constraint group the size and the rank of the product matrix ( $\Vek{A} \Vek{A}\trans$ in equation (11)) to verify linear independence.
      Constraint group, #con, rank           1           1           1
      Constraint group, #con, rank           2           1           1
      Constraint group, #con, rank           3           1           1
      Constraint group, #con, rank           4           1           1
    
     Rank of product matrix of constraints is           4  for           4  constraint equations
       
     grep  -B6 'Rank of product' pede.dump
    
    The results file millepede.res contains now several different sections:
  • Global parameter entries. One line for each global parameter with the three input values, the number of entries from the binary files, the index of the constraint group the parameter is part of (0 for none) and the parameter status (variable or fixed). For the first parameter of a parameter group (consecutive parameter (labels) appearing always togther) the line starts with '!>'.
     ! === global parameters ===
     ! fixed-1: by pre-sigma, -2: by entries cut, -3: by iterated entries cut
     !      Label       Value     Pre-sigma         Entries Cons. group  Status
     !>        50      0.0000        0.0000             185           0  variable
     !>        51      0.0000        0.0000             350           0  variable
     !>        52      0.0000        0.0000             366           0  variable
     !>        53      0.0000        0.0000             331           0  variable
     !>        54      0.0000        0.0000             357           1  variable
     !>        55      0.0000        0.0000             373           0  variable
     !>        56      0.0000        0.0000             376           0  variable
     ...
     !>       150      0.0000        0.0000              18           0  fixed-2
     ...
       
     grep -e"\! " -e "\!>" millepede.res
    
  • Global parameter appearance. One line for each global parameter with the file and record number of the first and last appearance in the binary files, the number of binary files with appearance and the number of paired global parameters (appearing together in measurements).
     !.
     !.Appearance statistics
     !.     Label  First file and record  Last file and record   #files  #paired-par
     !.        50          1        169          1       9931          1          1
     !.        51          1         63          1       9976          1          1
     !.        52          1         15          1       9914          1          1
     !.        53          1         64          1       9992          1          1
     !.        54          1        132          1       9968          1          1
     !.        55          1         22          1       9984          1          1
     !.        56          1          2          1       9984          1          1
     ...
     !.       100          1        550          1       9931          1          0
     !.       101          1         63          1       9925          1          0
     !.       102          1         58          1       9976          1          0
     ...
       
    The first plane (labels 50..99) contains a stereo layer. There each alignment offset in X or Y is paired with one in the other (Y or X) direction.
     grep -e"\!\." millepede.res
    
  • Constraint group entries. One line for each constraint group (described by contributing (variable global) parameters) with number of constraints, number of entries and label range. Additional lines with paired label ranges.
     * === constraint groups ===
     *  Group  #Cons.     Entries First label  Last label   Paired label range
     *      1       1        2231          54          94
     *:                                                        1054 ..        1054
     *:                                                        1064 ..        1064
     *:                                                        1074 ..        1074
     *:                                                        1084 ..        1084
     *:                                                        1094 ..        1094
     *      2       1        2049         504         544
     *:                                                        1504 ..        1504
     *:                                                        1514 ..        1514
     *:                                                        1524 ..        1524
     *:                                                        1534 ..        1534
     *:                                                        1544 ..        1544
     ...
       
     grep -e"* " -e"*:" millepede.res
    
  • Constraint group appearance. One line for each constraint group (described by contributing (variable global) parameters) with the file and record number of the first and last appearance in the binary files and the number of binary files with appearance.
     *.
     *.Appearance statistics
     *.     Group  First file and record  Last file and record   files
     *.         1          1         13          1       9998          1
     *.         2          1          5          1       9963          1
     *.         3          1         13          1       9998          1
     *.         4          1          5          1       9963          1
       
     grep -e"*\." millepede.res
    

Save steering and output files:

 mkdir run1; cp mp2???.txt run1; cp pede.dump run1; cp millepede.??? run1

2. Reject bad tracks

Depending on the quality of the local (track) fit very bad tracks are rejected ( $ \chi^2 > f_i\cdot \chi^2_\mathrm{cut} $) automatically and bad tracks can be removed by the user with the "chisqcut f1 f2" option. The scaling factor $f_i$ should start with a large value $f_1$ to allow for initial misalignment and be further reduced (to 1.) in the internal iterations of the global fit.

Modify steering:

Uncomment "chisqcut" option in the steering file.

 chisqcut 30.0 6.0        ! cut factor in iterations 1 and 2
 

Rerun pede:

   ./pede mp2str.txt > pede.dump

Check output:

Watch the pede exit code in millepede.end.

     2   Ended with severe warnings (insufficient measurements)
 

The severe warnings indicate a major problem. More details are in pede.dump.

  • During iterations:
      3  4 0.10086E+06 0.17E+02                0  0           1.6      33  00000 F
                                                                           00000   S
    
      ... warning ...
      global parameters with too few (< MREQENA) accepted entries:            2
      minimum entries:            7  for label          209
    
      4  5 0.10013E+06 0.53E+01 0.781 0.062    0  0 -1 1.000  1.0      79  00000 F
     ...
       
    After rejection of bad tracks there are only 7 entries left for the global parameter with the label "209". This is less than the minumum required by MREQENA (default 10, second parameter of entries cut).
     grep -C3 '... warning ...' pede.dump
    
  • After iterations:
           ngWarningWarningWarningWarningWarningWarningWarningWarningW
           gWarningWarningWarningWarningWarningWarningWarningWarningWa
    
             Possible bad elements =           6  in global vector
             (too few accepted entries)
             (indicated in millepede.res by counts<0)
             => please check mille data and ENTRIES cut
    
           WarningWarningWarningWarningWarningWarningWarningWarningWar
       
    During the iterations six times a global parameter had too few accepted entries. All those parameters are flagged in millepede.res by an negative count value (-(#entries+1)).
     grep -C4 'too few accepted' pede.dump
    

Get affected global parameter from result file millepede.res.

Select lines with negative count value, e.g.:

       209     0.76209E-02    0.0000       0.76209E-02   0.54649E-03          -8
       249    -0.11705E-01    0.0000      -0.11705E-01   0.51044E-03          -9
 

The global parameters with labels 209/249 have finally less than 10 entries.

 awk 'NF==6 && $6<0' millepede.res

Recheck input:

  • Get entries in binary files for affected global parameters from check input results, e.g.:
     !>       209      0.0000        0.0000              35           0  variable
     !>       249      0.0000        0.0000              40           0  variable
       
    The global parameters with labels 209/249 have initialy more than 25 entries.
     awk '($1=="!>" || $1=="! ") && ($2==209 || $2==249)' run1/millepede.res
    
  • Get number of paired parameters for affected global parameters from check input results, e.g.:
     !.       209          1         39          1       9867          1          1
     !.       249          1        637          1       9547          1          1
       
    The global parameters with labels 209/249 describe offsets in a stereo layer. Therefore the X and Y offsets for each sensor are paired (appearing together for each measurement). The number of paired parameters (last column) is 1.
     awk '$1=="!." && ($2==209 || $2==249)' run1/millepede.res
    

Options to resolve severe warnings

Options to resolve the severe warnings for problematic labels 209/249:

  1. Increase the statistics (more/larger binary files) to get the minimal number of accepted entries (7) above the entries cut value mreqena (10).
  2. Decrease the entries cut value mreqena (10) below the minimal number of accepted entries (7). It must still be larger than the maximal number of paired parameters (1).
  3. Increase the entries cut value mreqenf (25) above the maximum number of entries (40) to fix those parameters.
  4. Fix those parameters in the steering file:
        Parameter
        209  0.  -1.
        249  0.  -1.
    

Save steering and output files:

mkdir run2; cp mp2???.txt run2; cp pede.dump run2; cp millepede.??? run2

3. Reject bad tracks with decreased (accepted) entries cut (Option 2)

Modify steering further:

Uncomment "entries 25" option in the steering file and add mreqena with value 5.

 entries  25 5 ! default  limit on number of entries/parameter, reduced accepted entries
 

Rerun pede:

   ./pede mp2str.txt > pede.dump

Check output:

Watch the pede exit code in millepede.end.

     0   Ended normally
 

No warnings are reported. More details are in pede.dump.

Quality of global fit.

 Data rejected in last loop:
               0  (rank deficit/NaN)            0  (Ndf=0)              0  (huge)             79  (large)

 Sum(Chi^2)/Sum(Ndf) =   100125.94658222546
                     / (       99992 -         679 )
                     =   1.0081857015921929


 _____________________________Iteration-end________________________________
 

The final Sum(Chi^2)/Sum(Ndf) is very close to 1, the rejection rate less than 1% (79/10000) and 679 (of maximal 696) parameters have been fitted.

   grep -A5 -B3 "/Sum(Ndf)" pede.dump

The quality is fine.

Save steering and output files:

   mkdir run3; cp mp2???.txt run3; cp pede.dump run3; cp millepede.??? run3

4. Reject bad tracks with increased entries cut (Option 3)

Modify steering further:

Comment "entries 25" again and uncomment "entries 50" option in the steering file.

 *entries  25 5 ! default  limit on number of entries/parameter, reduced accepted entries
 entries  50 ! higher  limit on number of entries/parameter
 

Rerun pede:

   ./pede mp2str.txt > pede.dump

Check output:

Watch the pede exit code in millepede.end.

     1   Ended with warnings (bad measurements)
 

The warnings indicate some problem. More details are in pede.dump.

  • Quality of global fit:
     Data rejected in last loop:
                   0  (rank deficit/NaN)            0  (Ndf=0)              0  (huge)            276  (large)
    
     Sum(Chi^2)/Sum(Ndf) =   104495.47920560479
                         / (       99992 -         645 )
                         =   1.0518231975359578
    
    
    
       
    The final Sum(Chi^2)/Sum(Ndf) is close to 1, but the rejection rate is allmost 3% (276/10000) and only 645 (of maximal 696) parameters have been fitted.
     grep -A5 -B3 "/Sum(Ndf)" pede.dump
    
  • After iterations:
           ngWarningWarningWarningWarningWarningWarningWarningWarningW
           gWarningWarningWarningWarningWarningWarningWarningWarningWa
    
             Fraction of rejects =   2.76 %  (should be far below 1 %)
             => please provide correct mille data
    
           WarningWarningWarningWarningWarningWarningWarningWarningWar
       
     grep -C3 'Fraction of rejects' pede.dump
    

The result is not optimal. Too many globals parameters have been fixed at the misaligned values spoiling the local/track fit.

Save steering and output files:

   mkdir run4; cp mp2???.txt run4; cp pede.dump run4; cp millepede.??? run4

5. Reject bad tracks with increased entries cut (from truth)

Modify steering further:

Add text file "mp2param.txt" with start values to the steering file (e.g. after constraints file).

 mp2con.txt   ! constraints text file
 mp2param.txt ! start values (from truth)
 mp2tst.bin   ! binary data file
 

This defines a different starting geometry (in this case the "truth").

Rerun pede:

   ./pede mp2str.txt > pede.dump

Check output:

Watch the pede exit code in millepede.end.

     0   Ended normally
 

No warnings are reported. More details are in pede.dump.

Quality of global fit.

 Data rejected in last loop:
               0  (rank deficit/NaN)            0  (Ndf=0)              0  (huge)             24  (large)

 Sum(Chi^2)/Sum(Ndf) =   99050.817948768963
                     / (       99992 -         645 )
                     =  0.99701871167492695


 _____________________________Iteration-end________________________________
 

The final Sum(Chi^2)/Sum(Ndf) is very close to 1 and the rejection rate less than 1% (24/10000) despite only 645 (of maximal 696) parameters have been fitted.

   grep -A5 -B3 "/Sum(Ndf)" pede.dump

The quality is fine. Starting values matter for fixed parameters.

Save steering and output files:

   mkdir run5; cp mp2???.txt run5; cp pede.dump run5; cp millepede.??? run5

6. Compare millepede results

With the python2 script compareResults.py two millepede result text files (millepede.res) can be compared. As only the first three columns are evaluated a steering file (fragment) defining start values can be used too, e.g:

   python2 compareResults.py mp2param.txt run4/millepede.res

This compares the results from example 4 with the "truth".

 output file         mp2compare.txt
 total parameters    700
 matched parameters  700
 mean1, mean2        -0.000189471428571 -0.000542676439286
 rms1, rms2          0.00994723953843 0.00994581992618
 correlation         0.945525591556
 

For both sets of parameters the mean and RMS and the correlation of the sets is calculated.

A text file (default mp2compare.txt) is created.

It contains the parameter labels and values and can be directly used as ROOT input (with "tree->ReadFile()"):

 label/I:value_1/F:value_2/F
        50      0.01339     0.012557
        51     -0.00037  -0.00082443
        52      0.00403    0.0037749
        53      0.01061     0.010278
        54           -0   0.00020854
        55     -0.00712   -0.0075308
 ...
      1545     -0.01218    -0.010553
      1546      -0.0029    0.0032485
      1547     -0.01605   -0.0082543
      1548     -0.01081   -0.0027929
      1549      0.00236     0.008585
 

For all labels appearing in boths files ("matched"): label value1 value2. For all labels appearing only in one file: -label, value, value

Other examples

Other examples (C++, Python2 or Python3) are included in the GeneralBrokenLines package to explore track fitting and track based alignment with MP2. Those will produce binary files of type Cfiles.