Changeset 4706
- Timestamp:
- 08/23/04 13:37:54 (20 years ago)
- Location:
- trunk/MagicSoft/Mars/mtemp
- Files:
-
- 3 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/MagicSoft/Mars/mtemp/MFindStars.cc
r4685 r4706 120 120 UInt_t numPixels = geom->GetNumPixels(); 121 121 122 // fix the correlation if requested123 if ( find->GetUseCorrelatedGauss() && find->GetFixCorrelation() )124 {125 Double_t tandel;126 if (par[1] != 0.0)127 tandel = par[3]/par[1];128 else129 tandel = 1.e10;130 Double_t cxx = par[2]*par[2];131 Double_t cyy = par[4]*par[4];132 par[5] = tandel/(tandel*tandel-1.0) * (cyy-cxx)/sqrt(cxx*cyy);133 }134 122 135 123 //calculate chisquare … … 205 193 fNumVar = 6; 206 194 fUseCorrelatedGauss = kTRUE; 207 fFixCorrelation = kFALSE;208 195 209 196 fNumIntegratedEvents=0; … … 287 274 //------------------------------------------------------------------------- 288 275 // 289 // Set 'fUseCorrelatedGauss' and 'fFixCorrelation'flag276 // Set 'fUseCorrelatedGauss' flag 290 277 // 291 278 // if 'fUseCorrelatedGauss' is TRUE a 2dim Gaussian with correlation 292 279 // will be fitted 293 // if 'fFixCorrelation' is TRUE the orientation of the 2dim Gaussian 294 // will be kept fixed at the angle delta, 295 // where tan(delta) = ymean/xmean 296 297 void MFindStars::SetUseCorrelatedGauss(Bool_t usecorrgauss, Bool_t fixcorr) 280 // 281 282 void MFindStars::SetUseCorrelatedGauss(Bool_t usecorrgauss) 298 283 { 299 284 fUseCorrelatedGauss = usecorrgauss; 300 fFixCorrelation = fixcorr;301 285 302 286 if (usecorrgauss) … … 445 429 TIter Next(fStars->GetList()); 446 430 while ((starpos=(MStarPos*)Next())) { 447 starpos->SetCalcValues(40,40,starpos->GetXExp(),starpos->GetYExp(),fRingInterest/2,fRingInterest/2); 448 starpos->SetFitValues (40,40,starpos->GetXExp(),starpos->GetYExp(),fRingInterest/2,fRingInterest/2,0.,-1); 431 //starpos->SetCalcValues(40,40,starpos->GetXExp(), starpos->GetYExp(), 432 // fRingInterest/2, fRingInterest/2, 433 // 0., 0., 0., 0.); 434 //starpos->SetFitValues (40,40,starpos->GetXExp(), starpos->GetYExp(), 435 // fRingInterest/2, fRingInterest/2, 0., 436 // 0., 0., 0., 0., -1); 449 437 } 450 438 … … 483 471 484 472 while(FindPixelWithMaxDC(maxPixelDC, maxPixel)) { 485 486 473 MStarPos *starpos = new MStarPos; 487 starpos->SetExpValues(maxPixelDC,maxPixel.GetX(),maxPixel.GetY()); 488 starpos->SetCalcValues(maxPixelDC,maxPixelDC,maxPixel.GetX(),maxPixel.GetY(),fRingInterest/2,fRingInterest/2); 489 starpos->SetFitValues(maxPixelDC,maxPixelDC,maxPixel.GetX(),maxPixel.GetY(),fRingInterest/2,fRingInterest/2,0.,-1); 474 475 starpos->SetExpValues(maxPixelDC,maxPixel.GetX(),maxPixel.GetY()); 476 477 //starpos->SetCalcValues(maxPixelDC, maxPixelDC, 478 // maxPixel.GetX(), maxPixel.GetY(), 479 // fRingInterest/2, fRingInterest/2, 480 // 0., 0., 0., 0.); 481 482 //starpos->SetFitValues(maxPixelDC, maxPixelDC, 483 // maxPixel.GetX(), maxPixel.GetY(), 484 // fRingInterest/2, fRingInterest/2, 0., 485 // 0., 0., 0., 0., -1); 486 490 487 fStars->GetList()->Add(starpos); 491 492 488 ShadowStar(starpos); 493 489 } … … 792 788 Float_t expY = star->GetYExp(); 793 789 794 Float_t max=0; 795 UInt_t pixmax=0; 796 Float_t meanX=0; 797 Float_t meanY=0; 798 Float_t meanSqX=0; 799 Float_t meanSqY=0; 800 Float_t sumCharge=0; 801 UInt_t usedInnerPx=0; 802 UInt_t usedOuterPx=0; 790 Float_t max =0; 791 UInt_t pixmax =0; 792 Float_t meanX =0; 793 Float_t meanY =0; 794 Float_t meanSqX =0; 795 Float_t meanSqY =0; 796 Float_t sumCharge =0; 797 Float_t sumSqCharge =0; 798 UInt_t usedInnerPx =0; 799 UInt_t usedOuterPx =0; 800 Float_t factor = 1.0; 803 801 804 802 Float_t meanXold = 1.e10; … … 882 880 meanSqY = 0.0; 883 881 sumCharge = 0.0; 882 sumSqCharge = 0.0; 884 883 usedInnerPx = 0; 885 884 usedOuterPx = 0; … … 905 904 meanSqX += charge*pixXpos*pixXpos; 906 905 meanSqY += charge*pixYpos*pixYpos; 907 sumCharge += charge; 906 sumCharge += charge; 907 sumSqCharge += charge*charge; 908 908 } 909 909 } … … 915 915 meanSqX /= sumCharge; 916 916 meanSqY /= sumCharge; 917 factor = sqrt( sumSqCharge/(sumCharge*sumCharge) ); 917 918 918 919 // stop iteration when (meanX, meanY) becomes stable … … 962 963 } 963 964 965 Float_t deltameanX = factor * rmsX; 966 Float_t deltameanY = factor * rmsY; 967 968 964 969 // Substrack pedestal DC 965 970 sumCharge-= (usedInnerPx*innerped+usedOuterPx*outerped)/(usedInnerPx+usedOuterPx); … … 967 972 968 973 969 star->SetCalcValues( sumCharge, max, meanX, meanY, rmsX, rmsY); 974 star->SetCalcValues( sumCharge, max, meanX, meanY, rmsX, rmsY, 975 0., deltameanX*deltameanX, 0., deltameanY*deltameanY); 970 976 971 977 if (rmsX <= 0. || rmsY <= 0.) … … 1004 1010 } 1005 1011 1006 // set the correlation fixed if requested1007 if (fUseCorrelatedGauss && fFixCorrelation)1008 {1009 fStep[5] = 0.0;1010 fFix[5] = 1;1011 }1012 1013 1012 1014 1013 // ------------------------------------------- … … 1053 1052 Double_t sigmaY, sigmaYError; 1054 1053 Float_t chisquare = GetChisquare(); 1055 Int_t d regrees = GetDegreesofFreedom()-fNumVar;1054 Int_t degrees = GetDegreesofFreedom()-fNumVar; 1056 1055 1057 1056 if (!ierflg) … … 1078 1077 } 1079 1078 1080 // rwagner:get error matrix1079 // get error matrix 1081 1080 Double_t matrix[5][5]; 1082 1081 gMinuit->mnemat(&matrix[0][0],5); 1083 1082 1084 1083 star->SetFitValues(integratedCharge,maxFit, meanXFit, meanYFit, 1085 sigmaX, sigmaY, chisquare, dregrees, 1086 matrix[1][1], matrix[1][3], matrix[3][3]); 1087 1088 // set the results from the correlated Gauss fit to zero 1089 star->SetCGFitValues(100.0, 100.0, 0.0, 0.0, 0.0, 0.0, 0.0, 1090 0.0, 0.0, 0.0, 0.0, -1); 1084 sigmaX, sigmaY, 0.0, 1085 matrix[1][1], matrix[1][3], matrix[3][3], 1086 chisquare, degrees); 1091 1087 } 1092 1088 //---------- for the uncorrelated Gauss fit (end) --------- … … 1180 1176 if (fMinuitPrintOutLevel>=0) 1181 1177 { 1182 *fLog << "Before calling Set CGFitValues : " << endl;1178 *fLog << "Before calling SetFitValues : " << endl; 1183 1179 *fLog << "fValues = " << fValues[0] << ", " << fValues[1] << ", " 1184 1180 << fValues[2] << ", " << fValues[3] << ", " << fValues[4] … … 1215 1211 } 1216 1212 1217 star->SetCGFitValues(charge, fValues[0], fValues[1], fValues[3], 1218 fValues[2], fValues[4], fValues[5], 1219 matrix[1][1], matrix[1][3],matrix[3][3], 1220 chi2, ndof); 1221 // set the results from the uncorrelated Gauss fit to zero 1222 star->SetFitValues(100.0, 100.0, 0.0, 0.0, 0.0, 0.0, 0.0, -1, 1223 0.0, 0.0, 0.0); 1213 star->SetFitValues(charge, fValues[0], fValues[1], fValues[3], 1214 fValues[2], fValues[4], fValues[5], 1215 matrix[1][1], matrix[1][3],matrix[3][3], 1216 chi2, ndof); 1224 1217 } 1225 1218 … … 1248 1241 Float_t starYpos = star->GetMeanY(); 1249 1242 1250 Float_t starSize = 3*star->GetSigmaMajorAxis(); 1243 Float_t starSize = 3*sqrt( star->GetSigmaX()*star->GetSigmaX() 1244 + star->GetSigmaY()*star->GetSigmaY() ); 1251 1245 1252 1246 Float_t dist = sqrt((pixXpos-starXpos)*(pixXpos-starXpos)+ -
trunk/MagicSoft/Mars/mtemp/MFindStars.h
r4685 r4706 60 60 61 61 Bool_t fUseCorrelatedGauss; 62 Bool_t fFixCorrelation;63 62 64 63 TString *fVname; … … 88 87 Int_t PostProcess(); 89 88 90 void SetUseCorrelatedGauss(Bool_t usecorrgauss = kTRUE, 91 Bool_t fixcorr = kFALSE); 89 void SetUseCorrelatedGauss(Bool_t usecorrgauss = kTRUE); 92 90 Bool_t GetUseCorrelatedGauss() {return fUseCorrelatedGauss;} 93 Bool_t GetFixCorrelation() {return fFixCorrelation;}94 91 95 92 // setters -
trunk/MagicSoft/Mars/mtemp/MSourceDirections.cc
r4697 r4706 90 90 } 91 91 92 92 93 fStars = (MStarCam*)pList->FindCreateObj(AddSerialNumber("MStarCam"),"MSourceCam"); 93 94 if (!fStars) { 94 *fLog << err << AddSerialNumber("MStarCam") << " cannot be created ... aborting" << endl; 95 *fLog << err << AddSerialNumber("MStarCam") << " with name '" 96 << "MSourceCam" << " cannot be created ... aborting" << endl; 95 97 return kFALSE; 96 98 } 97 99 100 101 98 102 MObservatory magic1; 99 103 … … 110 114 fAstro.SetGeom(*geom); 111 115 fAstro.SetObservatory(magic1); 112 116 113 117 return kTRUE; 114 118 } … … 116 120 Int_t MSourceDirections::AddDirection(Float_t ra, Float_t dec, Float_t mag, TString name) 117 121 { 122 *fLog << "MSourceDirections::AddDirection; add the direction : ra, dec, mag, name = " 123 << ra << ", " << dec << ", " << mag << ", " << name << endl; 124 118 125 Int_t rc = fAstro.AddObject(ra,dec,1,name); 119 if (rc) { 120 MStarPos *starpos = new MStarPos; 121 starpos->SetName(name); 122 fStars->GetList()->Add(starpos); 123 } 126 124 127 return rc; 125 128 } … … 127 130 Int_t MSourceDirections::Process() 128 131 { 129 //Fi st delete the previous directions in the list132 //First delete the previous directions in the list 130 133 fStars->GetList()->Delete(); 131 134 … … 137 140 TIter Next(fStars->GetList()); 138 141 while ((starpos=(MStarPos*)Next())) { 139 starpos->SetCalcValues(40,40,starpos->GetXExp(),starpos->GetYExp(),0.,0.); 140 starpos->SetFitValues (40,40,starpos->GetXExp(),starpos->GetYExp(),0.,0.,0.,1); 142 //starpos->SetCalcValues(40,40,starpos->GetXExp(),starpos->GetYExp(), 143 // 0.,0.,0., 0.,0.,0.); 144 //starpos->SetFitValues (40,40,starpos->GetXExp(),starpos->GetYExp(), 145 // 0.,0.,0., 0.,0.,0., 0., -1); 141 146 } 142 147 143 148 if (fStars->GetList()->GetSize() == 0) { 144 149 *fLog << err << GetName() << "No directions inside the chosen FOV." << endl; … … 154 159 return kTRUE; 155 160 } 161
Note:
See TracChangeset
for help on using the changeset viewer.