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

Last change on this file since 1838 was 1838, checked in by ansari, 24 years ago

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