MUON/AliMUONSegment.cxx

   1 /**************************************************************************
   2  * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
   3  *                                                                        *
   4  * Author: The ALICE Off-line Project.                                    *
   5  * Contributors are mentioned in the code where appropriate.              *
   6  *                                                                        *
   7  * Permission to use, copy, modify and distribute this software and its   *
   8  * documentation strictly for non-commercial purposes is hereby granted   *
   9  * without fee, provided that the above copyright notice appears in all   *
  10  * copies and that both the copyright notice and this permission notice   *
  11  * appear in the supporting documentation. The authors make no claims     *
  12  * about the suitability of this software for any purpose. It is          *
  13  * provided "as is" without express or implied warranty.                  *
  14  **************************************************************************/
  15
  16 /*
  17 $Log$
  18 Revision 1.8.10.1  2002/10/11 06:56:48  hristov
  19 Updating VirtualMC to v3-09-02
  20
  21 Revision 1.9  2002/09/20 13:32:26  cussonno
  22 Minor bugs in the definition of the bending impact parameter corrected (thanks to A. Zinchenko)
  23
  24 Revision 1.8  2001/04/25 14:50:42  gosset
  25 Corrections to violations of coding conventions
  26
  27 Revision 1.7  2001/02/08 10:34:41  gosset
  28 Add a "real" default constructor.
  29
  30 Revision 1.6  2001/02/05 14:49:29  hristov
  31 Compare() declared const (R.Brun)
  32
  33 Revision 1.5  2001/01/08 11:01:02  gosset
  34 Modifications used for addendum to Dimuon TDR (JP Cussonneau):
  35 *. MaxBendingMomentum to make both a segment and a track (default 500)
  36 *. MaxChi2 per degree of freedom to make a track (default 100)
  37 *. MinBendingMomentum used also to make a track
  38    and not only a segment (default 3)
  39 *. wider roads for track search in stations 1 to 3
  40 *. extrapolation to actual Z instead of Z(chamber) in FollowTracks
  41 *. in track fit:
  42    - limits on parameters X and Y (+/-500)
  43    - covariance matrices in double precision
  44    - normalization of covariance matrices before inversion
  45    - suppression of Minuit printouts
  46 *. correction against memory leak (delete extrapHit) in FollowTracks
  47 *. RMax to 10 degrees with Z(chamber) instead of fixed values;
  48    RMin and Rmax cuts suppressed in NewHitForRecFromGEANT,
  49    because useless with realistic geometry
  50
  51 Revision 1.4  2000/06/30 10:15:48  gosset
  52 Changes to EventReconstructor...:
  53 precision fit with multiple Coulomb scattering;
  54 extrapolation to vertex with Branson correction in absorber (JPC)
  55
  56 Revision 1.3  2000/06/25 13:06:39  hristov
  57 Inline functions moved from *.cxx to *.h files instead of forward declarations
  58
  59 Revision 1.2  2000/06/15 07:58:48  morsch
  60 Code from MUON-dev joined
  61
  62 Revision 1.1.2.4  2000/06/12 10:10:21  morsch
  63 Dummy copy constructor and assignment operator added
  64
  65 Revision 1.1.2.3  2000/06/09 21:01:16  morsch
  66 Make includes consistent with new file structure.
  67
  68 Revision 1.1.2.2  2000/06/09 12:58:05  gosset
  69 Removed comment beginnings in Log sections of .cxx files
  70 Suppressed most violations of coding rules
  71
  72 Revision 1.1.2.1  2000/06/07 14:44:53  gosset
  73 Addition of files for track reconstruction in C++
  74 */
  75
  76 ///////////////////////////////////////////////////////////
  77 //
  78 // Segment for reconstruction
  79 // in
  80 // ALICE
  81 // dimuon
  82 // spectrometer:
  83 // two hits for reconstruction in the two chambers of one station
  84 //
  85 ///////////////////////////////////////////////////////////
  86
  87 #include "AliMUON.h"
  88 #include "AliMUONChamber.h"
  89 #include "AliMUONHitForRec.h"
  90 #include "AliMUONSegment.h"
  91 #include "AliMUONTrackParam.h"
  92 #include "AliRun.h" // for gAlice
  93
  94 ClassImp(AliMUONSegment) // Class implementation in ROOT context
  95
  96   //__________________________________________________________________________
  97 AliMUONSegment::AliMUONSegment()
  98 {
  99   // Default constructor
 100   fHitForRecPtr1 = 0; // pointer to HitForRec in first chamber
 101   fHitForRecPtr2 = 0; // pointer to HitForRec in second chamber
 102   // Bending plane:
 103   fBendingCoor = 0.0; // Coordinate in bending plane
 104   fBendingSlope = 0.0; // Slope in bending plane
 105   // Covariance in bending plane:
 106   fBendingCoorReso2 = 0.0; // Covariance(coordinate C1 in first chamber)
 107   fBendingSlopeReso2 = 0.0; // Covariance(slope)
 108   fBendingCoorSlopeReso2 = 0.0; // Covariance(C1,slope)
 109   fBendingImpact = 0.0; // Impact parameter in bending plane
 110   // Non Bending plane:
 111   fNonBendingCoor = 0.0; // Coordinate in non bending plane
 112   fNonBendingSlope = 0.0; // Slope in non bending plane
 113   // Covariance in non bending plane:
 114   fNonBendingCoorReso2 = 0.0; // Covariance(coordinate C1 in first chamber)
 115   fNonBendingSlopeReso2 = 0.0; // Covariance(slope)
 116   fNonBendingCoorSlopeReso2 = 0.0; // Covariance(C1,slope)
 117   fNonBendingImpact = 0.0; // Impact parameter in non bending plane
 118   fInTrack = kFALSE; // TRUE if segment belongs to one track
 119 }
 120
 121   //__________________________________________________________________________
 122 AliMUONSegment::AliMUONSegment(AliMUONHitForRec* Hit1, AliMUONHitForRec* Hit2)
 123 {
 124   // Constructor for AliMUONSegment from two HitForRec's,
 125   // one, in the first chamber of the station, pointed to by "Hit1",
 126   // the other one, in the second chamber of the station, pointed to by "Hit1".
 127   // Fills the pointers to both hits,
 128   // the slope, the covariance for (coordinate in first chamber, slope),
 129   // and the impact parameter at vertex (Z=0),
 130   // in bending and non bending planes.
 131   // Puts the "fInTrack" flag to "kFALSE".
 132   Double_t dz;
 133   // pointers to HitForRec's
 134   fHitForRecPtr1 = Hit1;
 135   fHitForRecPtr2 = Hit2;
 136   dz = Hit1->GetZ() - Hit2->GetZ();
 137   // bending plane
 138   fBendingCoor = Hit1->GetBendingCoor();
 139   fBendingSlope = (fBendingCoor - Hit2->GetBendingCoor()) / dz;
 140   fBendingImpact = fBendingCoor - Hit1->GetZ() * fBendingSlope;
 141   fBendingCoorReso2 = Hit1->GetBendingReso2();
 142   fBendingSlopeReso2 = ( Hit1->GetBendingReso2() +
 143                          Hit2->GetBendingReso2() ) / dz / dz;
 144   fBendingCoorSlopeReso2 = Hit1->GetBendingReso2() / dz;
 145   // non bending plane
 146   fNonBendingCoor = Hit1->GetNonBendingCoor();
 147   fNonBendingSlope = (fNonBendingCoor - Hit2->GetNonBendingCoor()) / dz;
 148   fNonBendingImpact = fNonBendingCoor - Hit1->GetZ() * fNonBendingSlope;
 149   fNonBendingCoorReso2 = Hit1->GetNonBendingReso2();
 150   fNonBendingSlopeReso2 = ( Hit1->GetNonBendingReso2() +
 151                             Hit2->GetNonBendingReso2() ) / dz / dz;
 152   fNonBendingCoorSlopeReso2 = Hit1->GetNonBendingReso2() / dz;
 153   // "fInTrack" flag to "kFALSE"
 154   fInTrack = kFALSE;
 155   return;
 156 }
 157
 158 AliMUONSegment::AliMUONSegment (const AliMUONSegment& MUONSegment)
 159 {
 160 // Dummy copy constructor
 161 }
 162
 163 AliMUONSegment & AliMUONSegment::operator=(const AliMUONSegment& MUONSegment)
 164 {
 165 // Dummy assignment operator
 166     return *this;
 167 }
 168
 169   //__________________________________________________________________________
 170 Int_t AliMUONSegment::Compare(const TObject* Segment) const
 171 {
 172   // "Compare" function to sort with increasing absolute value
 173   // of the "impact parameter" in bending plane.
 174   // Returns -1 (0, +1) if |impact parameter| of current Segment
 175   // is smaller than (equal to, larger than) |impact parameter| of Segment
 176   if (TMath::Abs(((AliMUONSegment*)this)->fBendingImpact)
 177       < TMath::Abs(((AliMUONSegment*)Segment)->fBendingImpact))
 178     return(-1);
 179   // continuous parameter, hence no need for testing equal case
 180   else return(+1);
 181 }
 182
 183   //__________________________________________________________________________
 184 Double_t AliMUONSegment::NormalizedChi2WithSegment(AliMUONSegment* Segment, Double_t Sigma2Cut)
 185 {
 186   // Calculate the normalized Chi2 between the current Segment (this)
 187   // and the Segment pointed to by "Segment",
 188   // i.e. the square deviations between the coordinates and the slopes,
 189   // in both the bending and the non bending plane,
 190   // divided by the variance of the same quantities and by "Sigma2Cut".
 191   // Returns 5 if none of the 4 quantities is OK,
 192   // something smaller than or equal to 4 otherwise.
 193   // Would it be more correct to use a real chi square
 194   // including the non diagonal term ????
 195   Double_t chi2, chi2Max, diff, normDiff;
 196   chi2 = 0.0;
 197   chi2Max = 5.0;
 198   // coordinate in bending plane
 199   diff = this->fBendingCoor - Segment->fBendingCoor;
 200   normDiff = diff * diff /
 201     (this->fBendingCoorReso2 + Segment->fBendingCoorReso2) / Sigma2Cut;
 202   if (normDiff > 1.0) return chi2Max;
 203   chi2 = chi2 + normDiff;
 204   // slope in bending plane
 205   diff = this->fBendingSlope - Segment->fBendingSlope;
 206   normDiff = diff * diff /
 207     (this->fBendingSlopeReso2 + Segment->fBendingSlopeReso2) / Sigma2Cut;
 208   if (normDiff > 1.0) return chi2Max;
 209   chi2 = chi2 + normDiff;
 210   // coordinate in non bending plane
 211   diff = this->fNonBendingCoor - Segment->fNonBendingCoor;
 212   normDiff = diff * diff /
 213     (this->fNonBendingCoorReso2 + Segment->fNonBendingCoorReso2) / Sigma2Cut;
 214   if (normDiff > 1.0) return chi2Max;
 215   chi2 = chi2 + normDiff;
 216   // slope in non bending plane
 217   diff = this->fNonBendingSlope - Segment->fNonBendingSlope;
 218   normDiff = diff * diff /
 219     (this->fNonBendingSlopeReso2 + Segment->fNonBendingSlopeReso2) / Sigma2Cut;
 220   if (normDiff > 1.0) return chi2Max;
 221   chi2 = chi2 + normDiff;
 222   return chi2;
 223 }
 224
 225   //__________________________________________________________________________
 226 AliMUONSegment* AliMUONSegment::CreateSegmentFromLinearExtrapToStation (Int_t Station, Double_t MCSfactor)
 227 {
 228   // Extrapolates linearly the current Segment (this) to station (0..) "Station".
 229   // Multiple Coulomb scattering calculated from "MCSfactor"
 230   // corresponding to one chamber,
 231   // with one chamber for the coordinate, two chambers for the angle,
 232   // due to the arrangement in stations.
 233   // Valid from station(1..) 4 to 5 or vice versa.
 234   // Returns the pointer to the created AliMUONSegment object
 235   // corresponding to this extrapolation.
 236   // The caller has the responsibility to delete this object.
 237   AliMUONSegment* extrapSegment = new AliMUONSegment(); // creates empty new segment
 238   // dZ from first hit of current Segment to first chamber of station "Station"
 239   AliMUON *pMUON = (AliMUON*) gAlice->GetModule("MUON"); // necessary ????
 240   Double_t dZ =
 241     (&(pMUON->Chamber(2 * Station)))->Z() - (this->fHitForRecPtr1)->GetZ();
 242   // Data in bending plane
 243   //  coordinate
 244   extrapSegment->fBendingCoor = this->fBendingCoor + this->fBendingSlope * dZ;
 245   //  slope
 246   extrapSegment->fBendingSlope = this->fBendingSlope;
 247   //  covariance, including multiple Coulomb scattering over dZ due to one chamber
 248   extrapSegment->fBendingCoorReso2 = this->fBendingCoorReso2 +
 249     (this->fBendingSlopeReso2 + MCSfactor) * dZ * dZ; // missing non diagonal term: "2.0 * this->fBendingCoorSlopeReso2 * dZ" !!!!
 250   extrapSegment->fBendingSlopeReso2 = this->fBendingSlopeReso2 + 2.0 * MCSfactor;
 251   extrapSegment->fBendingCoorSlopeReso2 =
 252     this->fBendingCoorSlopeReso2 + this->fBendingSlopeReso2 * dZ; // missing: contribution from multiple Coulomb scattering !!!!
 253   // Data in non bending plane
 254   //  coordinate
 255   extrapSegment->fNonBendingCoor =
 256     this->fNonBendingCoor + this->fNonBendingSlope * dZ;
 257   //  slope
 258   extrapSegment->fNonBendingSlope = this->fNonBendingSlope;
 259   //  covariance, including multiple Coulomb scattering over dZ due to one chamber
 260   extrapSegment->fNonBendingCoorReso2 = this->fNonBendingCoorReso2 +
 261     (this->fNonBendingSlopeReso2 + MCSfactor) *dZ * dZ; // missing non diagonal term: "2.0 * this->fNonBendingCoorSlopeReso2 * dZ" !!!!
 262   extrapSegment->fNonBendingSlopeReso2 =
 263     this->fNonBendingSlopeReso2 + 2.0 * MCSfactor;
 264   extrapSegment->fNonBendingCoorSlopeReso2 =
 265     this->fNonBendingCoorSlopeReso2 + this->fNonBendingSlopeReso2 * dZ; // missing: contribution from multiple Coulomb scattering !!!!
 266   return extrapSegment;
 267 }
 268
 269   //__________________________________________________________________________
 270 AliMUONHitForRec* AliMUONSegment::CreateHitForRecFromLinearExtrapToChamber (Int_t Chamber, Double_t MCSfactor)
 271 {
 272   // Extrapolates linearly the current Segment (this) to chamber(0..) "Chamber".
 273   // Multiple Coulomb scattering calculated from "MCSfactor"
 274   // corresponding to one chamber.
 275   // Valid from station(1..) 4 to 5 or vice versa.
 276   // Returns the pointer to the created AliMUONHitForRec object
 277   // corresponding to this extrapolation.
 278   // The caller has the responsibility to delete this object.
 279   AliMUONHitForRec* extrapHitForRec = new AliMUONHitForRec(); // creates empty new HitForRec
 280   // dZ from first hit of current Segment to chamber
 281   AliMUON *pMUON = (AliMUON*) gAlice->GetModule("MUON"); // necessary ????
 282   Double_t dZ =
 283     (&(pMUON->Chamber(Chamber)))->Z() - (this->fHitForRecPtr1)->GetZ();
 284   // Data in bending plane
 285   //  coordinate
 286   extrapHitForRec->SetBendingCoor(this->fBendingCoor + this->fBendingSlope * dZ);
 287   //  covariance, including multiple Coulomb scattering over dZ due to one chamber
 288   extrapHitForRec->SetBendingReso2(this->fBendingCoorReso2 +
 289                                    (this->fBendingSlopeReso2 + MCSfactor) * dZ * dZ); // missing non diagonal term: "2.0 * this->fBendingCoorSlopeReso2 * dZ" !!!!
 290   // Data in non bending plane
 291   //  coordinate
 292  extrapHitForRec ->SetNonBendingCoor(this->fNonBendingCoor +
 293                                      this->fNonBendingSlope * dZ);
 294   //  covariance, including multiple Coulomb scattering over dZ due to one chamber
 295   extrapHitForRec->
 296     SetNonBendingReso2(this->fNonBendingCoorReso2 +
 297                        (this->fNonBendingSlopeReso2 + MCSfactor) * dZ * dZ); // missing non diagonal term: "2.0 * this->fNonBendingCoorSlopeReso2 * dZ" !!!!
 298   return extrapHitForRec;
 299 }
 300
 301   //__________________________________________________________________________
 302 void AliMUONSegment::UpdateFromStationTrackParam(AliMUONTrackParam *TrackParam, Double_t MCSfactor, Double_t Dz1, Double_t Dz2, Double_t Dz3, Int_t Station, Double_t InverseMomentum)
 303 {
 304   // Fill data members with values calculated from the array of track parameters
 305   // pointed to by "TrackParam" (index = 0 and 1 for first and second chambers
 306   // of the station, respectively).
 307   // Multiple Coulomb scattering is taking into account with "MCSfactor"
 308   // corresponding to one chamber,
 309   // with one chamber for the coordinate, two chambers for the angle,
 310   // due to the arrangement in stations.
 311   // Resolution coming from:
 312   // coordinate in closest station at "Dz1" from current "Station",
 313   // slope between closest stations, with "Dz2" interval between them,
 314   // interval "Dz3" between chambers of closest station,
 315   // extrapolation over "Dz1" from closest station,
 316   // "InverseMomentum".
 317   // When called, "fBendingCoorReso2" and "fNonBendingCoorReso2"
 318   // are assumed to be filled
 319   // with the variance on bending and non bending coordinates.
 320   // The "road" is parametrized from the old reco_muon.F
 321   // with 8 cm between stations.
 322   AliMUONTrackParam *param0;
 323 //   Double_t cReso2, sReso2;
 324   // parameters to define the widths of the searching roads in station 0,1,2
 325   // width = p0 + p1/ (momentum)^2
 326   //                  station number:        0         1          2
 327 //   static Double_t p0BendingCoor[3] =     { 6.43e-2, 1.64e-2,   0.034 };
 328 //   static Double_t p1BendingCoor[3] =     {    986.,    821.,    446. };
 329 //   static Double_t p0BendingSlope[3] =    { 3.54e-6, 3.63e-6,  3.6e-6 };
 330 //   static Double_t p1BendingSlope[3] =    { 4.49e-3,  4.8e-3,   0.011 };
 331 //   static Double_t p0NonBendingCoor[3] =  { 4.66e-2, 4.83e-2,   0.049 };
 332 //   static Double_t p1NonBendingCoor[3] =  {   1444.,    866.,    354. };
 333 //   static Double_t p0NonBendingSlope[3] = { 6.14e-4, 6.49e-4, 6.85e-4 };
 334 //   static Double_t p1NonBendingSlope[3] = {      0.,      0.,      0. };
 335
 336   static Double_t p0BendingCoor[3] =     { 6.43e-2, 6.43e-2,   6.43e-2  };
 337   static Double_t p1BendingCoor[3] =     {    986.,    986.,       986. };
 338   static Double_t p0BendingSlope[3] =    {   3.6e-6,   3.6e-6,     3.6e-6  };
 339   static Double_t p1BendingSlope[3] =    {  1.1e-2,  1.1e-2,    1.1e-2  };
 340   static Double_t p0NonBendingCoor[3] =  {   0.049,   0.049,     0.049  };
 341   static Double_t p1NonBendingCoor[3] =  {   1444.,   1444.,      1444. };
 342   static Double_t p0NonBendingSlope[3] = {   6.8e-4,   6.8e-4,     6.8e-4  };
 343   static Double_t p1NonBendingSlope[3] = {      0.,      0.,         0. };
 344   param0 = &(TrackParam[0]);
 345
 346 // OLD version
 347 //   // Bending plane
 348 //   fBendingCoor = param0->GetBendingCoor(); // coordinate
 349 //   fBendingSlope = param0->GetBendingSlope(); // slope
 350 //   cReso2 = fBendingCoorReso2;
 351 //   sReso2 = 2.0 * cReso2 / Dz2 / Dz2;
 352 //   fBendingCoorReso2 = cReso2 + (sReso2 + MCSfactor) * Dz1 * Dz1;
 353 //   fBendingSlopeReso2 = sReso2 + 2.0 * MCSfactor;
 354 //   // Non bending plane
 355 //   fNonBendingCoor = param0->GetNonBendingCoor(); // coordinate
 356 //   fNonBendingSlope = param0->GetNonBendingSlope(); // slope
 357 //   cReso2 = fNonBendingCoorReso2;
 358 //   sReso2 = 2.0 * cReso2 / Dz2 / Dz2;
 359 //   fNonBendingCoorReso2 = cReso2 + (sReso2 + MCSfactor) * Dz1 * Dz1;
 360 //   fNonBendingSlopeReso2 = sReso2 + 2.0 * MCSfactor;
 361
 362   // Coordinate and slope
 363   // Bending plane
 364   fBendingCoor = param0->GetBendingCoor(); // coordinate
 365   fBendingSlope = param0->GetBendingSlope(); // slope
 366   // Non bending plane
 367   fNonBendingCoor = param0->GetNonBendingCoor(); // coordinate
 368   fNonBendingSlope = param0->GetNonBendingSlope(); // slope
 369
 370   // Resolutions
 371   // cReso2 and sReso2 have to be subtracted here from the parametrization
 372   // because they are added in the functions "NormalizedChi2WithSegment"
 373   // and "NormalizedChi2WithHitForRec"
 374   // Bending plane
 375 //   cReso2 = fBendingCoorReso2;
 376 //   sReso2 = (2. * cReso2 )/ (Dz3*Dz3) ;
 377   fBendingCoorReso2 = p0BendingCoor[Station] + p1BendingCoor[Station]*InverseMomentum*InverseMomentum ;  // - cReso2
 378   fBendingSlopeReso2 = p0BendingSlope[Station] + p1BendingSlope[Station]*InverseMomentum*InverseMomentum; //  - sReso2;
 379   // Non bending plane
 380 //   cReso2 = fNonBendingCoorReso2;
 381 //   sReso2 =  (2. * cReso2 )/ (Dz3*Dz3) ;
 382   fNonBendingCoorReso2 = p0NonBendingCoor[Station] + p1NonBendingCoor[Station]*InverseMomentum*InverseMomentum; // - cReso2;
 383   fNonBendingSlopeReso2 = p0NonBendingSlope[Station] + p1NonBendingSlope[Station]*InverseMomentum*InverseMomentum; //  - sReso2;
 384   return;
 385 }
 386
 387 // OLD function, with roads automatically calculated instead from being parametrized
 388 // kept because it would be a better solution,
 389 // if one can really find the right values.
 390 //   //__________________________________________________________________________
 391 // void AliMUONSegment::UpdateFromStationTrackParam(AliMUONTrackParam *TrackParam, Double_t MCSfactor, Double_t Dz1, Double_t Dz2)
 392 // {
 393 //   // Fill data members with values calculated from the array of track parameters
 394 //   // pointed to by "TrackParam" (index = 0 and 1 for first and second chambers
 395 //   // of the station, respectively).
 396 //   // Multiple Coulomb scattering is taking into account with "MCSfactor"
 397 //   // corresponding to one chamber,
 398 //   // with one chamber for the coordinate, two chambers for the angle,
 399 //   // due to the arrangement in stations.
 400 //   // Resolution coming from:
 401 //   // coordinate in closest station at "Dz1",
 402 //   // slope between closest stations, with "Dz2" interval between them,
 403 //   // extrapolation over "Dz" from closest station.
 404 //   // When called, "fBendingCoorReso2" and "fNonBendingCoorReso2"
 405 //   // are assumed to be filled
 406 //   // with the variance on bending and non bending coordinates.
 407 //   AliMUONTrackParam *param0;
 408 //   Double_t cReso2, sReso2;
 409 //   param0 = &(TrackParam[0]);
 410 //   // Bending plane
 411 //   fBendingCoor = param0->GetBendingCoor(); // coordinate
 412 //   fBendingSlope = param0->GetBendingSlope(); // slope
 413 //   cReso2 = fBendingCoorReso2;
 414 //   sReso2 = 2.0 * cReso2 / Dz2 / Dz2;
 415 //   fBendingCoorReso2 = cReso2 + (sReso2 + MCSfactor) * Dz1 * Dz1;
 416 //   fBendingSlopeReso2 = sReso2 + 2.0 * MCSfactor;
 417 //   // Non bending plane
 418 //   fNonBendingCoor = param0->GetNonBendingCoor(); // coordinate
 419 //   fNonBendingSlope = param0->GetNonBendingSlope(); // slope
 420 //   cReso2 = fNonBendingCoorReso2;
 421 //   sReso2 = 2.0 * cReso2 / Dz2 / Dz2;
 422 //   fNonBendingCoorReso2 = cReso2 + (sReso2 + MCSfactor) * Dz1 * Dz1;
 423 //   fNonBendingSlopeReso2 = sReso2 + 2.0 * MCSfactor;
 424 //   return;
 425 // }