Subversion Repositories f9daq

Rev

Blame | Last modification | View Log | RSS feed

  1. typedef struct  {
  2.    
  3.     unsigned short offset;
  4.     unsigned short fineadj_unipol;
  5.     unsigned short kill;
  6.    
  7. } channel_register;
  8.  
  9. channel_register creg[144];
  10. double gsigma=0;
  11. double gmean=0;
  12.  
  13.  
  14. double sa02calculateoffset(double x, int *offset, int * fineadj_unipol ){
  15.  
  16.     const double fscale= 2500/1023;      // 2500 mV / 1023 ch / dacch
  17.     const double foffset = 59.0/fscale;  // 59.0 mV / dacch (coarse)
  18.     const double ffine    = 5.4/fscale; // 5.4mV / dacch (fine)
  19.      
  20.     *offset        =  int(x/ foffset);
  21.     x=x-(*offset)*foffset;
  22.     *fineadj_unipol= int(x/ffine);
  23.         x=x-(*fineadj_unipol)*ffine;
  24.         return x;
  25.  
  26. }
  27.  
  28. int sa02fithisto(double mean0=512 , double sigma0=3){
  29. gmean=mean0;
  30. gsigma=sigma0;
  31.  
  32.  TH2F *h=NULL;
  33.  if (!h) h=((TH2F * ) gDirectory->Get("h2d"));
  34. double max= h->GetMaximum();
  35. TF1 *f = new TF1("merfc","[2]*TMath::Erfc((x-[0])*[1])", 0 , 1000);
  36. f->FixParameter(2,max*0.5);
  37.  
  38.  TNtuple *nt=new TNtuple("nt","nt","ch:mean:sigma");
  39.  
  40.  
  41.  h->FitSlicesY(f);
  42.  TAxis *ax= h->GetYaxis();
  43.  double xmin=  ax->GetXmin();
  44.  double xmax=  ax->GetXmax();
  45.  char n0[100];
  46.  char n1[100];
  47.  char n2[100];
  48.  sprintf(n0,"%s_0",h->GetName());
  49.  sprintf(n1,"%s_1",h->GetName());
  50.  sprintf(n2,"%s_2",h->GetName());
  51.  TH1F *h0 = ((TH1F * ) gDirectory->Get(n0));
  52.  TH1F *h1 = ((TH1F * ) gDirectory->Get(n1));
  53.  TH1F *h2 = ((TH1F * ) gDirectory->Get(n2));
  54.  //h0->Draw();
  55.  TCanvas *c= new TCanvas("c","c",1000,1000);
  56.  c->Divide(2,2);
  57.  
  58.  c->cd(1);
  59.  h0->SetMaximum(xmax);
  60.  h0->SetMinimum(xmin);
  61.  h0->Draw();
  62.  c->cd(3);
  63.  h1->SetMaximum(200);
  64.  h1->SetMinimum(0);
  65.  sprintf(n0,"%sc",h1->GetName());
  66.  TH1F *h1c = h1->Clone(n0);
  67.  
  68.  float fx[10];
  69.  for (int i=1;i<h1->GetNbinsX()+1;i++){
  70.     double first = h->GetBinContent(h->FindBin(i-1,xmin+1));
  71.     double last  = h->GetBinContent(h->FindBin(i-1,xmax-1));
  72.    
  73.     double sigma= TMath::Abs( h1->GetBinContent(i) );
  74.     double mean = TMath::Abs( h0->GetBinContent(i) );
  75.     if (sigma>0) sigma=1/sigma; else sigma=0;
  76.     fx[0]=i;
  77.     fx[1]=mean;
  78.     fx[2]=sigma;
  79.         nt->Fill(fx);
  80.         double            x = ((mean+sigma0*sigma)-mean0) ;
  81.     int offset        = 0;
  82.     int fineadj_unipol= 0;
  83.  
  84.         sa02calculateoffset(x, &offset, &fineadj_unipol );
  85.    
  86.    
  87.        
  88.        
  89.         int valid = 1;
  90.         if (last>0.5*h2d->GetMaximum())  valid = 0;
  91.         if (first<0.5*h2d->GetMaximum()) valid = 0;
  92.         if (valid) {
  93.           creg[i-1].offset         = offset;
  94.           creg[i-1].fineadj_unipol = fineadj_unipol;
  95.           creg[i-1].kill           = 0;
  96.         } else {
  97.           creg[i-1].offset         = 0;
  98.           creg[i-1].fineadj_unipol = 0;
  99.           creg[i-1].kill           = 1;
  100.         }
  101.    
  102.     printf("%d %f %f %f %f %d %d %f\n",i, first, last, mean,sigma,offset,fineadj_unipol ,x );
  103.    
  104.     h1c->SetBinContent(i,sigma);
  105.      
  106.     // of << ch << "\t" << (int)entries << "\t" << (((mean*2.4)-1245.2)) << "\t" << sigma*2.4 << endl;
  107.    
  108.  }
  109.  h1c->Draw();
  110.  c->cd(4);
  111.  nt->Draw("sigma","sigma>0&&sigma<100");
  112.  //h2->Draw();
  113.  c->cd(2);
  114.  
  115.  //h->Draw("colz");
  116.  /*
  117.  h0->SetLineColor(kRed);
  118.  h0->SetLineWidth(10);
  119.  h0->Draw("same");
  120.  */
  121.  nt->Draw("mean","mean>0&&mean<1000");
  122.  
  123.  c->Modified();
  124.  c->Update();
  125.  
  126. }
  127.  
  128. int sa02param( char *fname="../params/offset_orig.param", char *hname="offsettmp.param" ){
  129.  
  130.    
  131.     int ndim=400;
  132.     char line[400];
  133.     char cmd[400];
  134.     char sasic[400];
  135.     char v0[400];
  136.     char v1[400];
  137.     int asic=0, ch=0;
  138.     int  cval=0;
  139.     int i;
  140.  
  141.     FILE *fout = fopen(hname,"w");
  142.     if (!fout) {
  143.         printf("Error! Cannot open file %s\n",hname);
  144.         return -1;
  145.     }
  146.    
  147.    
  148.     FILE *fp = fopen(fname,"r");
  149.     if (!fp) {
  150.         printf("Error! Cannot open file %s\n",fname);
  151.         return -1;
  152.     }
  153.    
  154.         int write=1;
  155.         fprintf(fout,"### Created from sa02fit mean=%f sigma=%f\n", gmean,gsigma);
  156.         fprintf(fout,"### date tbi\n");
  157.         char *null=NULL;
  158.     while (fgets(line,ndim,fp)!=NULL) {
  159.                 write=1;
  160.         int nb = sscanf(line,"%s%s%s%s",cmd,sasic,v0,v1);
  161.         if (nb<1 || cmd[0]=='#') {
  162.            fprintf(fout,line);
  163.            continue;
  164.         }  
  165.  
  166.         asic =   atoi (sasic);
  167.         ch   =   atoi (v0);
  168.         cval =   strtoul (v1,null,0);
  169.         if (strcmp(cmd,"offset")==0)     {
  170.                 fprintf(fout,"offset %d %d %d\n", asic,ch, creg[asic*36+ch].offset);
  171.                 write=0;
  172.         }
  173.         if (strcmp(cmd,"fineadj_unipol")==0) {
  174.                 fprintf(fout,"fineadj_unipol %d %d %d\n", asic,ch, creg[asic*36+ch].fineadj_unipol);
  175.                 write=0;
  176.         }
  177.         if (strcmp(cmd,"kill")==0)       {
  178.                 fprintf(fout,"kill %d %d %d\n", asic,ch, creg[asic*36+ch].kill);
  179.                 write=0;
  180.         }
  181.         if (write){
  182.            fprintf(fout,line);
  183.         }
  184.        
  185.     }
  186.     fclose(fp);
  187.     fclose(fout);
  188. return 0;
  189. }
  190.  
  191. int sa02fit(char *hname,double mean0=400 , double sigma0=3, char *fnameparam="../../params/offset_orig.param"){
  192.  gROOT->SetMacroPath("../macros");
  193.  gROOT->ProcessLine(".L H2Dload.cxx");
  194.  
  195.  H2Dload(hname);
  196.  //H2Dload(gSystem->UnixPathName(hname));
  197.  sa02fithisto(mean0,sigma0);
  198.  char foutparname[400];
  199.  sprintf(foutparname,"%s.param",hname);
  200.  sa02param(fnameparam,foutparname);
  201. return 0;      
  202. }
  203.