TPC/Upgrade/macros/AnaDelta.C

   1 #include <TStyle.h>
   2 #include <TROOT.h>
   3 #include <TAxis.h>
   4 #include <TF1.h>
   5 #include <TFile.h>
   6 #include <TH1.h>
   7 #include <THn.h>
   8 #include <TObjArray.h>
   9 #include <TObject.h>
  10 #include <TString.h>
  11 #include <TVectorT.h>
  12 #include <TCanvas.h>
  13 #include <TProfile2D.h>
  14 #include <TGraphErrors.h>
  15 #include <TTreeStream.h>
  16
  17 #include <AliExternalTrackParam.h>
  18 #include <AliTPCComposedCorrection.h>
  19 #include <AliTPCCorrectionLookupTable.h>
  20
  21 #include <AliToyMCEventGenerator.h>
  22
  23 /*
  24
  25 .L $ALICE_ROOT/TPC/Upgrade/macros/AnaDelta.C+g
  26
  27
  28 */
  29 TVectorD* MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax);
  30 TVectorD* MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax);
  31 TVectorD* MakeArbitraryBinning(const char* bins);
  32
  33 void DumpHn(THn *hn, TTreeSRedirector &stream);
  34
  35 void AnaDelta(TString file, TString outDir=".")
  36 {
  37   //
  38   //
  39   //
  40
  41   gStyle->SetOptFit();
  42
  43   TTreeSRedirector stream(Form("%s/deltas.root",outDir.Data()));
  44   gROOT->cd();
  45
  46   TFile f(file);
  47   THn *hn=(THn*)f.Get("hn");
  48
  49   DumpHn(hn, stream);
  50
  51   delete hn;
  52 }
  53
  54
  55 void AnaDeltaTree(TString file, TString outFile="deltas_tree.root")
  56 {
  57   TFile f(file);
  58   gROOT->cd();
  59   TTree *t = (TTree*)f.Get("delta");
  60   Float_t soneOverPt=0.;
  61   Float_t radius=0.;
  62   Float_t trackPhi=0.;
  63   Float_t trackY=0.;
  64   Float_t trackZ=0.;
  65   Float_t resRphi=0.;
  66   Double_t trackRes=0.;
  67   Float_t pointY=0.;
  68   Float_t pointZ=0.;
  69   Short_t npTRD=0.;
  70   Short_t event=0.;
  71
  72   t->SetBranchAddress("soneOverPt" , &soneOverPt);
  73   t->SetBranchAddress("r"          , &radius);
  74   t->SetBranchAddress("trackPhi"   , &trackPhi);
  75   t->SetBranchAddress("trackY"     , &trackY);
  76   t->SetBranchAddress("trackZ"     , &trackZ);
  77   t->SetBranchAddress("resRphi"    , &resRphi);
  78   t->SetBranchAddress("trackRes"   , &trackRes);
  79   t->SetBranchAddress("pointY"     , &pointY);
  80   t->SetBranchAddress("pointZ"     , &pointZ);
  81   t->SetBranchAddress("npTRD"      , &npTRD);
  82   t->SetBranchAddress("event"      , &event);
  83
  84   // make binning
  85   TVectorD *vR   = MakeLinBinning(10,86.,250.);
  86   TVectorD *vPhi = MakeLinBinning(18*8,0.,2*TMath::Pi());
  87   TVectorD *vZ   = MakeLinBinning(50,-250.,250.);
  88
  89   const Int_t nbins=4;
  90   Int_t bins[nbins]    = {vR->GetNrows()-1, vPhi->GetNrows()-1, vZ->GetNrows()-1, 80};
  91 //   Int_t bins[nbins]    = {16, 18*5, 50, 80};
  92   Double_t xmin[nbins] = {86. , 0.,           -250., -2.};
  93   Double_t xmax[nbins] = {250., 2*TMath::Pi(), 250.,  2.};
  94   THnF *hn = new THnF("hn", "hn", nbins, bins, xmin, xmax);
  95
  96   hn->GetAxis(0)->Set(vR  ->GetNrows()-1, vR  ->GetMatrixArray());
  97   hn->GetAxis(1)->Set(vPhi->GetNrows()-1, vPhi->GetMatrixArray());
  98   hn->GetAxis(2)->Set(vZ  ->GetNrows()-1, vZ  ->GetMatrixArray());
  99
 100   hn->GetAxis(0)->SetNameTitle("r","r (cm)");
 101   hn->GetAxis(1)->SetNameTitle("phi","#varphi");
 102   hn->GetAxis(2)->SetNameTitle("z","z (cm)");
 103   hn->GetAxis(3)->SetNameTitle("drphi","#Delta(r#varphi)");
 104
 105   for (Int_t iev=0; iev<t->GetEntries(); ++iev) {
 106     t->GetEntry(iev);
 107
 108     // cuts
 109     // -- on trd
 110     if (npTRD<2) continue;
 111     Double_t pt=1./TMath::Abs(soneOverPt);
 112     if (pt<0.8) continue;
 113
 114     Float_t resRphiRandom = resRphi*trackRes;
 115     Float_t deviation     = pointY-(trackY+resRphiRandom);
 116
 117     Double_t xx[4]={radius, trackPhi, trackZ ,deviation};
 118     hn->Fill(xx);
 119   }
 120
 121   // do fits and fill tree
 122   TTreeSRedirector stream(outFile.Data());
 123   gROOT->cd();
 124
 125   DumpHn(hn, stream);
 126
 127   stream.GetFile()->cd();
 128   hn->Write();
 129
 130   delete hn;
 131   delete vR;
 132   delete vPhi;
 133   delete vZ;
 134 }
 135
 136
 137 void AnaDeltaTree2(TString file/*, TString outDir="."*/)
 138 {
 139   //
 140   // NOTE: not finished
 141   //
 142   TFile f(file);
 143   gROOT->cd();
 144   TTree *t = (TTree*)f.Get("delta");
 145   Float_t soneOverPt=0.;
 146   Float_t radius=0.;
 147   Float_t trackPhi=0.;
 148   Float_t trackY=0.;
 149   Float_t trackZ=0.;
 150   Float_t resRphi=0.;
 151   Float_t trackRes=0.;
 152   Float_t pointY=0.;
 153   Float_t pointZ=0.;
 154   Float_t npTRD=0.;
 155   Float_t event=0.;
 156
 157   t->SetBranchAddress("soneOverPt" , &soneOverPt);
 158   t->SetBranchAddress("r"          , &radius);
 159   t->SetBranchAddress("trackPhi"   , &trackPhi);
 160   t->SetBranchAddress("trackY"     , &trackY);
 161   t->SetBranchAddress("trackZ"     , &trackZ);
 162   t->SetBranchAddress("resRphi"    , &resRphi);
 163   t->SetBranchAddress("trackRes"   , &trackRes);
 164   t->SetBranchAddress("pointY"     , &pointY);
 165   t->SetBranchAddress("pointZ"     , &pointZ);
 166   t->SetBranchAddress("npTRD"      , &npTRD);
 167   t->SetBranchAddress("event"      , &event);
 168
 169   // make binning
 170   TVectorD *vZ   = MakeLinBinning(50,-250.,250.);
 171   TVectorD *vPhi = MakeLinBinning(18*8,0.,TMath::Pi());
 172   TVectorD *vR   = MakeLinBinning(16,86.,250.);
 173
 174   TObjArray arrZ(vZ->GetNrows()-1);
 175   arrZ.SetOwner();
 176
 177   for (Int_t iev=0; iev<t->GetEntries(); ++iev) {
 178     t->GetEntry(iev);
 179
 180     // cuts
 181     // -- on trd
 182     if (npTRD<2) continue;
 183
 184     Float_t resRphiRandom=resRphi*trackRes;
 185
 186     Int_t binZ   = TMath::BinarySearch(vZ->GetNrows(),vZ->GetMatrixArray(),(Double_t)trackZ);
 187     Int_t binPhi = TMath::BinarySearch(vPhi->GetNrows(),vPhi->GetMatrixArray(),(Double_t)trackPhi);
 188     Int_t binR   = TMath::BinarySearch(vR->GetNrows(),vR->GetMatrixArray(),(Double_t)radius);
 189
 190     if (binZ<0)   binZ=0;
 191     if (binPhi<0) binPhi=0;
 192     if (binR<0)   binR=0;
 193
 194     TObjArray *arrPhi=(TObjArray*)arrZ.UncheckedAt(binZ);
 195     if (!arrPhi) {
 196       arrPhi=new TObjArray(vPhi->GetNrows()-1);
 197       arrZ.AddAt(arrPhi,binZ);
 198     }
 199
 200     TObjArray *arrR=(TObjArray*)arrPhi->UncheckedAt(binPhi);
 201     if (!arrR) {
 202       arrR=new TObjArray(vR->GetNrows()-1);
 203       arrPhi->AddAt(arrR,binPhi);
 204     }
 205
 206     TH1S *h = (TH1S*)arrR->UncheckedAt(binR);
 207     if (!h) {
 208       h = new TH1S(Form("h_%02d_%02d_%d02",binZ, binPhi, binR),
 209                    Form("z,phi,r: %02d,%02d,%d02; #Delta r#phi (cm)",binZ, binPhi, binR),
 210                    80, -2., 2.);
 211       arrR->AddAt(h, binR);
 212     }
 213
 214     h->Fill(trackY+resRphiRandom-pointY);
 215   }
 216
 217   // do fits and fill tree
 218 }
 219
 220 void AnaDeltaResiduals(TString fluctuationMap, TString averageMap, TString outFile="deltas_residuals.root")
 221 {
 222   //
 223   //
 224   //
 225
 226   TFile fFluct(fluctuationMap);
 227   AliTPCCorrectionLookupTable *corrFluct = (AliTPCCorrectionLookupTable*)fFluct.Get("map");
 228   fFluct.Close();
 229
 230   TFile fAverage(averageMap);
 231   AliTPCCorrectionLookupTable *corrAverage = (AliTPCCorrectionLookupTable*)fAverage.Get("map");
 232   fAverage.Close();
 233
 234 //   TObjArray *arrMaps = new TObjArray(2);
 235 //   arrMaps->Add(corrAverage); // correction with the average Map
 236 //   arrMaps->Add(corrFluct);   // distortion with the fluctuation Map
 237
 238   // create the composed correction
 239   // if the weight are set to +1 and -1, the first map will be responsible for the distortions
 240   // The second map for the corrections
 241   // !!!!! In AliTPCComposedCorrection::GetDistortion MakeInverseIterator is called !!!!
 242   // for this reason we have to add the maps in the wrong order
 243
 244 //   AliTPCComposedCorrection *residualDistortion = new AliTPCComposedCorrection(arrMaps, AliTPCComposedCorrection::kQueueResidual);
 245   Float_t dummy=0;
 246 //   TVectorD weights(2);
 247 //   weights(0)=+1.;
 248 //   weights(1)=-AliToyMCEventGenerator::GetSCScalingFactor(corrFluct, corrAverage,dummy);
 249 //   residualDistortion->SetWeights(&weights);
 250
 251   corrAverage->SetCorrScaleFactor(AliToyMCEventGenerator::GetSCScalingFactor(corrFluct, corrAverage,dummy));
 252
 253   TVectorD *vR   = MakeLinBinning(10,86.,250.);
 254   TVectorD *vPhi = MakeLinBinning(18*8,0.,2*TMath::Pi());
 255   TVectorD *vZ   = MakeLinBinning(50,-250.,250.);
 256
 257   const Int_t nbins=4;
 258   Int_t bins[nbins]    = {vR->GetNrows()-1, vPhi->GetNrows()-1, vZ->GetNrows()-1, 80};
 259   //   Int_t bins[nbins]    = {16, 18*5, 50, 80};
 260   Double_t xmin[nbins] = {86. , 0.,           -250., -2.};
 261   Double_t xmax[nbins] = {250., 2*TMath::Pi(), 250.,  2.};
 262   THnF *hn = new THnF("hn", "hn", nbins, bins, xmin, xmax);
 263
 264   hn->GetAxis(0)->Set(vR  ->GetNrows()-1, vR  ->GetMatrixArray());
 265   hn->GetAxis(1)->Set(vPhi->GetNrows()-1, vPhi->GetMatrixArray());
 266   hn->GetAxis(2)->Set(vZ  ->GetNrows()-1, vZ  ->GetMatrixArray());
 267
 268   AliExternalTrackParam vv;
 269
 270   for (Float_t iz=-245; iz<=245; iz+=2) {
 271     Short_t roc=(iz>=0)?0:18;
 272     for (Float_t ir=86; ir<250; ir+=1) {
 273       for (Float_t iphi=0; iphi<TMath::TwoPi(); iphi+=0.5*TMath::DegToRad()){
 274         Float_t x=ir*(Float_t)TMath::Cos(iphi);
 275         Float_t y=ir*(Float_t)TMath::Sin(iphi);
 276         Float_t x3[3]    = {x,y,iz};
 277         Float_t x3dc[3]    = {x,y,iz};
 278         Float_t dx3[3]   = {0.,0.,0.};
 279 //         residualDistortion->GetDistortion(x3,roc,dx3);
 280         corrFluct->DistortPoint(x3dc,roc);
 281         corrAverage->CorrectPoint(x3dc,roc);
 282         dx3[0]=x3dc[0]-x3[0];
 283         dx3[1]=x3dc[1]-x3[1];
 284         dx3[2]=x3dc[2]-x3[2];
 285
 286         Double_t ddx3[3]={dx3[0], dx3[1], dx3[2]};
 287         vv.Global2LocalPosition(ddx3,iphi);
 288
 289         Double_t xx[4]={ir, iphi, iz ,ddx3[1]};
 290         hn->Fill(xx);
 291
 292       }
 293     }
 294   }
 295
 296   TTreeSRedirector stream(outFile.Data());
 297   gROOT->cd();
 298
 299   DumpHn(hn, stream);
 300
 301   stream.GetFile()->cd();
 302   hn->Write();
 303
 304   delete hn;
 305   delete vR;
 306   delete vPhi;
 307   delete vZ;
 308
 309 //   delete residualDistortion;
 310 }
 311
 312 void DumpHn(THn *hn, TTreeSRedirector &stream)
 313 {
 314   TAxis *ar   = hn->GetAxis(0);
 315   TAxis *aphi = hn->GetAxis(1);
 316   TAxis *az   = hn->GetAxis(2);
 317
 318   // output Hn
 319   const Int_t nbins=3;
 320   Int_t bins[nbins]    = {1,1,1};
 321   Double_t xmin[nbins] = {0.,0.,0.};
 322   Double_t xmax[nbins] = {1.,1.,1.};
 323   THnF hnRes("hnRes", "hnRes", nbins, bins, xmin, xmax);
 324
 325   ar  ->Copy(*hnRes.GetAxis(0));
 326   aphi->Copy(*hnRes.GetAxis(1));
 327   az  ->Copy(*hnRes.GetAxis(2));
 328
 329
 330   for (Int_t iz=0; iz<az->GetNbins(); ++iz) {
 331     az->SetRange(iz+1,iz+1);
 332     TObjArray arrFits;
 333     arrFits.SetName(Form("z_%02d",iz));
 334     arrFits.SetOwner();
 335
 336     for (Int_t ir=0; ir<ar->GetNbins(); ++ir) {
 337       ar->SetRange(ir+1,ir+1);
 338       for (Int_t iphi=0; iphi<aphi->GetNbins(); ++iphi) {
 339         aphi->SetRange(iphi+1,iphi+1);
 340
 341         Float_t cr       = 0.;
 342         Float_t cphi     = 0.;
 343         Float_t cz       = 0.;
 344         Float_t mean     = 0.;
 345         Float_t meanErr  = 0.;
 346         Float_t sigma    = 0.;
 347         Float_t sigmaErr = 0.;
 348         Int_t   entries  = 0.;
 349         Float_t chi2ndf  = 0.;
 350         Float_t mean2    = 0.;
 351         Float_t meanErr2 = 0.;
 352         Float_t rms2     = 0.;
 353         Float_t rmsErr2  = 0.;
 354
 355         TH1 *hProj = hn->Projection(3);
 356         if (hProj->GetEntries()>1) {
 357           TF1 fg("fg","gaus",-2,2);
 358           cr   = ar->GetBinCenter(ir+1);
 359           cphi = aphi->GetBinCenter(iphi+1);
 360           cz   = az->GetBinCenter(iz+1);
 361           hProj->SetNameTitle(Form("h_%02d_%02d_%02d",iz, iphi, ir),
 362           Form("z,phi,r: %02d,%02d,%02d (%.2f, %.2f, %.2f)",iz,iphi,ir, cz, cphi, cr )
 363           );
 364           hProj->Fit(&fg,"LMQR");
 365           arrFits.Add(hProj);
 366
 367           mean     = fg.GetParameter(1);
 368           meanErr  = fg.GetParError(1);
 369           sigma    = fg.GetParameter(2);
 370           sigmaErr = fg.GetParError(2);
 371           entries  = hProj->GetEntries();
 372           chi2ndf  = fg.GetChisquare()/fg.GetNDF();
 373           mean2    = hProj->GetMean();
 374           meanErr2 = hProj->GetMeanError();
 375           rms2     = hProj->GetRMS();
 376           rmsErr2  = hProj->GetRMSError();
 377         } else {
 378           delete hProj;
 379         }
 380
 381         stream << "d" <<
 382         "ir="          << ir       <<
 383         "iphi="        << iphi     <<
 384         "iz="          << iz       <<
 385         "cr="          << cr       <<
 386         "cphi="        << cphi     <<
 387         "cz="          << cz       <<
 388         "mean="        << mean     <<
 389         "meanErr="     << meanErr  <<
 390         "sigma="       << sigma    <<
 391         "sigmaErr="    << sigmaErr <<
 392         "histMean="    << mean2    <<
 393         "histMeanErr=" << meanErr2 <<
 394         "histRMS="     << rms2    <<
 395         "histRMSErr="  << rmsErr2 <<
 396         "entries="     << entries  <<
 397         "chi2ndf="     << chi2ndf  <<
 398         "\n";
 399
 400 //        Double_t x[nbins]={cr, cphi, cz};
 401 //        if (meanErr<0.3) hnRes.Fill(x,mean);
 402       }
 403     }
 404     stream.GetFile()->cd();
 405     arrFits.Write(0x0,TObject::kSingleKey);
 406     gROOT->cd();
 407   }
 408
 409   stream.GetFile()->cd();
 410   hnRes.Write();
 411   gROOT->cd();
 412 }
 413
 414 //______________________________________________________________________________
 415 TVectorD* MakeLogBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
 416 {
 417   //
 418   // Make logarithmic binning
 419   // the user has to delete the array afterwards!!!
 420   //
 421
 422   //check limits
 423   if (xmin<1e-20 || xmax<1e-20){
 424     printf("For Log binning xmin and xmax must be > 1e-20. Using linear binning instead!");
 425     return MakeLinBinning(nbinsX, xmin, xmax);
 426   }
 427   if (xmax<xmin){
 428     Double_t tmp=xmin;
 429     xmin=xmax;
 430     xmax=tmp;
 431   }
 432   TVectorD *binLim=new TVectorD(nbinsX+1);
 433   Double_t first=xmin;
 434   Double_t last=xmax;
 435   Double_t expMax=TMath::Log(last/first);
 436   for (Int_t i=0; i<nbinsX+1; ++i){
 437     (*binLim)[i]=first*TMath::Exp(expMax/nbinsX*(Double_t)i);
 438   }
 439   return binLim;
 440 }
 441
 442 //______________________________________________________________________________
 443 TVectorD* MakeLinBinning(Int_t nbinsX, Double_t xmin, Double_t xmax)
 444 {
 445   //
 446   // Make linear binning
 447   // the user has to delete the array afterwards!!!
 448   //
 449   if (xmax<xmin){
 450     Double_t tmp=xmin;
 451     xmin=xmax;
 452     xmax=tmp;
 453   }
 454   TVectorD *binLim=new TVectorD(nbinsX+1);
 455   Double_t first=xmin;
 456   Double_t last=xmax;
 457   Double_t binWidth=(last-first)/nbinsX;
 458   for (Int_t i=0; i<nbinsX+1; ++i){
 459     (*binLim)[i]=first+binWidth*(Double_t)i;
 460   }
 461   return binLim;
 462 }
 463
 464 //_____________________________________________________________________________
 465 TVectorD* MakeArbitraryBinning(const char* bins)
 466 {
 467   //
 468   // Make arbitrary binning, bins separated by a ','
 469   //
 470   TString limits(bins);
 471   if (limits.IsNull()){
 472     printf("Bin Limit string is empty, cannot add the variable");
 473     return 0x0;
 474   }
 475
 476   TObjArray *arr=limits.Tokenize(",");
 477   Int_t nLimits=arr->GetEntries();
 478   if (nLimits<2){
 479     printf("Need at leas 2 bin limits, cannot add the variable");
 480     delete arr;
 481     return 0x0;
 482   }
 483
 484   TVectorD *binLimits=new TVectorD(nLimits);
 485   for (Int_t iLim=0; iLim<nLimits; ++iLim){
 486     (*binLimits)[iLim]=(static_cast<TObjString*>(arr->At(iLim)))->GetString().Atof();
 487   }
 488
 489   delete arr;
 490   return binLimits;
 491 }
 492
 493 void PlotFromTree(TTree *d, TString outDir=".")
 494 {
 495   TCanvas *c=new TCanvas;
 496   gStyle->SetOptStat(0);
 497   d->SetMarkerStyle(20);
 498   d->SetMarkerSize(1);
 499
 500   TProfile2D pRZ("pRZ",";z (cm); r(cm)",50,-250,250,10,85,250);
 501   d->Draw("entries:cr:cz>>pRZ","","profcolz");
 502   pRZ.GetZaxis()->UnZoom();
 503   c->SaveAs(Form("%s/entries_average.png",outDir.Data()));
 504   d->Draw("entries:cr:cz>>pRZ","iphi==2","profcolz");
 505   c->SaveAs(Form("%s/entries_onePhi.png",outDir.Data()));
 506
 507   pRZ.SetMaximum(0.04);
 508   d->Draw("meanErr:cr:cz>>pRZ","","profcolz");
 509   c->SaveAs(Form("%s/meanErr_average.png",outDir.Data()));
 510   d->Draw("meanErr:cr:cz>>pRZ","iphi==2","profcolz");
 511   c->SaveAs(Form("%s/meanErr_onePhi.png",outDir.Data()));
 512
 513
 514   d->Draw("mean:cphi:cr","iz==25","colz");
 515   c->SaveAs(Form("%s/mean_oneZ_phi_allR.png",outDir.Data()));
 516   d->Draw("mean:meanErr:cphi","iz==25&&ir==2","goff");
 517   TGraphErrors *grmean_phi=new TGraphErrors(d->GetSelectedRows(),d->GetV3(),d->GetV1(),0,d->GetV2());
 518   grmean_phi->SetTitle(";#varphi;#LT#Delta r#varphi#GT");
 519   grmean_phi->SetMarkerStyle(20);
 520   grmean_phi->SetMarkerSize(1);
 521   grmean_phi->Draw("ap");
 522   c->SaveAs(Form("%s/mean_oneZ_phi_oneR.png",outDir.Data()));
 523
 524   d->Draw("mean:cr:cphi","iz==25","colz");
 525   c->SaveAs(Form("%s/mean_oneZ_r_allPhi.png",outDir.Data()));
 526
 527   d->Draw("mean:meanErr:cr","iz==25&&iphi==2","goff");
 528   TGraphErrors *grmean_r=new TGraphErrors(d->GetSelectedRows(),d->GetV3(),d->GetV1(),0,d->GetV2());
 529   grmean_r->SetTitle(";r (cm);#LT#Delta r#varphi#GT");
 530   grmean_r->SetMarkerStyle(20);
 531   grmean_r->SetMarkerSize(1);
 532   grmean_r->Draw("ap");
 533   c->SaveAs(Form("%s/mean_oneZ_r_onePhi.png",outDir.Data()));
 534
 535
 536   d->Draw("meanErr:cphi:cr","iz==25","colz");
 537   c->SaveAs(Form("%s/meanErr_oneZ_phi_allR.png",outDir.Data()));
 538   d->Draw("meanErr:cphi","iz==25&&ir==2");
 539   c->SaveAs(Form("%s/meanErr_oneZ_phi_oneR.png",outDir.Data()));
 540
 541   d->Draw("meanErr:cr:cphi","iz==25","colz");
 542   c->SaveAs(Form("%s/meanErr_oneZ_r_allPhi.png",outDir.Data()));
 543
 544   d->Draw("meanErr:cr","iz==25&&iphi==2");
 545   c->SaveAs(Form("%s/meanErr_oneZ_r_onePhi.png",outDir.Data()));
 546
 547 }
 548
 549