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

Last change on this file since 1837 was 1829, checked in by perderos, 24 years ago

Creation du module PIGcont, Drawer de trace de contour / a partir
de l'adaptation du code de trace de contour de GNUplot, interfacage
avec la librairie PI .

Olivier Perdereau 19/12/2001

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