// This may look like C code, but it is really -*- C++ -*- // Classe d'interface avec le code de GNUplot // pour calculer les contours & classe de trace de contours // O. Perdereau // LAL (Orsay) / IN2P3-CNRS DAPNIA/SPP (Saclay) / CEA #include "sopnamsp.h" #include "machdefs.h" #include #include #include #include #include #include "histos.h" #include "ntuple.h" #include "nbtri.h" #include "picntools.h" #include "pigncont.h" #include #include //++ // Class GNUPlotContour // Lib PIGcont // include pigncont.h // // Classe d'interface avec le code de GNUplot qui calcule les courbes de niveaux // GNUplot part d'une grille (structure iso_curve) qui est remplie par les // createurs // Les parametres du calcul sont modifiable a partir de la fenetre d'options // du traceur de contour (classe PIContourDrawer) // Options (reconnues dans un enum t_contour_kind dans GNUplot) // // 1) nombre de courbes et et niveaux determines automatiquement par le programme (defaut) // En general, 5 courbes sont tracees LEVELS_AUTO // // 2) nombre de niveaux fixe par l'utilisateur. Les niveaux sont equidistants entre // le min et le max de la surface LEVELS_NUM // // 3) nombre et hauteurs des niveaux fixes par l'utilisateur. Les niveaux sont passes // dans un tableau LEVELS_DISCRETE // // 4) nombre niveau min et ecart entre niveaux fixe par l'utilisateur LEVELS_INCREMENTAL // // Les contours sont traces en interpolant entre les points de la grille. // On peut choisir 3 algorithmes d'interpolation a travers p.ex. les options // de PIContourDrawer // (enum t_contour_levels_kind ds GNUplot) // // a) interpolation simple CONTOUR_KIND_LINEAR // // b) spline cubique CONTOUR_KIND_CUBIC_SPL // // c) B-spline CONTOUR_KIND_BSPLINE // //-- //++ // Links Voir aussi // PIContourDrawer // PICnTools //-- //++ // Titre Constructeurs et Méthodes //-- //++ // GNUPlotContour(struct iso_curve * iso) : // constructeur a partir d'une structure iso_curve // (structure de GNUplot) //-- GNUPlotContour::GNUPlotContour(struct iso_curve * iso){ int k=0; iso1 = iso; struct iso_curve *iso_cur = iso ; struct iso_curve *iso_old; _xmin = 1.e37; _ymin = 1.e37; _xmax = -1.e37; _ymax = -1.e37; _myiso = NULL; _zmin = 1.e37; _zmax = -1.e37; while(iso_cur){ k++; iso_cur = iso_cur->next; _ny = iso_cur->p_count; for(int k=0 ; k<_ny ; k++){ double xx=iso->points[k].x; double yy=iso->points[k].y; double zz=iso->points[k].z; if(xx<_xmin) _xmin = xx; if(xx>_xmax) _xmax = xx; if(yy<_ymin) _ymin = yy; if(yy>_ymax) _ymax = yy; if(zz<_zmin) _zmin = zz; if(zz>_zmax) _zmax = zz; } } _nx = k; _contours = NULL; My_Levels = NULL; _npolys=-1; _Dxp =_Dyp =1.; _Exy =_Invy =_Invy = FALSE ; } //++ // GNUPlotContour(NTupleInterface *nt,int *nz=NULL,double *z=NULL,double *dz=NULL): // constructeur a partir d'un NTupleInterface // la structure iso_curve est remplie par le constructeur //Inputs //| NTupleInterface *nt //| double *nz nombre de contours (optionnel) //| double *z tabelau avec les niveaux des contours (optionnel) //| double *dz increment entre contours (optionnel) //-- GNUPlotContour::GNUPlotContour(NTupleInterface *nt,int *nz=NULL,double *z=NULL,double *dz=NULL){ _xmin = 1.e37; _ymin = 1.e37; _xmax = -1.e37; _ymax = -1.e37; _zmin = 1.e37; _zmax = -1.e37; const char *nom[4]={"x","y","z","niso"}; _myiso = new NTuple(4,nom); //_myiso->Show(); float xnt[4]; double mn,mx; nt->GetMinMax(0,_xmin,_xmax); nt->GetMinMax(1,_ymin,_ymax); nt->GetMinMax(2,_zmin,_zmax); struct iso_curve *istmp = NULL ; struct iso_curve *iso = NULL ; struct iso_curve *iso_old = NULL ; double x,y; int nen = nt->NbLines(); // histo des y for(int k= 0 ; k<1000 ; k++){ istmp = iso_alloc(100); iso_free(istmp); } #define NUMPTMIN 25 #define MAXISO 100 #define MINISO 5 Histo prep(_ymin,_ymax,nen); for(int ix=0 ; ixGetCell(ix,1),1.); Histo Intg = prep; Intg.HInteg(0); //prep.Print(); //Intg.Print(); float *bin; int *cnt; int niso = sqrt( (double)nen); bin = new float[niso+1]; cnt = new int[niso+1]; int i; for(i=0 ; i0){ numx = numx>numy ? numx : numy ; numi = numy>numi ? numi : numy ; } } if (numi<=0) { return; } _nx = 0; for(i=0;iGetCell(ix,1)>=ymi&&nt->GetCell(ix,1)points[ip].type = INRANGE; istmp->points[ip].x = nt->GetCell(ix,0); istmp->points[ip].y = nt->GetCell(ix,1); istmp->points[ip].z = nt->GetCell(ix,2); indx[ip] = ip; xp[ip++] = nt->GetCell(ix,0); } } // tri des points de l'iso_curve // istmp->p_count = ip; tri_double(xp,indx,ip); // remplissage structure finale iso = iso_alloc(numx) ; if( _nx == 0 ) iso1 = iso; //premiere structure remplie if (iso_old) iso_old->next = iso; for(ix=0 ; ixpoints[ix].type = INRANGE ; iso->points[ix].x = istmp->points[jx].x ; iso->points[ix].y = istmp->points[jx].y ; iso->points[ix].z = istmp->points[jx].z ; xnt[0] = iso->points[ix].x; xnt[1] = iso->points[ix].y; xnt[2] = iso->points[ix].z; xnt[3] = _nx; //cout << " fill avec "<Show(); _myiso->Fill(xnt); } //cout << " completion "<> "<points[ix].type = INRANGE ; iso->points[ix].x = istmp->points[jmx].x ; iso->points[ix].y = istmp->points[jmx].y ; iso->points[ix].z = istmp->points[jmx].z ; xnt[0] = iso->points[ix].x; xnt[1] = iso->points[ix].y; xnt[2] = iso->points[ix].z; xnt[3] = _nx; //cout << " fill avec "<Fill(xnt); //_myiso->Show(); } iso->p_count = numx; iso_old = iso; _nx++; //cout << " destruction structure tempo "<<_nx <Write("myiso.ppf"); //cout << " fin du remplissage des iso_curves xmin,max "<< _xmin<<","<<_xmax<<" ymin,max "<<_ymin<<","<<_ymax<XSize()-1; // _ny = arr->YSize()-1; _nx = arr->XSize() ; _ny = arr->YSize() ; //cout << " size x,y "<<_nx<<" , "<<_ny <XYfromxy(0,0,x,y); arr->XYfromxy(1,1,x1,y1); _Dxp = (x1-x); _Dyp = (y1-y); double zz; for(int ix=0 ; ix<_nx ; ix++){ iso = iso_alloc(_ny) ; if(ix==0) iso1 = iso; if (iso_old) iso_old->next = iso; for(int iy = 0 ; iy<_ny ; iy++){ /*i = ix*100+iy;*/ arr->XYfromxy(ix,iy,x,y); iso->points[iy].type = INRANGE; iso->points[iy].x = x + _Dxp/2.; iso->points[iy].y = y + _Dyp/2.; zz=arr->Value(ix,iy); iso->points[iy].z = zz; //iso->points[iy].z = ((x-25)*(x-25) + (y-25)*(y-25))*.05 ; if(x <_xmin) _xmin = x; if(x >_xmax) _xmax = x; if(y <_ymin) _ymin = y; if(y >_ymax) _ymax = y; if(zz <_zmin) _zmin = zz; if(zz >_zmax) _zmax = zz; //printf("ix=%d iy=%d x=%f y=%f z=%f \n",ix,iy,iso->points[iy].x,iso->points[iy].y,iso->points[iy].z); xnt[0] = iso->points[iy].x; xnt[1] = iso->points[iy].y; xnt[2] = iso->points[iy].z; xnt[3] = _nx; _myiso->Fill(xnt); } iso->p_count = _ny; iso_old = iso; /*printf(" remplissage %lx old %lx old next %lx \n",iso,iso_old,iso->next,iso_old->next);*/ /*iso->next = iso ;*/ } //cout << " fin du remplissage des iso_curves xmin,max "<< _xmin<<","<<_xmax<<" ymin,max "<<_ymin<<","<<_ymax<next; iso_free(iso_cur);iso_cur = NULL; iso_cur = iso_nx; } iso1 = NULL; struct gnuplot_contours *cntcur = _contours; /***** struct gnuplot_contours *cntold; while(cntcur) { //cout << " ~GNUPlotContour() : destruction de _contours" << _contours <next; gp_free(cntold); cntold=NULL; } *******/ contour_free(cntcur); _contours = NULL; if(My_Levels) { //cout << " ~GNUPlotContour() : destruction de MyLevels "< "<< r_usage.ru_maxrss <<" , "<< r_usage.ru_ixrss <<" , "<< r_usage.ru_ixrss <next; gp_free(cntold); cntold = NULL; } ******/ contour_free(cntcur); _contours = NULL; } //struct gnuplot_contours *cntcur; //struct gnuplot_contours *cntold; //rcus = getrusage( RUSAGE_SELF , &r_usage); //if(rcus==0) // cout << " rusage -> "<< r_usage.ru_maxrss <<" , "<< r_usage.ru_ixrss <<" , "<< r_usage.ru_ixrss < "<< r_usage.ru_maxrss <<" , "<< r_usage.ru_ixrss <<" , "<< r_usage.ru_ixrss <vec) // Setting du tableau des niveaux des contours //Inputs //| //| vector vec : vecteur des niveaux //-- void GNUPlotContour::SetMyLevels(vector vec){ if(My_Levels!=NULL) delete My_Levels; int sz = vec.size(); assert (sz>0); My_Levels = new double[sz]; for(int i=0 ; iSaveGraphicAtt(); struct gnuplot_contours *cntcur; struct gnuplot_contours *cntold; cntcur = _contours; double a,b; int tot=0; bool Invx, Invy, Exy; PIColorMap * cmap = NULL; CMapId mcmapid = GetGraphicAtt().GetColMapId(); if(mcmapid !=CMAP_OTHER ){ cmap = new PIColorMap(mcmapid); } int numdr = 0; while(cntcur){ int npts = cntcur->num_pts ; if(npts<0) continue; xx = new PIGrCoord[npts]; yy = new PIGrCoord[npts]; double lev = (cntcur->coords)[0].z ; for(int l=0 ; l< npts ; l++){ a = (cntcur->coords)[l].x; b = (cntcur->coords)[l].y; xx[l] = a ; yy[l] = b ; } // choix de la couleur g->SelForeground(GetGraphicAtt().GetFgColor()); if(mcmapid !=CMAP_OTHER ){ if(_zmax-_zmin!=0){ int kc = ((lev - _zmin)/(_zmax-_zmin))*cmap->NCol(); g->SelForeground(*cmap,kc); } } // traitement des des markers if(IsMarkOn()==true){ g->SelMarker(GetGraphicAtt().GetMarkerSz(), GetGraphicAtt().GetMarker()); g->DrawMarkers(xx,yy,npts); } // traitement des lignes if(mLineOn==true){ g->SelLine(GetGraphicAtt().GetLineAtt()); g->DrawPolygon(xx,yy,npts,false); } // traitement des labels // SIMPLISTE POUR L'INSTANT : UN AFFICHAGE if(mLabelOn==true){ // char strg[10]; sprintf(strg,"%g",lev); PIFont myfont(GetGraphicAtt().GetFont()); g->SelFont(myfont); double px,py; px = (cntcur->coords)[0].x+2*(Xmax()-Xmin())/100.; py = (cntcur->coords)[0].y; g->DrawString(px,py,strg, PI_HorizontalCenter&&PI_VerticalCenter ); } numdr++; delete[] xx; xx=NULL; delete[] yy; yy=NULL; cntcur = cntcur->next; tot++; } //cout << "PIContourDrawer::Draw fin de trace des "<RestoreGraphicAtt(); } /* --Methode-- */ void PIContourDrawer::UpdateLimits() { // Doit calculer les limites double xmn = _xmin; double xmx = _xmax; double ymn = _ymin; double ymx = _ymax; SetLimits(xmn, xmx, ymn, ymx); } /* --Methode-- */ void PIContourDrawer::ShowControlWindow(PIBaseWdgGen* wdg) { PICnTools::SetCurrentBaseWdg(wdg); PICnTools::SetCurrentCnDrw(this); PICnTools::ShowPICnTools(); } //++ void PIContourDrawer::DeactivateControlWindow(PIBaseWdgGen* wdg) // // Desactivation de la fenetre de controle specialisee //-- { // si wdg != NULL, c'est un Detach (Drawer detache du PIBaseWdg // si wdg == NULL, c'est un delete du PIHisto2D (du PIDrawer) PICnTools::SetCurrentBaseWdg(NULL); PICnTools::SetCurrentCnDrw(NULL); PICnTools::HidePICnTools(); return; } /* --Methode-- */ bool PIContourDrawer::IsLabelOn() const { return (mLabelOn); } bool PIContourDrawer::IsLineOn() const { return (mLineOn); } bool PIContourDrawer::IsMarkOn() const { return (mMarkOn); } void PIContourDrawer::SetLabelOn(bool state) { mLabelOn = state; } void PIContourDrawer::SetMarkOn(bool state){ mMarkOn = state; } void PIContourDrawer::SetLineOn(bool state){ mLineOn = state; } void PIContourDrawer::InitAtts(){ mLineOn = true; mMarkOn = false; mLabelOn = false; } int PIContourDrawer::DecodeOptionString(vector & opt, bool rmdecopt) { // inspire de int PINTuple::DecodeOptionString // OP 11/2002 LAL Orsay int ndec = opt.size(); if (ndec < 1) return(0); // cout << " PIDrawer::DecodeOptionString : ndec = "<< ndec <" << ndec < udopt; // On gardera ici les options non decodees int iopt = 0; ndec = opt.size(); int recalc = 0; for (iopt = 0 ; ioptSetNLevel(nlev); this->SetCntLevelKind(LEVELS_NUM); recalc =1; }else if(deb=="niv"||deb=="lev"){ // hauteur des niveaux string fin = sopt.substr(sopt.find_first_of("=")+1 , string::npos); char * buff = strdup(fin.c_str()); char *tmp; tmp = strtok(buff,","); if(tmp==NULL){ cerr<< " PICOntourDrawer::DecodeOptionString ERREUR decodage des niveaux impossible " << buff < ztmp; while(tmp!=NULL){ ztmp.push_back(atof(tmp)); tmp = strtok(NULL,","); } int nlev = ztmp.size(); // cout << " PICOntourDrawer::DecodeOptionString "<SetCntLevelKind(LEVELS_DISCRETE); this->SetNLevel(nlev); this->SetMyLevels(ztmp); recalc =1; } }else if(deb=="lstep" ){ // niveaux incrementaux : args = nombre,depart,pas string fin = sopt.substr(sopt.find_first_of("=")+1 , string::npos); char * buff = strdup(fin.c_str()); char *tmp; tmp = strtok(buff,","); if(tmp==NULL){ cerr<< " PICOntourDrawer::DecodeOptionString ERREUR decodage nvx/incr. impossible " << buff < ztmp; while(tmp!=NULL){ ztmp.push_back(atof(tmp)); tmp = strtok(NULL,","); } if(ztmp.size() !=2) { cerr<< " PICOntourDrawer::DecodeOptionString ERREUR nb params incorrect(incr) " << buff < zlev; zlev.push_back(ztmp[1]); zlev.push_back(ztmp[2]); this->SetCntLevelKind(LEVELS_INCREMENTAL); this->SetNLevel( (int) ztmp[0]); this->SetMyLevels(zlev); recalc=1; } } }else{ cout<<"PICOntourDrawer::DecodeOptionString ERREUR option "<SetLabelOn(); else if ( sopt == "laboff" ) this->SetLabelOn(false); else if ( sopt == "autolevels" ) this->SetCntLevelKind(LEVELS_AUTO); else if (sopt == "linear" ){ this->SetCntKind(CONTOUR_KIND_BSPLINE); recalc = 1; } else if (sopt == "bspline" ){ this->SetCntKind(CONTOUR_KIND_BSPLINE); recalc = 1; }else if ((sopt == "3spl" )||(sopt == "cubicspl")) { this->SetCntKind(CONTOUR_KIND_CUBIC_SPL); recalc = 1; } else{ cout<<"PICOntourDrawer::DecodeOptionString ERREUR option "<CalcContour(); } this->Refresh(); // // S'il faut supprimer les options decodees, on remplace l'argument opt // par le vecteur des options non decodees. if (rmdecopt) opt = udopt; return(ndec+ndec1); } // Methode ajoute en Oct 2007 - Reza int PIContourDrawer::OptionToString(vector & opt) const { PIDrawer::OptionToString(opt); if (GetCntKind()==CONTOUR_KIND_LINEAR) opt.push_back("linear"); else if (GetCntKind()==CONTOUR_KIND_BSPLINE) opt.push_back("bspline"); else if (GetCntKind()==CONTOUR_KIND_CUBIC_SPL) opt.push_back("cubicspl"); if (IsLabelOn()) opt.push_back("labon"); else opt.push_back("laboff"); return 1; } void PIContourDrawer::GetOptionsHelpInfo(string& info) { info += " ---- PIContourDrawer options help info : \n" ; info += " autolevels : automatic selection of levels and number of contours \n"; info += " ncont=nLevel (or nc=NLevel) : sets the number of contour\n"; info += " lev=v1,v2,v3... (or niv=v1,v2,v3...) set the number and levels of contours\n"; info += " lstep=nLev,start,step : define incremental levels \n"; info += " labon/laboff : display of contour level values on/off \n"; info += " linear/bspline/cubicspl=3spl : select contour kind \n"; // On recupere ensuite la chaine info de la classe de base PIDrawer::GetOptionsHelpInfo(info); return; }