← Back to branch summary

~ubuntu-branches/debian/jessie/eso-midas/jessie

« back to all changes in this revision

Viewing changes to prim/plot/src/plotcon.c

  • Committer: Package Import Robot
  • Author(s): Ole Streicher
  • Date: 2014-04-22 14:44:58 UTC
  • Revision ID: package-import@ubuntu.com-20140422144458-okiwi1assxkkiz39
Tags: upstream-13.09pl1.2+dfsg
Import upstream version 13.09pl1.2+dfsg

Show diffs side-by-side

added added

removed removed

Lines of Context:
prim/plot/src/plotcon.c
 
1
/*===========================================================================
 
2
  Copyright (C) 1995-2006 European Southern Observatory (ESO)
 
3
 
 
4
  This program is free software; you can redistribute it and/or 
 
5
  modify it under the terms of the GNU General Public License as 
 
6
  published by the Free Software Foundation; either version 2 of 
 
7
  the License, or (at your option) any later version.
 
8
 
 
9
  This program is distributed in the hope that it will be useful,
 
10
  but WITHOUT ANY WARRANTY; without even the implied warranty of
 
11
  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 
12
  GNU General Public License for more details.
 
13
 
 
14
  You should have received a copy of the GNU General Public 
 
15
  License along with this program; if not, write to the Free 
 
16
  Software Foundation, Inc., 675 Massachusetts Ave, Cambridge, 
 
17
  MA 02139, USA.
 
18
 
 
19
  Correspondence concerning ESO-MIDAS should be addressed as follows:
 
20
        Internet e-mail: midas@eso.org
 
21
        Postal address: European Southern Observatory
 
22
                        Data Management Division 
 
23
                        Karl-Schwarzschild-Strasse 2
 
24
                        D 85748 Garching bei Muenchen 
 
25
                        GERMANY
 
26
===========================================================================*/
 
27
 
 
28
/*++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
 
29
.IDENTifer   PLOTCON
 
30
.AUTHOR      R.M. van Hees IPG-ESO Garching
 
31
.KEYWORDS    Graphics, bulk data frame, one-dimensional plotting
 
32
.LANGUAGE    C
 
33
.PURPOSE     Plot or overplot a contour map of two-dimensional FRAME
 
34
  in/output: IN_A/C/1/60   = input frame
 
35
             IN_B/C/1/60   = coord_string
 
36
             P3/R/1/4      = scales in X and Y and offset X and Y
 
37
                             (default is auto scaling to device filling)
 
38
             INPUTC/C/1/72 = contour levels
 
39
             INPUTC/C/73/4 = contour type
 
40
             INPUTI/I/1/1  = smoothing parameter 
 
41
 
 
42
.COMMENTS    none
 
43
.ENVIRONment MIDAS
 
44
             #include <midas_def.h>     Prototypes for MIDAS interfaces
 
45
             #include <plot_def.h>      General symbols for Plot routines
 
46
 
 
47
.VERSION     1.2     25-Nov-1993   Removed restrictions on frame size, RvH
 
48
 
 
49
 060524         last modif
 
50
 
 
51
------------------------------------------------------------*/
 
52
  
 
53
#define _POSIX_SOURCE 1
 
54
 
 
55
/*
 
56
 * definition of the used functions
 
57
 */
 
58
#include <stdio.h>
 
59
#include <string.h>
 
60
#include <stdlib.h>
 
61
#include <math.h>
 
62
 
 
63
/*
 
64
 * define some macros and constants
 
65
 */
 
66
 
 
67
#include <midas_def.h>
 
68
#include <plot_def.h>
 
69
 
 
70
#define MAXLEV          50                  /* maximum num of contour levels */
 
71
#define MAXPIX          512    /* max frame dimension (X,Y) accessed at once */
 
72
#define MAXSIZ          (MAXPIX * MAXPIX)
 
73
 
 
74
/* 
 
75
 
 
76
*/
 
77
 
 
78
/*++++++++++++++++++++
 
79
.IDENTifer   GET_CONTOURS
 
80
.PURPOSE     fill two arrays with contour levels and line types
 
81
     output:   int    *nlevl  number of requested contour levels
 
82
               float  *clevl  array with contour levels
 
83
               int    *ctype  line type for each contour
 
84
.COMMENTS    static function
 
85
--------------------*/
 
86
#ifdef __STDC__
 
87
      static void GET_CONTOURS( int *nlevl, float *clevl, int *ctype )
 
88
#else
 
89
      static void GET_CONTOURS( nlevl, clevl, ctype )
 
90
      int   *nlevl, *ctype;
 
91
      float *clevl;
 
92
#endif
 
93
{
 
94
int  actvals, ic, ltype;
 
95
char *cbuff, option[5], input[72];
 
96
 
 
97
char  *err_usrin = "*** FATAL: error detected in USRINP";
 
98
 
 
99
/*
 
100
 * get contour levels as character string
 
101
 */
 
102
(void) SCKGETC( "INPUTC", 1, 72, &actvals, input );
 
103
 
 
104
cbuff  = (char *) clevl;
 
105
if ( USRINP( 'r', input, MAXLEV, cbuff, nlevl ) != ERR_NORMAL )
 
106
   SCETER( 1, err_usrin );
 
107
/*
 
108
 * sort the contour levels
 
109
 */
 
110
SORLEV( *nlevl-1, clevl );
 
111
 
 
112
/*
 
113
 * get C_TYPE option for the plotting of the contours
 
114
 */
 
115
(void) SCKGETC( "INPUTC", 73, 4, &actvals, option );
 
116
if ( *option == 'n' || *option == 'N' )          /* dashed negative contours */
 
117
   { for ( ic = 0; ic < *nlevl; ic++ )
 
118
         { if ( fabs( (double) clevl[ic] ) < PLT_EPS )
 
119
              ctype[ic] = 2;                 /*dashed*/
 
120
           else if ( clevl[ic] < 0.0 )
 
121
              ctype[ic] = 1;                 /*dotted*/
 
122
           else
 
123
              ctype[ic] = 0;                 /*solid*/
 
124
         }
 
125
   }
 
126
else if ( *option == 'l' || *option == 'L' )                    /* use LTYPE */
 
127
   { PCKRDI( "LTYPE", 1, &actvals, &ltype );
 
128
     ltype = MYMAX( ltype-1, 0 );
 
129
     for ( ic = 0; ic < *nlevl; ic++ ) ctype[ic] = ltype;
 
130
   }
 
131
else                                             /* draw odd contours dashed */
 
132
   { for ( ic = 0; ic < *nlevl; ic++ )
 
133
         { if ( ic % 2 == 0 )
 
134
              ctype[ic] = 1;                 /*dotted*/
 
135
           else
 
136
              ctype[ic] = 0;                 /*solid*/
 
137
         }
 
138
   }
 
139
return;
 
140
}
 
141
/* 
 
142
 
 
143
*/
 
144
 
 
145
/*++++++++++++++++++++++++++++++++++++++++++++++++++
 
146
 *
 
147
 * here starts the code of the function
 
148
 */
 
149
 
 
150
int main()
 
151
 
 
152
{
 
153
register int  ii;
 
154
int     actvals, chunks, imf, indx, knul, naxis, nr, size, stat, unit, 
 
155
        ctype[MAXLEV], ndum[PLDIM2], npix[PLDIM2], nrpix[PLDIM2], 
 
156
        sublo[PLDIM2], subhi[PLDIM2];
 
157
 
 
158
float   amin, amax, *p_img, last_row, zmin, zmax, area[4], image[4], 
 
159
        wcfram[12], clevl[MAXLEV];
 
160
 
 
161
double  start[PLDIM2], step[PLDIM2];
 
162
 
 
163
char    cmnd[21], ident[33], cunit[49], name[61], input[61], buff[81],
 
164
        *label[4];
 
165
/*
 
166
 * initialised variables
 
167
 */
 
168
char  *err_1dim  = "*** FATAL: Frame has only one dimension",
 
169
      *err_coord = "*** FATAL: invalid coordinate input ...",
 
170
      *err_xcuts = "*** FATAL: range in x has no overlap with current graph abscissa - NO PLOT",
 
171
      *err_ycuts = "*** FATAL: range in y has no overlap with current graph abscissa - NO PLOT",
 
172
      *info_flat = "*** WARNING: zero dynamic range in data at %13.8g",
 
173
      *info_levl = "*** WARNING: no contour level falls within the dynamic range of your data";
 
174
 
 
175
static char *axis[PLDIM2] = { "MANU", "MANU" };
 
176
 
 
177
int  access =  0,                      /* parameter for PCOPEN: plot mode */
 
178
     plmode = -1,                   /* plot mode taken from keyword PMODE */
 
179
     nlevl  =  0,                   /* Default number of contours is zero */
 
180
     ismoot =  1,                      /* Default no smooting of the data */
 
181
     more_in = 0,                   /* for cursor input, we may want more */
 
182
     nrdrawn = 0;                                /* no contours drawn yet */
 
183
 
 
184
/*
 
185
 * allocate memory for different character pointers and initialise a few
 
186
 */
 
187
for ( ii = 0; ii < 4; ii++ ) label[ii] = osmmget(81);
 
188
 
 
189
(void) strcpy( label[0], "Position (" );
 
190
(void) strcpy( label[1], "Position (" );
 
191
(void) strcpy( label[2], "Frame: " );
 
192
(void) strcpy( label[3], "Ident: " );
 
193
 
 
194
/*
 
195
 * start of executable code
 
196
 */
 
197
(void) SCSPRO( "PLTCON" );                   /*contact with the MIDAS monitor*/
 
198
 
 
199
/*
 
200
 * plot or overplot mode
 
201
 */
 
202
(void) SCKGETC( "MID$CMND", 1, 20, &actvals, cmnd );
 
203
if ( *cmnd == 'O' ) access = 1;                         
 
204
 
 
205
/*
 
206
 * find file name and read header information
 
207
 */
 
208
(void) SCKGETC( "IN_A", 1, 60, &actvals, name );
 
209
(void) SCFOPN( name, D_R4_FORMAT, 0, F_IMA_TYPE, &imf );
 
210
(void) SCDRDI( imf, "NAXIS", 1, 1, &actvals, &naxis, &unit, &knul );
 
211
(void) SCDRDI( imf, "NPIX" , 1, PLDIM2, &actvals, npix , &unit, &knul );
 
212
/*
 
213
 * check frame parameters
 
214
 */
 
215
if ( naxis < 2 || (npix[0] == 1 || npix[1] == 1) ) SCETER( 1, err_1dim );
 
216
 
 
217
/*
 
218
 * read the descriptor
 
219
 */
 
220
(void) SCDRDD( imf, "START", 1, PLDIM2, &actvals, start, &unit, &knul );
 
221
(void) SCDRDD( imf, "STEP" , 1, PLDIM2, &actvals, step , &unit, &knul );
 
222
(void) SCDGETC( imf, "IDENT", 1, 32, &actvals, ident );
 
223
(void) SCDGETC( imf, "CUNIT", 1, 48, &actvals, cunit );
 
224
 
 
225
/*
 
226
 * Get the manual setting for the axes
 
227
 */
 
228
PCKRDR( "XAXIS", 4, &actvals, wcfram );
 
229
PCKRDR( "YAXIS", 4, &actvals, wcfram+FOR_Y );
 
230
 
 
231
/*
 
232
 * read window coordinates and take action
 
233
 */
 
234
(void) SCKGETC( "IN_B", 1, 60, &actvals, input );
 
235
 
 
236
if ( *input == 'm' || *input == 'M' )                      /* manual scaling */
 
237
   {
 
238
   BOXWTP( wcfram,       npix[0], start[0], step[0], image );
 
239
   BOXWTP( wcfram+FOR_Y, npix[1], start[1], step[1], image+2 );
 
240
   }
 
241
else
 
242
   {
 
243
   if ( *input == 'c' || *input == 'C' )                  /* display input */
 
244
      {
 
245
      (void) SCKRDR( "OUTPUTR", 10, 1, &actvals, image, &unit, &knul);
 
246
      (void) SCKRDR( "OUTPUTR", 11, 1, &actvals, image+2, &unit, &knul); 
 
247
      (void) SCKRDR( "OUTPUTR", 15, 1, &actvals, image+1, &unit, &knul);
 
248
      (void) SCKRDR( "OUTPUTR", 16, 1, &actvals, image+3, &unit, &knul);
 
249
      more_in = 1;
 
250
      }
 
251
   else                                               /* automatic scaling */
 
252
      {
 
253
      stat = Convcoo(1,imf,input,PLDIM2,ndum,sublo,subhi);
 
254
      if ( stat != 0 ) SCETER( 2, err_coord );
 
255
 
 
256
      image[0] = sublo[0] + 1;
 
257
      image[1] = subhi[0] + 1;
 
258
      image[2] = sublo[1] + 1;
 
259
      image[3] = subhi[1] + 1;
 
260
      }
 
261
   }
 
262
 
 
263
BOXPTW( image,   npix[0], start[0], step[0], area );
 
264
BOXPTW( image+2, npix[1], start[1], step[1], area+2 );
 
265
 
 
266
PCKWRR( "PIXEL", 4, image );
 
267
 
 
268
if ( access == 0 )
 
269
   { 
 
270
/*
 
271
 * get size of frame along X-axis
 
272
 */
 
273
     if ( fabs( *wcfram ) < PLT_EPS && fabs( *(wcfram+1) ) < PLT_EPS )
 
274
        { axis[0] = "AUTO";
 
275
          wcfram[0] = area[0];
 
276
          wcfram[1] = area[1]; 
 
277
          wcfram[2] = wcfram[3] = 0.0;
 
278
        }
 
279
/*
 
280
 * get size of frame along Y-axis
 
281
 */
 
282
     if ( fabs( *(wcfram+FOR_Y) ) < PLT_EPS 
 
283
          && fabs( *(wcfram+FOR_Y+1) ) < PLT_EPS )
 
284
        { axis[1] = "AUTO";
 
285
          wcfram[FOR_Y]   = area[2];
 
286
          wcfram[FOR_Y+1] = area[3]; 
 
287
          wcfram[FOR_Y+2] = wcfram[FOR_Y+3] = 0.0;
 
288
        }
 
289
     GETFRM( axis[0], wcfram );
 
290
     GETFRM( axis[1], wcfram + FOR_Y );
 
291
     PCKWRR( "XWNDL", 4, wcfram );
 
292
     PCKWRR( "YWNDL", 4, wcfram+FOR_Y );
 
293
   }
 
294
else                                                         /* overplot mode*/
 
295
   { PCKRDR( "XWNDL", 4, &actvals, wcfram );
 
296
     PCKRDR( "YWNDL", 4, &actvals, wcfram+FOR_Y );
 
297
/*
 
298
 * does overplot data  fall within plotted frame? 
 
299
 */
 
300
     amin = MYMIN( *wcfram, *(wcfram + 1) );
 
301
     amax = MYMAX( *wcfram, *(wcfram + 1) );
 
302
     if ( ( MYMAX( area[0], area[1] ) < amin ) ||
 
303
          ( MYMIN( area[0], area[1] ) > amax ) )
 
304
        SCETER( 3, err_xcuts );
 
305
 
 
306
     amin = MYMIN( *(wcfram + FOR_Y), *(wcfram + FOR_Y + 1) );
 
307
     amax = MYMAX( *(wcfram + FOR_Y), *(wcfram + FOR_Y + 1) );
 
308
     if ( ( MYMAX( area[2], area[3] ) < amin ) ||
 
309
          ( MYMIN( area[2], area[3] ) > amax ) )
 
310
        SCETER( 4, err_ycuts );
 
311
   }
 
312
/*
 
313
 * get contours levels
 
314
 */
 
315
GET_CONTOURS( &nlevl, clevl, ctype );
 
316
 
 
317
/*
 
318
 * setup graphic device according to MIDAS settings
 
319
 */
 
320
PCOPEN( " ", " ", access, &plmode );
 
321
 
 
322
/*
 
323
 * get the smooting parameter
 
324
 */
 
325
(void) SCKRDI( "INPUTI", 1, 1, &actvals, &ismoot, &unit, &knul);
 
326
 
 
327
/*
 
328
 * determine number of pixels and number of chunks
 
329
 */
 
330
nrpix[0] = (int) fabs( image[1] - image[0] ) + 1;
 
331
nrpix[1] = (int) fabs( image[3] - image[2] ) + 1;
 
332
last_row = MYMAX( image[2], image[3]);
 
333
chunks   = (int) ceil( (double) nrpix[0] * nrpix[1] / MAXSIZ );
 
334
nrpix[1] = (int) ceil( (double) nrpix[1] / chunks );
 
335
 
 
336
/*
 
337
 * allocate virtual memory and scratch space
 
338
 */
 
339
size  = nrpix[0] * nrpix[1];
 
340
p_img = (float *) osmmget( size * sizeof( float ));
 
341
 
 
342
wcfram[FOR_Z]   =  1e+12,                                /* minimum in frame */
 
343
wcfram[FOR_Z+1] = -1e+12;                                /* maximum in frame */
 
344
for ( ii = 0; ii < chunks; ii++ )
 
345
    { if ( image[3] > image[2] )                        /* size of chunk p.c.*/
 
346
         { if ( ii > 0 ) image[2] += nrpix[1] - 1.0;
 
347
           image[3] = MYMIN( last_row, image[2] + nrpix[1] - 1 );
 
348
         }
 
349
      else
 
350
         { if ( ii > 0 ) image[3] += nrpix[1] - 1.0;
 
351
           image[2] = MYMIN( last_row, image[3] + nrpix[1] - 1 );
 
352
         }
 
353
/*
 
354
 * get size of chunk in w.c.
 
355
 */
 
356
      BOXPTW( image+2, npix[1], start[1], step[1], area+2 );      
 
357
/*
 
358
 * extract chunck from original frame
 
359
 */
 
360
      GETDAT( imf, MAXSIZ, npix, image, ismoot, p_img );
 
361
/*
 
362
 * get dynamic range of the data
 
363
 */
 
364
      MINMAX( p_img, size, &zmin, &zmax );
 
365
      wcfram[FOR_Z]   = MYMIN( wcfram[FOR_Z], zmin );
 
366
      wcfram[FOR_Z+1] = MYMAX( wcfram[FOR_Z+1], zmax );
 
367
/*
 
368
 * determine the actual number of contours to be plotted
 
369
 */
 
370
      indx = 0;
 
371
      while ( clevl[indx] < zmin ) indx++;
 
372
      nr = nlevl;
 
373
      while ( clevl[nr-1] > zmax ) nr--;
 
374
      nr -= indx;
 
375
      nrdrawn = MYMAX( nr, nrdrawn );
 
376
/*
 
377
 * draw the contours
 
378
 */
 
379
      PLCON( p_img, image, area, step, nr, clevl+indx, ctype+indx );
 
380
/*
 
381
 * updating for the next loop
 
382
 */
 
383
      nrpix[1] = MYMIN( last_row - ii * nrpix[1], nrpix[1] );
 
384
      size  = nrpix[0] * nrpix[1];
 
385
    }
 
386
(void) SCFCLO( imf );
 
387
 
 
388
/*
 
389
 * give a warning if no contour is drawn
 
390
 */
 
391
if ( nrdrawn == 0 ) SCTPUT( info_levl );
 
392
 
 
393
/*
 
394
 * store min and max found in the whole frame
 
395
 */
 
396
if ( wcfram[FOR_Z] == wcfram[FOR_Z+1] )
 
397
   { (void) sprintf( buff, info_flat, wcfram+FOR_Z );
 
398
     SCTPUT( buff );
 
399
   }
 
400
PCKWRR( "ZWNDL", 2, wcfram+FOR_Z );
 
401
 
 
402
/*
 
403
 * draw the axes and the label
 
404
 */
 
405
if ( plmode >= 0 && access == 0 )
 
406
   { if ( strlen( cunit ) > (size_t) 32 )
 
407
        { (void) strcat( label[1], cunit+32 );
 
408
          *(cunit+32) = '\0';
 
409
        }
 
410
     if ( strlen( cunit ) > (size_t) 16 ) (void) strcat( label[0], cunit+16 );
 
411
 
 
412
     for ( ii = 0; ii < PLDIM2; ii++ ) 
 
413
         { (void) strcat( label[ii], ")" );
 
414
           LABSTR( label[ii] );
 
415
         }
 
416
/*
 
417
 * get format for the axes
 
418
 */
 
419
/*
 
420
 * plot axes and labels
 
421
 */
 
422
     PCFRAM( wcfram, wcfram+FOR_Y, label[0], label[1] );
 
423
 
 
424
     if ( plmode == 1 )
 
425
        { (void) strcat( label[2], name );
 
426
          (void) strcat( label[3], ident );
 
427
          PLIDEN( plmode, label[2], label[3] );
 
428
        }
 
429
     else if ( plmode == 2 )
 
430
        PLCONI( plmode, name, ident, clevl, ctype, nlevl );
 
431
   }
 
432
/*
 
433
 * close plot file and terminate graphic operations
 
434
 */
 
435
PCCLOS();
 
436
 
 
437
if (more_in == 1)
 
438
   {                    /* thus, we'll loop on cursor input */
 
439
   (void) SCKWRI("MORE_IN",&more_in,1,1,&unit);
 
440
   }
 
441
 
 
442
return SCSEPI();
 
443
}