source: Sophya/trunk/SophyaPI/PIGcont/pigncont.cc@ 1886

Last change on this file since 1886 was 1857, checked in by perderos, 24 years ago

1) correction de fuites de memoire
2) integration de glpc_set.* ds gp_contour (variables en commun)
File size: 15.8 KB
RevLine 
[1829]1#include "machdefs.h"
2#include <stdio.h>
3#include <stdlib.h>
4#include <iostream.h>
5#include <math.h>
6
7#include "histos.h"
8#include "ntuple.h"
9
[1838]10#include "nbtri.h"
[1829]11
12#include "picntools.h"
13#include "pigncont.h"
14
[1857]15#include <sys/time.h>
16#include <sys/resource.h>
17
[1838]18/* A virer ?? Reza 21/12/2001
19static TBOOLEAN multiplot;
20static double zero =0;
21static TBOOLEAN fg_polar;
22static int xleft, xright, ybot, ytop;
23static float xsize;
24static float ysize;
25static float xoffset, yoffset;
26*/
[1829]27
28GNUPlotContour::GNUPlotContour(struct iso_curve * iso){
29 int k=0;
30 iso1 = iso;
31 struct iso_curve *iso_cur = iso ;
32 struct iso_curve *iso_old;
33 _xmin = 1.e37;
34 _ymin = 1.e37;
35 _xmax = -1.e37;
36 _ymax = -1.e37;
37 _myiso = NULL;
38 _zmin = 1.e37;
39 _zmax = -1.e37;
40 while(iso_cur){
41 k++;
42 iso_cur = iso_cur->next;
43 _ny = iso_cur->p_count;
44 for(int k=0 ; k<_ny ; k++){
45 double xx=iso->points[k].x;
46 double yy=iso->points[k].y;
47 double zz=iso->points[k].z;
48 if(xx<_xmin) _xmin = xx;
49 if(xx>_xmax) _xmax = xx;
50 if(yy<_ymin) _ymin = yy;
51 if(yy>_ymax) _ymax = yy;
52 if(zz<_zmin) _zmin = zz;
53 if(zz>_zmax) _zmax = zz;
54
55 }
56 }
57
58 _nx = k;
59 _contours = NULL;
60 My_Levels = NULL;
61
62 _npolys=-1;
63 _Dxp =_Dyp =1.;
64 _Exy =_Invy =_Invy = FALSE ;
65
66}
67
68
69GNUPlotContour::GNUPlotContour(NTupleInterface *nt,int *nz=NULL,double *z=NULL,double *dz=NULL){
70 _xmin = 1.e37;
71 _ymin = 1.e37;
72 _xmax = -1.e37;
73 _ymax = -1.e37;
74 _zmin = 1.e37;
75 _zmax = -1.e37;
76 char *nom[4]={"x","y","z","niso"};
77 _myiso = new NTuple(4,nom);
78 //_myiso->Show();
79 float xnt[4];
80 double mn,mx;
81
82 nt->GetMinMax(0,_xmin,_xmax);
83 nt->GetMinMax(1,_ymin,_ymax);
84 nt->GetMinMax(2,_zmin,_zmax);
85
86
87 struct iso_curve *istmp = NULL ;
88 struct iso_curve *iso = NULL ;
89 struct iso_curve *iso_old = NULL ;
90 double x,y;
91 int nen = nt->NbLines();
92
93 // histo des y
94 for(int k= 0 ; k<1000 ; k++){
95 istmp = iso_alloc(100);
96 iso_free(istmp);
97 }
98
99#define NUMPTMIN 25
100#define MAXISO 100
101#define MINISO 5
102
103 Histo prep(_ymin,_ymax,nen);
104 for(int ix=0 ; ix<nen ; ix++) prep.Add(nt->GetCell(ix,1),1.);
105 Histo Intg = prep;
106 Intg.HInteg(0);
107 //prep.Print();
108 //Intg.Print();
109 float *bin;
110 int *cnt;
111
112 int niso = sqrt( (double)nen);
113 bin = new float[niso+1];
114 cnt = new int[niso+1];
115
[1842]116 int i;
117 for(i=0 ; i<niso ; i++){
[1829]118 bin[i] = prep.BinHighEdge(prep.BinPercent( ((double)i)/(double)niso));
119 cnt[i] = Intg(prep.BinPercent( ((double)i)/(double)niso));
120 }
121
122 bin[niso]=_ymax;
123
124 cnt[niso] = prep.NEntries();
125 int numx = 0;
126 int numi = 9999;
[1842]127 for(i=0;i<niso ; i++){
[1829]128 int numy = cnt[i+1]-cnt[i];
129
130 if(numy>0){
131 numx = numx>numy ? numx : numy ;
132 numi = numy>numi ? numi : numy ;
133 }
134 }
135
136
137 if (numi<=0) {
138
139 return;
140 }
141 _nx = 0;
142
[1842]143 for(i=0;i<niso ; i++){
[1829]144 double ymi = bin[i];
145 double ymx = bin[i+1];
146 int numy = cnt[i+1]-cnt[i];
147
148 if(numy<=0) continue;
149 double *xp;
150 xp = new double[numy];
151
152 int *indx;
153 indx = new int[numy];
154
155 istmp = iso_alloc(numy) ; // structure temporaire
156
157 int ip=0;
[1842]158 int ix;
159 for(ix=0 ; ix<nen ; ix++){
[1829]160
161 if(nt->GetCell(ix,1)>=ymi&&nt->GetCell(ix,1)<ymx){
162
163
164 istmp->points[ip].type = INRANGE;
165 istmp->points[ip].x = nt->GetCell(ix,0);
166 istmp->points[ip].y = nt->GetCell(ix,1);
167 istmp->points[ip].z = nt->GetCell(ix,2);
168 indx[ip] = ip;
169 xp[ip++] = nt->GetCell(ix,0);
170 }
171 }
172
173 // tri des points de l'iso_curve
174 //
175 istmp->p_count = ip;
176
177 tri_double(xp,indx,ip);
178
179 // remplissage structure finale
180
181 iso = iso_alloc(numx) ;
182 if( _nx == 0 ) iso1 = iso; //premiere structure remplie
183 if (iso_old) iso_old->next = iso;
[1842]184 for(ix=0 ; ix<numy ; ix++){
[1829]185 int jx = indx[ix];
186
187
188 iso->points[ix].type = INRANGE ;
189 iso->points[ix].x = istmp->points[jx].x ;
190 iso->points[ix].y = istmp->points[jx].y ;
191 iso->points[ix].z = istmp->points[jx].z ;
192 xnt[0] = iso->points[ix].x;
193 xnt[1] = iso->points[ix].y;
194 xnt[2] = iso->points[ix].z;
195 xnt[3] = _nx;
196 //cout << " fill avec "<<xnt[0]<<" "<<xnt[1] <<" "<<xnt[2] << " "<<xnt[3]<<endl;
197 // _myiso->Show();
198
199 _myiso->Fill(xnt);
200 }
201 //cout << " completion "<<ip<<" =>> "<<numx<<endl;
202 // on complete
[1842]203 for(ix=ip; ix<numx ; ix++ ){
[1829]204 int jmx = indx[ip-1];
205 iso->points[ix].type = INRANGE ;
206 iso->points[ix].x = istmp->points[jmx].x ;
207 iso->points[ix].y = istmp->points[jmx].y ;
208 iso->points[ix].z = istmp->points[jmx].z ;
209 xnt[0] = iso->points[ix].x;
210 xnt[1] = iso->points[ix].y;
211 xnt[2] = iso->points[ix].z;
212 xnt[3] = _nx;
213 //cout << " fill avec "<<xnt[0]<<" "<<xnt[1] <<" "<<xnt[2] << " "<<xnt[3]<<endl;
214 _myiso->Fill(xnt);
215 //_myiso->Show();
216 }
217
218 iso->p_count = numx;
219
220 iso_old = iso;
221 _nx++;
222 //cout << " destruction structure tempo "<<_nx <<endl;
223 if(istmp){
224 iso_free( istmp );// destruction de la structure temporaire
225 istmp = NULL;
226 }
227 if(xp){delete[] xp ; xp=NULL;}
228 if(indx) {delete[] indx ; indx = NULL;}
229
230 }
231 //if(_myiso!=NULL)_myiso->Write("myiso.ppf");
232 //cout << " fin du remplissage des iso_curves xmin,max "<< _xmin<<","<<_xmax<<" ymin,max "<<_ymin<<","<<_ymax<<endl;
233 _contours = NULL;
234 My_Levels = NULL;
235 _npolys=-1;
236 Contour_Levels_kind = LEVELS_AUTO ;
237 Contour_kind = CONTOUR_KIND_LINEAR ;
238 set_contour_kind(Contour_kind);
239 Contour_Levels = 5;
240 set_contour_levels(Contour_Levels);
241 set_contour_levels_kind(Contour_Levels_kind);
242
243
244 if(dz!=NULL){
245 My_Levels = new double[2];
246 Contour_Levels = *nz;
247 My_Levels[0] = *z ;
248 My_Levels[1] = *dz;
249 Contour_Levels_kind = LEVELS_INCREMENTAL ;
250 set_contour_levels(Contour_Levels);
251 set_contour_levels_kind(Contour_Levels_kind);
252 set_contour_levels_list(My_Levels);//,2);
253
254 }else if(nz!=NULL){
255
256 Contour_Levels = *nz;
257 Contour_Levels_kind = LEVELS_NUM ;
258 set_contour_levels(Contour_Levels);
259 set_contour_levels_kind(Contour_Levels_kind);
260 }
261
262}
263
264GNUPlotContour::GNUPlotContour(P2DArrayAdapter*arr,int *nz=NULL,double *z=NULL,double *dz=NULL){
265 _xmin = 1.e37;
266 _ymin = 1.e37;
267 _xmax = -1.e37;
268 _ymax = -1.e37;
269 _zmin = 1.e37;
270 _zmax = -1.e37;
271 // _nx = arr->XSize()-1;
272 // _ny = arr->YSize()-1;
273 _nx = arr->XSize() ;
274 _ny = arr->YSize() ;
275 //cout << " size x,y "<<_nx<<" , "<<_ny <<endl;
276 char *nom[4]={"x","y","z","niso"};
277
278 _myiso = new NTuple(4,nom);
279 float xnt[4];
280 struct iso_curve *iso =NULL ;
281 struct iso_curve *iso_old =NULL ;
282 double x,y;
283 // Calcul de la taille en X et en Y d'un pixel
284 double dxp, dyp;
285 double x1,y1;
286 arr->XYfromxy(0,0,x,y);
287 arr->XYfromxy(1,1,x1,y1);
288 _Dxp = (x1-x);
289 _Dyp = (y1-y);
290 double zz;
291 for(int ix=0 ; ix<_nx ; ix++){
292 iso = iso_alloc(_ny) ;
293 if(ix==0) iso1 = iso;
294 if (iso_old) iso_old->next = iso;
295 for(int iy = 0 ; iy<_ny ; iy++){
296 /*i = ix*100+iy;*/
297 arr->XYfromxy(ix,iy,x,y);
298
299 iso->points[iy].type = INRANGE;
300 iso->points[iy].x = x + _Dxp/2.;
301 iso->points[iy].y = y + _Dyp/2.;
302 zz=arr->Value(ix,iy);
303 iso->points[iy].z = zz;
304 //iso->points[iy].z = ((x-25)*(x-25) + (y-25)*(y-25))*.05 ;
305 if(x <_xmin) _xmin = x;
306 if(x >_xmax) _xmax = x;
307 if(y <_ymin) _ymin = y;
308 if(y >_ymax) _ymax = y;
309 if(zz <_zmin) _zmin = zz;
310 if(zz >_zmax) _zmax = zz;
311 //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);
312 xnt[0] = iso->points[iy].x;
313 xnt[1] = iso->points[iy].y;
314 xnt[2] = iso->points[iy].z;
315 xnt[3] = _nx;
316
317 _myiso->Fill(xnt);
318 }
319 iso->p_count = _ny;
320 iso_old = iso;
321 /*printf(" remplissage %lx old %lx old next %lx \n",iso,iso_old,iso->next,iso_old->next);*/
322 /*iso->next = iso ;*/
323 }
324 //cout << " fin du remplissage des iso_curves xmin,max "<< _xmin<<","<<_xmax<<" ymin,max "<<_ymin<<","<<_ymax<<endl;
325 _contours = NULL;
326 My_Levels = NULL;
327 _npolys=-1;
328 Contour_kind = CONTOUR_KIND_LINEAR ;
329 set_contour_kind(Contour_kind);
330 Contour_Levels_kind = LEVELS_AUTO ;
331 Contour_Levels = 5;
332 set_contour_levels(Contour_Levels);
333 set_contour_levels_kind(Contour_Levels_kind);
334
335
336 if(dz!=NULL){
337 My_Levels = new double[2];
338 Contour_Levels = *nz;
339 My_Levels[0] = *z ;
340 My_Levels[1] = *dz;
341 Contour_Levels_kind = LEVELS_INCREMENTAL ;
342 set_contour_levels(Contour_Levels);
343 set_contour_levels_kind(Contour_Levels_kind);
344 set_contour_levels_list(My_Levels);//,2);
345
346 }else if(nz!=NULL){
347 Contour_Levels = *nz;
348 Contour_Levels_kind = LEVELS_NUM ;
349 set_contour_levels(Contour_Levels);
350 set_contour_levels_kind(Contour_Levels_kind);
351 }
352
353}
354
355GNUPlotContour::~GNUPlotContour(){
356 cout << " destructeur de GNUPlotContour "<<endl;
357
358 struct iso_curve *iso_cur;
359 struct iso_curve *iso_nx;
360 iso_cur = iso1;
361
362 while(iso_cur){
363 // cout << " ~GNUPlotContour() : destruction de iso1 "<< iso_cur << endl;
364 iso_nx = iso_cur->next;
365 iso_free(iso_cur);iso_cur = NULL;
366
367 iso_cur = iso_nx;
368
369 }
370
371 iso1 = NULL;
372
373 struct gnuplot_contours *cntcur = _contours;
[1857]374 /*****
[1829]375 struct gnuplot_contours *cntold;
376
[1857]377 while(cntcur) {
[1829]378 //cout << " ~GNUPlotContour() : destruction de _contours" << _contours <<endl;
[1857]379 cntold = cntcur;
380 cntcur = cntold->next;
381 gp_free(cntold);
382 cntold=NULL;
383 }
384 *******/
385 contour_free(cntcur);
[1829]386 _contours = NULL;
387
388 if(My_Levels) {
389 //cout << " ~GNUPlotContour() : destruction de MyLevels "<<My_Levels <<endl;
390 delete[] My_Levels; My_Levels=NULL;}
391
392 if(_myiso!=NULL){ delete _myiso; _myiso=NULL;}
393
394}
395
396void GNUPlotContour::CalcContour(){
[1857]397 //cout << " GNUPlotContour::CalcContour(): determination des contours "<<endl;
[1829]398
[1857]399 struct rusage r_usage;
400 int rcus;
401
[1829]402 set_contour_kind(Contour_kind);
403 if(Contour_Levels_kind == LEVELS_INCREMENTAL){
404
405 set_contour_levels(Contour_Levels);
406 set_contour_levels_kind(Contour_Levels_kind);
407 //free_contour_levels_list();
408 set_contour_levels_list(My_Levels);//,2);
409
410 }else if(Contour_Levels_kind == LEVELS_DISCRETE ){
411 set_contour_levels(Contour_Levels);
412 set_contour_levels_kind(Contour_Levels_kind);
413 //free_contour_levels_list();
414 set_contour_levels_list(My_Levels);//,Contour_Levels);
415
416 }else if(Contour_Levels_kind == LEVELS_NUM ){
417
418 set_contour_levels(Contour_Levels);
419 set_contour_levels_kind(Contour_Levels_kind);
420 }else if(Contour_Levels_kind == LEVELS_AUTO){
421 Contour_Levels = 5;
422 set_contour_levels(Contour_Levels);
423 set_contour_levels_kind(Contour_Levels_kind);
424 }
425
[1857]426
427 rcus = getrusage( RUSAGE_SELF , &r_usage);
428
429 //if(rcus==0)
430 // cout << " rusage -> "<< r_usage.ru_maxrss <<" , "<< r_usage.ru_ixrss <<" , "<< r_usage.ru_ixrss <<endl;
431 //else
432 // perror("1er appel");
433
[1829]434 if(_contours) {
[1857]435 cout << " GNUPlotContour::CalcContour(): destruction des contours "<<_contours<<endl;
[1829]436 struct gnuplot_contours *cntcur = _contours;
[1857]437 /******
438 struct gnuplot_contours *cntold;
439 while(cntcur) {
440 cout << " GNUPlotContour::CalcContour(): destruction des contours "<< cntcur<<endl;
441 cntold = cntcur;
442 cntcur = cntold->next;
443 gp_free(cntold);
444 cntold = NULL;
445 }
446 ******/
447 contour_free(cntcur);
[1829]448 _contours = NULL;
449
450 }
451
[1857]452 //struct gnuplot_contours *cntcur;
[1829]453 //struct gnuplot_contours *cntold;
454
[1857]455
456 //rcus = getrusage( RUSAGE_SELF , &r_usage);
457 //if(rcus==0)
458 // cout << " rusage -> "<< r_usage.ru_maxrss <<" , "<< r_usage.ru_ixrss <<" , "<< r_usage.ru_ixrss <<endl;
459 //else
460 // perror("2d appel");
461
462
[1829]463 _contours = contour (_nx,iso1);
[1857]464
465
466 //rcus = getrusage( RUSAGE_SELF , &r_usage);
467 //if(rcus==0)
468 // cout << " rusage -> "<< r_usage.ru_maxrss <<" , "<< r_usage.ru_ixrss <<" , "<< r_usage.ru_ixrss <<endl;
469 //else
470 // perror("3ieme appel");
471
[1829]472 //free_contour_levels_list();
473
[1857]474
[1829]475}
476//_+_ Methode
477void
478GNUPlotContour::SetMyLevels(double *ptr,int k){
479
480 if(My_Levels!=NULL) delete My_Levels;
481 My_Levels = new double[k];
482 for(int i=0 ; i<k ; i++)My_Levels[i] = ptr[i];
483}
484
485// Class PIContourDrawer
486/* --Methode-- */
487PIContourDrawer::PIContourDrawer(P2DArrayAdapter* arr,bool autodel,int *nz,double *z0,double *dz )
488 :GNUPlotContour(arr,nz,z0,dz)
489
490{
491 _arr = arr;
492 _nti = NULL;
493 _nz=5;
494 if(nz!=NULL)_nz = *nz;
495 _autodel = autodel;
[1853]496
497 // This drawer has specific control tools
498 mFgSpecContWind = true;
[1829]499 InitAtts();
500}
501
502PIContourDrawer::PIContourDrawer(NTupleInterface* nti,bool autodel,int *nz,double *z0,double *dz )
503 :GNUPlotContour(nti,nz,z0,dz)
504
505{
506 _arr = NULL;
507 _nti = nti;
508 _nz=5;
509 if(nz!=NULL)_nz = *nz;
510 _autodel = autodel;
511 InitAtts();
512}
513
514
515
516/* --Methode-- */
517PIContourDrawer::~PIContourDrawer()
518{
519 if(_autodel){
520 if(_arr!=NULL){delete _arr; _arr=NULL;}
521 if(_nti!=NULL){delete _nti ; _nti=NULL;}
522 }
523 //if (mColorMap!=NULL){delete mColorMap; mColorMap=NULL;}
524}
525
526/* --Methode-- */
527void PIContourDrawer::Draw(PIGraphicUC* g, double xmin, double ymin, double xmax, double ymax)
528{
529 // double xmin, double ymin, double xmax, double ymax LIMITES DE LA ZONE A REFRESHER
530 //PIGrCoord * xx = new PIGrCoord[10];
531 //PIGrCoord * yy = new PIGrCoord[10];
532 PIGrCoord * xx =NULL;
533 PIGrCoord * yy =NULL;
534 PIGrCoord xa;
535 PIGrCoord ya;
536 char buff[64];
537 // On met les bons attributs graphiques
538 //SelGraAtt(g);
539 g->SaveGraphicAtt();
540 struct gnuplot_contours *cntcur;
541 struct gnuplot_contours *cntold;
542 cntcur = _contours;
543
544 double a,b;
545 int tot=0;
546 bool Invx, Invy, Exy;
547 PIColorMap * cmap = NULL;
548 if(mCmapid !=CMAP_OTHER ){
549 cmap = new PIColorMap(mCmapid);
550 }
551 int numdr = 0;
552 while(cntcur){
553 int npts = cntcur->num_pts ;
554 if(npts<0) continue;
555
556 xx = new PIGrCoord[npts];
557 yy = new PIGrCoord[npts];
558 double lev = (cntcur->coords)[0].z ;
559 for(int l=0 ; l< npts ; l++){
560 a = (cntcur->coords)[l].x;
561 b = (cntcur->coords)[l].y;
562
563 xx[l] = a ;
564 yy[l] = b ;
565
566 }
567
568 // choix de la couleur
569
570 g->SelForeground(mFCol);
571 if(mCmapid !=CMAP_OTHER ){
572 if(_zmax-_zmin!=0){
573 int kc = ((lev - _zmin)/(_zmax-_zmin))*cmap->NCol();
574
575 g->SelForeground(*cmap,kc);
576 }
577 }
578 // traitement des des markers
579 if(IsMarkOn()==true){
580
581 g->SelMarker(mMSz,mMrk);
582 g->DrawMarkers(xx,yy,npts);
583 }
584 // traitement des lignes
585 if(mLineOn==true){
586
587 g->SelLine(mLAtt);
588 g->DrawPolygon(xx,yy,npts,false);
589 }
590 // traitement des labels
591 // SIMPLISTE POUR L'INSTANT : UN AFFICHAGE
592 if(mLabelOn==true){
593 //
594
595 char strg[10];
596 sprintf(strg,"%g",lev);
597 PIFont myfont(mFName);
598 myfont.SetFontAtt(mFAtt);
599 myfont.SetFontSz(mFSz);
600
601 g->SelFont(myfont);
602 double px,py;
603 px = (cntcur->coords)[0].x+2*(Xmax()-Xmin())/100.;
604 py = (cntcur->coords)[0].y;
605 g->DrawString(px,py,strg, PI_HorizontalCenter&&PI_VerticalCenter );
606
607 }
608 numdr++;
609 delete[] xx;
610 xx=NULL;
611 delete[] yy;
612 yy=NULL;
613 cntcur = cntcur->next;
614 tot++;
615 }
616
617 //cout << "PIContourDrawer::Draw fin de trace des "<<tot<< " polygones "<<endl;
618
619 //
620 if(cmap !=NULL)delete cmap;
621
622 g->RestoreGraphicAtt();
623}
624
625/* --Methode-- */
626void PIContourDrawer::UpdateLimits()
627{
628 // Doit calculer les limites
629
630 double xmn = _xmin;
631 double xmx = _xmax;
632 double ymn = _ymin;
633 double ymx = _ymax;
634
635 SetLimits(xmn, xmx, ymn, ymx);
636}
637
638
639/* --Methode-- */
640void PIContourDrawer::ShowControlWindow(PIBaseWdgGen* wdg)
641{
642 PICnTools::SetCurrentBaseWdg(wdg);
643 PICnTools::SetCurrentCnDrw(this);
644 PICnTools::ShowPICnTools();
645}
646
647
648/* --Methode-- */
649bool PIContourDrawer::IsLabelOn(){
650 return (mLabelOn);
651}
652
653bool PIContourDrawer::IsLineOn(){
654 return (mLineOn);
655}
656
657bool PIContourDrawer::IsMarkOn(){
658 return (mMarkOn);
659}
660
661void PIContourDrawer::SetLabelOn(bool state){
662 mLabelOn = state;
663}
664
665void PIContourDrawer::SetMarkOn(bool state){
666 mMarkOn = state;
667}
668
669void PIContourDrawer::SetLineOn(bool state){
670 mLineOn = state;
671}
672
673void PIContourDrawer::InitAtts(){
674 mLineOn = false;
675 mMarkOn = true;
676 mLabelOn = false;
677}
678
679
680
Note: See TracBrowser for help on using the repository browser.