Changeset 362 in PSPA for Interface_Web

Timestamp:

Mar 5, 2013, 11:48:58 AM (12 years ago)

Author:

touze

Message:

histogramme en x, y ou z

Location:

Interface_Web/trunk/pspaWT/sources

Files:

• controler/include/particleBeam.h (modified) (2 diffs)
• controler/src/particleBeam.cc (modified) (1 diff)
• userInterface/include/GWt_pspaApplication.h (modified) (1 diff)
• userInterface/src/GWt_pspaApplication.cc (modified) (6 diffs)

Interface_Web/trunk/pspaWT/sources/controler/include/particleBeam.h

@@ -43 +43 @@
   unsigned indexFromPspaToTransport(unsigned index) const;
 
+  void preparationDesHistogrammes(int index,vector<double>& vshf,double out[3]);
+  void remplissageDesHistogrammes(double out[3],vector<double>vshf,vector<double>&xcor,vector<int>& hist,int& cnts);
+
  public:
 
@@ -73 +76 @@
   virtual bool FileInput(ifstream& ifs);
   void histogramme(vector<double>&xcor,vector<int>& hist,int& cnts,double out[3]);
+  void histogramme(int index,vector<double>&xcor,vector<int>& hist,int& cnts,double out[3]);
 };
 #endif

Interface_Web/trunk/pspaWT/sources/controler/src/particleBeam.cc

@@ -510 +510 @@
 }
 
+void particleBeam::histogramme(int iabs,vector<double>&xcor,vector<int>& hist,int& cnts,double out[3])
+{ 
+  vector<double> vshf;
+  preparationDesHistogrammes(iabs,vshf,out);
+  remplissageDesHistogrammes(out,vshf,xcor,hist,cnts);
+
+  // sortie pour la legende: out[0]= entries, out[1]= mean, out[2]= std deviation
+  if(iabs < 3) {
+    out[1] *= 1.0 ; // en m ?
+    out[2] *= 1000.0; //mm
+  } else {
+    out[1] *= EREST_MeV ; //MeV
+    out[2] *= EREST_keV; //KeV
+  }
+}
+
+void particleBeam::preparationDesHistogrammes(int index,vector<double>& vshf,double out[3])
+{
+  double vmin= GRAND;
+  double vmax= -vmin;
+  double vmoy= 0.0;
+  double ecatyp= 0.0;
+
+  double sum= (float)goodPartic_.size();
+  if(index < 3) 
+    {
+      // collecte d'une coordonnee x(index= 0), y(index= 1) ou z(index= 2) des particules
+      for (unsigned int k = 0; k < goodPartic_.size(); k++) {
+        double val = goodPartic_.at(k).getPosition().getComponent(index);
+        if (val < vmin) vmin = val;
+        else if (val > vmax) vmax = val;
+        vmoy += val;
+        ecatyp += val*val;
+      }
+
+      out[0]= sum;
+      vmoy /= sum;
+      out[1]= vmoy;
+      ecatyp /= sum;
+      ecatyp = sqrt(abs(ecatyp-vmoy*vmoy));
+      out[2]= ecatyp;
+    } 
+  else 
+    {
+      // collecte de l'energie acquise par les particules
+      for (unsigned int k = 0; k < goodPartic_.size(); k++) {
+        double gamma = goodPartic_.at(k).getGamma();
+        if (gamma < vmin) vmin = gamma;
+        else if (gamma > vmax) vmax = gamma;
+        double val = gamma-1.0;
+        vmoy += val;
+        ecatyp += val*val;
+      }
+      
+      out[0]= sum;
+      vmoy /= sum;
+      out[1]= vmoy;
+      ecatyp /= sum;
+      ecatyp = sqrt(abs(ecatyp-vmoy*vmoy));
+      out[2]= ecatyp;
+    }
+  
+  if(index == 0) {
+    cout << "position X -moyenne " << vmoy << " m " << "-mini " << vmin << " m " << "-maxi " << vmax << " m " << "ecart type " << ecatyp*1000.0 << " mm" << endl;
+  } else if(index == 1) {
+    cout << "position Y -moyenne " << vmoy << " m " << "-mini " << vmin << " m " << "-maxi " << vmax << " m " << "ecart type " << ecatyp*1000.0 << " mm" << endl;
+  } else if(index == 2) {
+    cout << "position Z -moyenne " << vmoy << " m " << "-mini " << vmin << " m " << "-maxi " << vmax << " m " << "ecart type " << ecatyp*1000.0 << " mm" << endl;
+  } else {
+    double gmin = vmin-1.0;
+    double gmax = vmax-1.0;
+    cout << "energie cinetique -moyenne " << vmoy*EREST_MeV << " Mev " << "-mini " << gmin*EREST_MeV << " Mev " << "-maxi " << gmax*EREST_MeV << " Mev " << "ecart type " << ecatyp*EREST_MeV << " Kev" << endl;
+  }
+
+  //pouVoir: autre maniere de calculer l'ecart-type
+  double stdev= 0.0;
+
+  if(!vshf.empty()) vshf.clear();
+  vshf.resize(goodPartic_.size());  
+  if(index < 3) 
+    {
+      for (unsigned int k = 0; k < goodPartic_.size(); k++) {
+        double val = goodPartic_.at(k).getPosition().getComponent(index);
+        vshf[k]= val-vmoy;
+        stdev += vshf[k]*vshf[k];
+      }
+    } 
+  else 
+    {
+      for (unsigned int k = 0; k < goodPartic_.size(); k++) {
+        double gamma = goodPartic_.at(k).getGamma();
+        double val = gamma-1.0;
+        vshf[k]= val-vmoy;
+        stdev += vshf[k]*vshf[k];
+      }
+    }
+
+  stdev /= sum;
+  stdev = sqrt(stdev);
+  cout << "stdev " << stdev << ", ecatyp " << ecatyp << ", diff " << stdev-ecatyp << endl;
+}
+
+void particleBeam::remplissageDesHistogrammes(double out[3],vector<double>vshf,vector<double>&xcor,vector<int>& hist,int& cnts)
+{
+  // recupere la moyenne et l'ecart-type
+  double vmoy= out[1];
+  double ecatyp= out[2];
+
+  // demi fenetre hfene= max(3.*ecatyp-vmoy,vmoy-3.*ecatyp);
+  double hfene= 3.*ecatyp;
+  // et pas de l'histogramme
+  double hpas = hfene/25.;
+  
+  cout << "demi fenetre " << hfene << ", hpas= " << hpas << endl;
+
+  double vmin = -hfene;
+  double dfen = 2.*hfene;
+  int ihist = dfen/hpas;
+  double phist = ihist*hpas;
+  double dpas = hpas-(dfen-phist);
+  if(dpas <= hpas*1.e-03) {
+    ihist++;
+    phist= ihist*hpas;
+  }
+  double vmax= vmin+hpas*ihist;
+  
+  cout << "'xAxisNumberOfBins= " << ihist <<", xAxisMinimum= " << vmin << ", xAxisMaximum= " << vmax << ", NParticules= " << vshf.size() << ", phist " << phist << endl;
+
+  if(!xcor.empty()) xcor.clear();
+  xcor.resize(ihist);  
+  for (int i = 0; i < ihist; ++i) {
+
+    // on gradue l'abscisse en pourcents
+    xcor[i]= 100.*( vmin+i*hpas )/vmoy;
+  }
+
+  if(!hist.empty()) hist.clear();
+  hist.resize(ihist,0);  
+  for (unsigned i = 0; i < vshf.size(); ++i) {
+    double var= vshf[i]-vmin;
+    if(var < 0 ) {
+      cout<<"lesser that the minimum "<<var<<", ("<< i<<")"<< endl;
+      hist[0]++;
+    }
+
+    if(var >= phist) {
+      cout<<"greater that the maximum "<<var<<", ("<< i<<")"<< endl;
+      hist[ihist-1]++;
+    }
+
+    int kk= (int)floor(var/hpas);
+    hist[kk]++;
+  }
+
+  cnts= 0;
+  for (int i = 0; i < ihist; ++i) {
+    if(hist.at(i) > 0) cnts++;
+  }
+  cout<< " ... cnts=  " << cnts << endl;
+}
+
 void particleBeam::histogramme(vector<double>&xcor,vector<int>& hist,int& cnts,double out[3])
 {

Interface_Web/trunk/pspaWT/sources/userInterface/include/GWt_pspaApplication.h

@@ -53 +53 @@
   WComboBox* choixEnveloppeDessin_;
   WComboBox* choixHistoDessin_;
+  WComboBox* choixVariableHisto_;
 
   WContainerWidget* toto_;

Interface_Web/trunk/pspaWT/sources/userInterface/src/GWt_pspaApplication.cc

@@ -273 +273 @@
   /////////////////////////////////////////////////////////////////////
 
+  choixVariableHisto_= new WComboBox();
+  choixVariableHisto_->addItem("x");
+  choixVariableHisto_->addItem("y");
+  choixVariableHisto_->addItem("z");
+  choixVariableHisto_->addItem("dE/E");
+  choixVariableHisto_->setCurrentIndex(3);
+
   choixHistoDessin_ = new WComboBox();
-  Wt::WRadioButton *button3= new WRadioButton(" histogram after element");
+  Wt::WRadioButton *button3= new WRadioButton(" histogram of");
   group_->addButton(button3,3);
+
   glayout->addWidget(button3,3,1);
-  glayout->addWidget(choixHistoDessin_,3,2);
+  glayout->addWidget(choixVariableHisto_,3,2);
+  glayout->addWidget(new WText("after element"),3,3);
+  glayout->addWidget(choixHistoDessin_,3,4);
   /////////////////////////////////////////////////////////////////////
   
@@ -805 +815 @@
   histoDialog->setMinimumSize(400,400);
   histoDialog->setClosable(true);
-  histoDialog->show();
   
   int index = choixHistoDessin_->currentIndex();    
@@ -821 +830 @@
   }
 
+  histoDialog->show();
   vector<double> xcor;
   vector<int> hist;
   int cnts;
   double out[3]= {0.0};
-  beam->histogramme(xcor,hist,cnts,out);
-  
+
+  int iabs= choixVariableHisto_->currentIndex();  
+  string xname= choixVariableHisto_->currentText().toUTF8();
+  beam->histogramme(iabs,xcor,hist,cnts,out);
+
   cout<<"xcor.size()= "<<xcor.size()<<", hist.size()= "<<hist.size()<<endl;
   //////////////////////////////////////////////////////////////////////////////////
@@ -833 +846 @@
   int nbpts= xcor.size()+cnts+1;
   WStandardItemModel *model = new WStandardItemModel(nbpts,2,w);
-  model->setHeaderData(1, WString("energie"));
+  //model->setHeaderData(1, WString("energie"));
   
   int j= 0;
@@ -868 +881 @@
   new WText(" entries : "+ mixedTools::intToString((int)out[0]),w);
   new WBreak(w);
-  new WText(" mean : "+ mixedTools::doubleToString(out[1])+" MeV",w);
-  new WBreak(w);
-  new WText(" sigma rms : "+ mixedTools::doubleToString(out[2])+" KeV",w);
-  
+  if(iabs < 3) {
+    new WText(" mean : "+ mixedTools::doubleToString(out[1])+" m",w);
+    new WBreak(w);
+    new WText(" sigma rms : "+ mixedTools::doubleToString(1000.0*out[2])+" mm",w);
+  } else {
+    new WText(" mean : "+ mixedTools::doubleToString(out[1])+" MeV",w);
+    new WBreak(w);
+    new WText(" sigma rms : "+ mixedTools::doubleToString(out[2])+" KeV",w);
+  }
+   
   WCartesianChart *chart = new WCartesianChart(w);
-  chart->setTitle("Histogram of kinetic energy");
+
+  if(iabs == 0) {
+    chart->setTitle("Histogram of particle x-coordinate");
+  } else if(iabs == 1) {
+    chart->setTitle("Histogram of particle y-coordinate");
+  } else if(iabs == 2) {
+    chart->setTitle("Histogram of particle z-coordinate");
+  } else {
+    chart->setTitle("Histogram of kinetic energy");
+  }
+
   chart->setModel(model);        // set the model
   chart->setXSeriesColumn(0);    // set the column that holds the categories
@@ -888 +917 @@
   axis.setLabelFormat("%.3f");
   axis.setGridLinesEnabled(true);
-  axis.setTitle(WString(" dEcin/Ecin (%)"));
+
+  if(iabs < 3) {
+    axis.setTitle(WString(" ??? (%)"));
+  } else {
+    axis.setTitle(WString(" dEcin/Ecin (%)"));
+  }
 
   chart->axis(Y1Axis).setGridLinesEnabled(true);
-  chart->axis(Y1Axis).setTitle(WString("legende y"));
+  //chart->axis(Y1Axis).setTitle(WString("legende y"));
     
   WDataSeries s(1, LineSeries);