00001 #define MAIN_COXPLOT_FILE
00002 #include "coxplot.h"
00003 
00004 /*****************************************************************************
00005   This software is copyrighted and owned by the Medical College of Wisconsin.
00006   See the file README.Copyright for details.
00007 ******************************************************************************/
00008 
00009 static int             num_plotar  = 0 ;
00010 static MEM_plotdata ** plotar      = NULL ;
00011 static int             active_plot = -1 ;
00012 
00013 static float           active_color = (float) RGB_TO_COL(1.0,1.0,1.0) ;
00014 static float           active_thick = 0.0 ;
00015 
00016 static float           active_opacity = 1.0 ;   /* 22 Jul 2004 */
00017 
00018 #define STATUS(str) fprintf(stderr,"** " str "\n")
00019 
00020 /*------------------------------------------------------------------------
00021    Function to return a pointer to an in-memory plot with the given
00022    ident.  A NULL return value indicates failure.
00023    A NULL input string will return the current "active" plot.
00024 --------------------------------------------------------------------------*/
00025 
00026 MEM_plotdata * find_memplot( char * id )
00027 {
00028    int ip ;
00029 
00030    if( num_plotar == 0 || plotar == NULL ) return NULL ;
00031 
00032    if( id == NULL || id[0] == '\0' ){
00033       if( active_plot < 0 || active_plot >= num_plotar ) return NULL ;
00034       return plotar[active_plot] ;
00035    }
00036 
00037    for( ip=0 ; ip < num_plotar ; ip++ )
00038       if( strcmp(plotar[ip]->ident,id) == 0 ) return plotar[ip] ;
00039 
00040    return NULL ;
00041 }
00042 
00043 /*------------------------------------------------------------------------
00044   Function to create an in-memory plot with the given ident and aspect
00045   ratio (length of x / length of y).  Also sets the "active" plot to be
00046   this new one.  Nonzero return value indicates an error.
00047 --------------------------------------------------------------------------*/
00048 
00049 int create_memplot( char * id , float aspect )
00050 {
00051    MEM_plotdata * pd ;
00052    static int plotpak_framed = 0 ;
00053    real asp ;
00054 
00055    if( find_memplot(id) != NULL ) return 1 ;
00056 
00057    INIT_MEMPLOT(pd,id) ;
00058 
00059    if( plotar == NULL ){
00060       plotar     = (MEM_plotdata **) malloc( sizeof(MEM_plotdata *) ) ;
00061       num_plotar = 0 ;
00062    } else {
00063       plotar = (MEM_plotdata **)
00064                   realloc( plotar , sizeof(MEM_plotdata *)*(num_plotar+1) ) ;
00065    }
00066 
00067    active_plot = num_plotar ;
00068    plotar[num_plotar++] = pd ;
00069 
00070    ADDTO_MEMPLOT( pd , 1.0,0.0,0.0,0.0 , 0.0 , -THCODE_OPAC ) ;  /* 22 Jul 2004 */
00071 
00072    if( aspect <= 0.0 ) aspect = 1.3 ;
00073    asp        = aspect ;
00074    pd->aspect = aspect ;
00075    memplt_( &asp ) ;                  /* setup PLOTPAK */
00076 
00077    return 0 ;
00078 }
00079 
00080 /*------------------------------------------------------------------
00081    20 Sep 2001: make a plot, fer shur
00082 --------------------------------------------------------------------*/
00083 
00084 int create_memplot_surely( char *id , float aspect )
00085 {
00086    int ii , jj ;
00087    char str[256] ;
00088 
00089    if( aspect <= 0.0 ) aspect = 1.0 ;  /* backup for stupid users */
00090 
00091    if( id != NULL && id[0] != '\0'){
00092       ii = create_memplot(id,aspect) ;
00093       if( ii == 0 ) return 0 ;
00094    } else {
00095       id = "ElvisWalksTheEarth" ;
00096    }
00097 
00098    for( jj=0 ; ; jj++ ){
00099       sprintf(str,"%.240s_%d",id,jj) ;
00100       ii = create_memplot(str,aspect) ;
00101       if( ii == 0 ) return 0 ;
00102    }
00103 
00104    return 1 ; /* actually, unreachable */
00105 }
00106 
00107 /*-------------------------------------------------------------------------
00108    Function to set the "active" plot to a given ident.
00109    A nonzero return value indicates an error.
00110 ---------------------------------------------------------------------------*/
00111 
00112 int set_active_memplot( char * id )
00113 {
00114    int ip ;
00115 
00116    if( id == NULL || id[0] == '\0' || num_plotar == 0 || plotar == NULL )
00117       return 1 ;
00118 
00119    for( ip=0 ; ip < num_plotar ; ip++ )
00120       if( strcmp(plotar[ip]->ident,id) == 0 ){
00121          real asp = plotar[ip]->aspect ;
00122          active_plot = ip ;
00123          memplt_( &asp ) ;    /* re-setup PLOTPAK */
00124          return 0 ;
00125       }
00126 
00127    return 1 ;
00128 }
00129 
00130 MEM_plotdata * get_active_memplot(void)
00131 {
00132    return find_memplot(NULL) ;
00133 }
00134 
00135 int nline_active_memplot(void)
00136 {
00137    MEM_plotdata * mp ;
00138    mp = find_memplot(NULL) ;
00139    if( mp == NULL ) return 0 ;
00140    return MEMPLOT_NLINE(mp) ;
00141 }
00142 
00143 /*-------------------------------------------------------------------------
00144   Functions to set line thickness and color.
00145   Color is given as an RGB triple from [0,1] x [0,1] x [0,1].
00146   Thickness is given in the same units as coordinates; zero means thin.
00147 ---------------------------------------------------------------------------*/
00148 
00149 void set_color_memplot( float r , float g , float b )
00150 {
00151    if( r > 1.0 || g > 1.0 || b > 1.0 ){        /* 22 Mar 2002:     */
00152       r /= 255.0 ; g /= 255.0 ; b /= 255.0 ;   /* allow for 0..255 */
00153    }
00154    if( r < 0.0 ) r = 0.0 ; else if ( r > 1.0 ) r = 1.0 ;
00155    if( g < 0.0 ) g = 0.0 ; else if ( g > 1.0 ) g = 1.0 ;
00156    if( b < 0.0 ) b = 0.0 ; else if ( b > 1.0 ) b = 1.0 ;
00157 
00158    active_color = (float) RGB_TO_COL(r,g,b) ;
00159    return ;
00160 }
00161 
00162 /*----- This routine is called from color.f -----*/
00163 
00164 void zzmpco_( float * r , float * g , float * b )
00165 {
00166    set_color_memplot( *r , *g , *b ) ;
00167    return ;
00168 }
00169 
00170 void set_thick_memplot( float th )
00171 {
00172    if( th < 0.0 ) th = 0.0 ;
00173    active_thick = th ;
00174    return ;
00175 }
00176 
00177 float get_thick_memplot( void )
00178 {
00179    return active_thick ;
00180 }
00181 
00182 void set_opacity_memplot( float th )  /* 22 Jul 2004 */
00183 {
00184    MEM_plotdata *mp ;
00185 
00186         if( th < 0.0 ) th = 0.0 ;
00187    else if( th > 1.0 ) th = 1.0 ;
00188    active_opacity = th ;
00189 
00190    /* Set opacity for further drawing [22 Jul 2004] */
00191 
00192    if( active_plot < 0 || active_plot >= num_plotar ||
00193        num_plotar == 0 || plotar == NULL            ||
00194        plotar[active_plot] == NULL                    ) return ;
00195 
00196    mp = plotar[active_plot] ;
00197    ADDTO_MEMPLOT( mp , th,0.0,0.0,0.0 , 0.0 , -THCODE_OPAC ) ;
00198    return ;
00199 }
00200 
00201 float get_opacity_memplot( void )
00202 {
00203    return active_opacity ;
00204 }
00205 
00206 /*------------------------------------------------------------------
00207   where the actual plotting into the memplot is done from
00208   the coxplot functions
00209 --------------------------------------------------------------------*/
00210 
00211 void plotline_memplot( float x1 , float y1 , float x2 , float y2 )
00212 {
00213    MEM_plotdata * mp ;
00214 
00215    if( active_plot < 0 || active_plot >= num_plotar ||
00216        num_plotar == 0 || plotar == NULL            ||
00217        plotar[active_plot] == NULL                    ) return ;
00218 
00219    mp = plotar[active_plot] ;
00220 
00221 #if 0
00222 fprintf(stderr,"** plotline_memplot %d: (%f,%f) to (%f,%f)\n",
00223         MEMPLOT_NLINE(mp) , x1,y1,x2,y2) ;
00224 #endif
00225 
00226    ADDTO_MEMPLOT( mp , x1,y1,x2,y2 , active_color , active_thick ) ;
00227    return ;
00228 }
00229 
00230 void plotrect_memplot( float x1 , float y1 , float x2 , float y2 ) /* 21 Mar 2001 */
00231 {
00232    MEM_plotdata * mp ;
00233 
00234    if( active_plot < 0 || active_plot >= num_plotar ||
00235        num_plotar == 0 || plotar == NULL            ||
00236        plotar[active_plot] == NULL                    ) return ;
00237 
00238    mp = plotar[active_plot] ;
00239 
00240    ADDTO_MEMPLOT( mp , x1,y1,x2,y2 , active_color , -THCODE_RECT ) ;
00241    return ;
00242 }
00243 
00244 void plotcirc_memplot( float x1 , float y1 , float rad ) /* 10 Mar 2002 */
00245 {
00246    MEM_plotdata * mp ;
00247 
00248    if( active_plot < 0 || active_plot >= num_plotar ||
00249        num_plotar == 0 || plotar == NULL            ||
00250        plotar[active_plot] == NULL                    ) return ;
00251 
00252    mp = plotar[active_plot] ;
00253 
00254    ADDTO_MEMPLOT( mp , x1,y1,rad,0.0 , active_color , -THCODE_CIRC ) ;
00255    return ;
00256 }
00257 
00258 /*----- This routine is called from zzphph.f to draw 1 actual line -----*/
00259 
00260 void zzmpli_( float * x1 , float * y1 , float * x2 , float * y2 )
00261 {
00262    plotline_memplot( *x1 , *y1 , *x2 , *y2 ) ;
00263    return ;
00264 }
00265 
00266 /*------------------------------------------------------------------------
00267    Delete the active in-memory plot.
00268    After this, there is no "active" plot.
00269 --------------------------------------------------------------------------*/
00270 
00271 void delete_active_memplot(void)
00272 {
00273    int ip ;
00274 
00275    if( active_plot < 0 || active_plot >= num_plotar ||
00276        num_plotar == 0 || plotar == NULL            ||
00277        plotar[active_plot] == NULL                    ) return ;
00278 
00279    DESTROY_MEMPLOT( plotar[active_plot] ) ;
00280 
00281    if( num_plotar == 1 ){
00282       free(plotar) ; plotar = NULL ; num_plotar = 0 ;
00283    } else {
00284       for( ip=active_plot+1 ; ip < num_plotar ; ip++ ) plotar[ip-1] = plotar[ip] ;
00285       num_plotar-- ; plotar[num_plotar] = NULL ;
00286    }
00287 
00288    active_plot = -1 ;
00289    return ;
00290 }
00291 
00292 /*------------------------------------------------------------------------*/
00293 
00294 void delete_memplot( MEM_plotdata * mp )
00295 {
00296    int ip ;
00297 
00298    if( num_plotar == 0 || plotar == NULL || mp == NULL ) return ;
00299 
00300    for( ip=0 ; ip < num_plotar ; ip++ ) if( plotar[ip] == mp ) break ;
00301 
00302    if( ip < num_plotar ){
00303            if( active_plot == ip ) active_plot = -1 ;
00304       else if( active_plot >  ip ) active_plot-- ;
00305 
00306       for( ip++ ; ip < num_plotar ; ip++ ) plotar[ip-1] = plotar[ip] ;
00307 
00308       num_plotar-- ; plotar[num_plotar] = NULL ;
00309    }
00310 
00311    DESTROY_MEMPLOT( mp ) ;
00312    return ;
00313 }
00314 
00315 /*-----------------------------------------------------------------------
00316    Scale data inside an memplot -- 26 Feb 2001
00317       x_new     = sx * x_old + tx
00318       y_new     = sy * y_old + ty
00319       thick_new = st * thick_old
00320 -------------------------------------------------------------------------*/
00321 
00322 void scale_memplot( float sx , float tx ,
00323                     float sy , float ty , float st , MEM_plotdata * mp )
00324 {
00325    int ii,nn ;
00326    if( mp == NULL ) return ;
00327 
00328    for( nn=ii=0 ; ii < mp->nxyline ; ii++ ){
00329       mp->xyline[nn] = mp->xyline[nn] * sx + tx ; nn++ ; /* x1 */
00330       mp->xyline[nn] = mp->xyline[nn] * sy + ty ; nn++ ; /* y1 */
00331       mp->xyline[nn] = mp->xyline[nn] * sx + tx ; nn++ ; /* x2 */
00332       mp->xyline[nn] = mp->xyline[nn] * sy + ty ; nn++ ; /* y2 */
00333                                                   nn++ ; /* color */
00334       if( mp->xyline[nn] > 0.0 )
00335         mp->xyline[nn] = mp->xyline[nn] * st    ; nn++ ; /* thick */
00336    }
00337    return ;
00338 }
00339 
00340 /*-----------------------------------------------------------------------
00341    Append data from one memplot to another -- 26 Feb 2001
00342 -------------------------------------------------------------------------*/
00343 
00344 void append_to_memplot( MEM_plotdata * mp , MEM_plotdata * ap )
00345 {
00346    int nn , nold ;
00347    if( mp == NULL || ap == NULL || ap->nxyline <= 0 ) return ;
00348 
00349    nn = mp->nxyline + ap->nxyline ;
00350    mp->xyline = (float *) realloc(mp->xyline,
00351                                   sizeof(float)*NXY_MEMPLOT*nn) ;
00352 
00353    memcpy( mp->xyline + NXY_MEMPLOT*mp->nxyline ,
00354            ap->xyline , sizeof(float)*NXY_MEMPLOT*ap->nxyline ) ;
00355 
00356    mp->nxyline = mp->nxyline_all = nn ;
00357    return ;
00358 }
00359 
00360 /*-----------------------------------------------------------------------
00361    Make a copy of a memplot; the new one will be the active memplot
00362    -- 26 Feb 2001 -- RWCox
00363 -------------------------------------------------------------------------*/
00364 
00365 MEM_plotdata * copy_memplot( MEM_plotdata * mp )
00366 {
00367    MEM_plotdata * np ;
00368    char str[256] ; int nn ;
00369 
00370    if( mp == NULL ) return NULL ;
00371 
00372    /* make a new ID string */
00373 
00374    for( nn=1 ; nn <= 9999 ; nn++ ){
00375       sprintf(str,"%.240sCopy%04d",mp->ident,nn) ;
00376       if( find_memplot(str) == NULL ) break ;
00377    }
00378    if( nn == 1000 ) return NULL ; /* this is bad (but unlikely) */
00379 
00380    /* make the new memplot */
00381 
00382    nn = create_memplot( str , mp->aspect ) ;
00383    if( nn ) return NULL ;         /* this is real bad */
00384 
00385    np = find_memplot(NULL) ;      /* is the new one */
00386    if( np == NULL ) return NULL ; /* shouldn't happen */
00387 
00388    /* copy data from old one into new one */
00389 
00390    nn = np->nxyline = np->nxyline_all = mp->nxyline ;
00391    np->xyline = (float *) realloc(np->xyline,
00392                                   sizeof(float)*NXY_MEMPLOT*nn) ;
00393    memcpy( np->xyline , mp->xyline , sizeof(float)*NXY_MEMPLOT*nn ) ;
00394 
00395    return np ;
00396 }
00397 
00398 /*----------------------------------------------------------------------
00399    Flip a memplot inplace - 30 Aug 2001 - RWCox
00400      rot    = one of the MRI_ROT_ codes (see coxplot.h)
00401      mirror = whether to left-right mirror after rotation
00402 ------------------------------------------------------------------------*/
00403 
00404 void flip_memplot( int rot , int mirror , MEM_plotdata *mp )
00405 {
00406    int fopt , ii,nn ;
00407    float xtop , ytop=1.0 , x1,y1,x2,y2 ;
00408 
00409    if( mp == NULL ) return ;                          /* nothing in */
00410    if( rot == MRI_ROT_0 && mirror == FALSE ) return ; /* do nothing */
00411 
00412    xtop = mp->aspect ;
00413 
00414    fopt = (mirror) ? (rot+MRI_FLMADD) : (rot) ;
00415    switch( fopt ){
00416 
00417       default: return ;  /* should never happen */
00418 
00419       case MRI_ROT_90:
00420        for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00421           x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00422           x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00423           mp->xyline[nn  ] = ytop - y1 ;
00424           mp->xyline[nn+1] = x1 ;
00425           mp->xyline[nn+2] = ytop - y2 ;
00426           mp->xyline[nn+3] = x2 ;
00427        }
00428       break ;
00429 
00430       case MRI_ROT_180:
00431        for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00432           x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00433           x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00434           mp->xyline[nn  ] = xtop - x1 ;
00435           mp->xyline[nn+1] = ytop - y1 ;
00436           mp->xyline[nn+2] = xtop - x2 ;
00437           mp->xyline[nn+3] = ytop - y2 ;
00438        }
00439       break ;
00440 
00441       case MRI_ROT_270:
00442        for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00443           x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00444           x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00445           mp->xyline[nn  ] = y1 ;
00446           mp->xyline[nn+1] = xtop - x1 ;
00447           mp->xyline[nn+2] = y2 ;
00448           mp->xyline[nn+3] = xtop - x2 ;
00449        }
00450       break ;
00451 
00452       case (MRI_ROT_0+MRI_FLMADD):
00453        for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00454           x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00455           x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00456           mp->xyline[nn  ] = xtop - x1 ;
00457           mp->xyline[nn+1] = y1 ;
00458           mp->xyline[nn+2] = xtop - x2 ;
00459           mp->xyline[nn+3] = y2 ;
00460        }
00461       break ;
00462 
00463       case (MRI_ROT_90+MRI_FLMADD):
00464        for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00465           x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00466           x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00467           mp->xyline[nn  ] = y1 ;
00468           mp->xyline[nn+1] = x1 ;
00469           mp->xyline[nn+2] = y2 ;
00470           mp->xyline[nn+3] = x2 ;
00471        }
00472       break ;
00473 
00474       case (MRI_ROT_180+MRI_FLMADD):
00475        for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00476           x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00477           x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00478           mp->xyline[nn  ] = x1 ;
00479           mp->xyline[nn+1] = ytop - y1 ;
00480           mp->xyline[nn+2] = x2 ;
00481           mp->xyline[nn+3] = ytop - y2 ;
00482        }
00483       break ;
00484 
00485       case (MRI_ROT_270+MRI_FLMADD):
00486        for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00487           x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00488           x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00489           mp->xyline[nn  ] = ytop - y1 ;
00490           mp->xyline[nn+1] = xtop - x1 ;
00491           mp->xyline[nn+2] = ytop - y2 ;
00492           mp->xyline[nn+3] = xtop - x2 ;
00493        }
00494       break ;
00495    }
00496 
00497    return ;
00498 }
00499 
00500 /*---------------------------------------------------------------------------
00501   Set the insertion point for over-writing -- 15 Nov 2001
00502 -----------------------------------------------------------------------------*/
00503 
00504 void insert_at_memplot( int ii , MEM_plotdata *mp )
00505 {
00506    if( mp != NULL ) mp->insert_at = ii ;
00507    return ;
00508 }
00509 
00510 /*---------------------------------------------------------------------------
00511   Cut lines nbot..ntop out of a memplot -- 15 Nov 2001
00512 -----------------------------------------------------------------------------*/
00513 
00514 void cutlines_memplot( int nbot , int ntop , MEM_plotdata *mp )
00515 {
00516    if( mp == NULL          ) return ;  /* bad or meaningless stuff */
00517    if( nbot <  0           ) return ;
00518    if( ntop >= mp->nxyline ) return ;
00519    if( nbot > ntop         ) return ;
00520 
00521    if( ntop == mp->nxyline-1 ){  /* just set num lines to nbot */
00522 
00523       mp->nxyline = nbot ;
00524 
00525    } else {                      /* must move things above ntop down */
00526 
00527       memmove( mp->xyline + NXY_MEMPLOT*nbot ,
00528                mp->xyline + NXY_MEMPLOT*(ntop+1) ,
00529                sizeof(float)*NXY_MEMPLOT*(mp->nxyline-1-ntop) ) ;
00530 
00531       mp->nxyline -= (ntop-nbot+1) ;
00532 
00533    }
00534    return ;
00535 }
00536 
00537 /*----------------------------------------------------------------------------*/
00538 
00539 #ifdef __GNUC__
00540 # define INLINE inline
00541 #else
00542 # define INLINE /*nada*/
00543 #endif
00544 
00545 /*----------------------------------------------------------------------------
00546   Clip a line to a rectangle.  Return is -1 if the line is totally outside.
00547   Otherwise, return is 0 and *x1in (etc.) is altered to the clipped line.
00548 ------------------------------------------------------------------------------*/
00549 
00550 static INLINE int clip_line_to_rect( float xclbot , float yclbot ,
00551                                      float xcltop , float ycltop ,
00552                                      float *x1in  , float *y1in  ,
00553                                      float *x2in  , float *y2in   )
00554 {
00555    float x1=*x1in , y1=*y1in , x2=*x2in , y2=*y2in , dx,dy,slope,temp ;
00556    int inter=0 ;
00557 
00558    /* Make sure that x1 < x2 by interchanging the points if necessary */
00559 
00560    if( x1 > x2 ){
00561      temp=x1 ; x1=x2 ; x2=temp;
00562      temp=y1 ; y1=y2 ; y2=temp; inter=1 ;
00563    }
00564 
00565    /* if outside entire region, throw line away */
00566 
00567    if( x2 < xclbot || x1 > xcltop ) return -1;
00568 
00569    if( y1 < y2 ){
00570      if( y2 < yclbot || y1 > ycltop ) return -1;
00571    } else {
00572      if( y1 < yclbot || y2 > ycltop ) return -1;
00573    }
00574 
00575    /* if inside entire region, then do nothing */
00576 
00577    if( x1 >= xclbot && x2 <= xcltop ){
00578      if( y1 < y2 ){
00579        if( y1 >= yclbot && y2 <= ycltop ) return 0 ;
00580      } else {
00581        if( y2 >= yclbot && y1 <= ycltop ) return 0 ;
00582      }
00583    }
00584 
00585    /* Clip line in X direction */
00586 
00587    dx = x2 - x1 ;
00588    if( dx > 0.0 ){  /* only clip if line has some x range */
00589      slope = (y2-y1)/dx ;
00590      if( x1 < xclbot ){  /* intercept of line at left side */
00591        y1 = y1 + slope*(xclbot-x1) ;
00592        x1 = xclbot ;
00593      }
00594      if( x2 > xcltop ){  /* intercept at right */
00595        y2 = y2 + slope*(xcltop-x2) ;
00596        x2 = xcltop ;
00597      }
00598    }
00599 
00600    /* Check line again to see if it falls outside of plot region */
00601 
00602    if( y1 < y2 ){
00603      if( y2 < yclbot || y1 > ycltop ) return -1;
00604    } else {
00605      if( y1 < yclbot || y2 > ycltop ) return -1;
00606 
00607      temp=x1 ; x1=x2 ; x2=temp;                 /* make sure y1 <= y2 */
00608      temp=y1 ; y1=y2 ; y2=temp; inter=!inter ;
00609    }
00610 
00611    /* Clip y-direction.  To do this, must have y1 <= y2 [supra] */
00612 
00613    dy = y2 - y1 ;
00614    if( dy > 0.0 ){  /* only clip if line has some Y range */
00615      slope = (x2-x1)/dy ;
00616      if( y1 < yclbot ){ /* intercept of line at bottom */
00617        x1 = x1 + slope*(yclbot-y1) ;
00618        y1 = yclbot ;
00619      }
00620      if( y2 > ycltop ){ /* intercept at top */
00621        x2 = x2 + slope*(ycltop-y2) ;
00622        y2 = ycltop ;
00623      }
00624    }
00625 
00626    /* Line is now guaranteed to be totally inside the plot region.
00627       Copy local clipped coordinates to output values and return.
00628       Note that we must restore points to original input order,
00629       if they were interchanged at some point above.              */
00630 
00631    if( inter ){
00632      *x1in = x2 ; *x2in = x1 ; *y1in = y2 ; *y2in = y1 ;
00633    } else {
00634      *x1in = x1 ; *y1in = y1 ; *x2in = x2 ; *y2in = y2 ;
00635    }
00636 
00637    return 0 ;
00638 }
00639 
00640 #undef INSIDE
00641 #define INSIDE(x,y)                                                    \
00642   ( (x) >= xclbot && (x) <= xcltop && (y) >= yclbot && (y) <= ycltop )
00643 
00644 /*---------------------------------------------------------------------------
00645   Clip a memplot to a rectangle, producing a new memplot.
00646 -----------------------------------------------------------------------------*/
00647 
00648 MEM_plotdata * clip_memplot( float xclbot, float yclbot,
00649                              float xcltop, float ycltop , MEM_plotdata *mp )
00650 {
00651    MEM_plotdata *np ;
00652    char str[256] ;
00653    int nn , ii , qq ;
00654    float x1,y1 , x2,y2 , col,th ;
00655 
00656    if( mp == NULL       ) return NULL ;  /* bad or meaningless stuff */
00657    if( xclbot >= xcltop ) return NULL ;
00658    if( yclbot >= ycltop ) return NULL ;
00659 
00660    sprintf(str,"%.240sCopy",mp->ident) ;
00661    nn = create_memplot_surely( str , mp->aspect ) ;
00662    np = find_memplot(NULL) ;
00663    if( np == NULL ) return NULL ; /* shouldn't happen */
00664 
00665    for( nn=ii=0 ; ii < mp->nxyline ; ii++,nn+=NXY_MEMPLOT ){
00666      x1 = mp->xyline[nn  ] ; y1 = mp->xyline[nn+1] ;
00667      x2 = mp->xyline[nn+2] ; y2 = mp->xyline[nn+3] ;
00668      col= mp->xyline[nn+4] ; th = mp->xyline[nn+5] ;
00669 
00670      if( th < 0.0 ){               /** Not a line! */
00671        int thc = (int)(-th) ;
00672        switch( thc ){
00673          case THCODE_RECT:         /* rectangle */
00674                                    /* both corners inside */
00675            if( INSIDE(x1,y1) && INSIDE(x2,y2) ){
00676              ADDTO_MEMPLOT(np,x1,y1,x2,y2,col,th) ;
00677            }
00678          break ;
00679 
00680          case THCODE_CIRC:{        /* circle */
00681                                    /* +/- 1 radius inside */
00682            float xx,yy , rr=x2 ;
00683            xx = x1+rr ; if( !INSIDE(xx,y1) ) break ;
00684            xx = x1-rr ; if( !INSIDE(xx,y1) ) break ;
00685            yy = y1+rr ; if( !INSIDE(x1,yy) ) break ;
00686            yy = y1-rr ; if( !INSIDE(x1,yy) ) break ;
00687            ADDTO_MEMPLOT(np,x1,y1,x2,y2,col,th) ;
00688          }
00689          break ;
00690        }
00691 
00692      } else {                      /** Truly a line! **/
00693 
00694        qq = clip_line_to_rect( xclbot,yclbot , xcltop,ycltop ,
00695                                &x1,&y1       , &x2,&y2        ) ;
00696        if( qq == 0 ){
00697          ADDTO_MEMPLOT(np,x1,y1,x2,y2,col,th) ;
00698        }
00699      }
00700    }
00701 
00702    if( np->nxyline == 0 ) DESTROY_MEMPLOT(np) ;
00703 
00704    return np ;
00705 }
00706 
00707 /****************************************************************************
00708   Functions to interface with PLOTPAK Fortran routines
00709 *****************************************************************************/
00710 
00711 /*-------------------------------
00712   Has no function at this time.
00713 ---------------------------------*/
00714 void plotpak_frame(void) { frame_() ; }
00715 
00716 /*-----------------------------------------------
00717    Draws a sequence of lines in one swell foop.
00718 -------------------------------------------------*/
00719 void plotpak_curve( float * x , float * y , int n )
00720 {
00721    integer nn = n ;
00722    curve_( (real *) x , (real *) y , &nn ) ;
00723 }
00724 
00725 /*---------------------------------------------------
00726   Establishes first point of a series of lines
00727   to be drawn one at time using plotpak_vector.
00728 -----------------------------------------------------*/
00729 void plotpak_frstpt( float x , float y )
00730 {
00731    real xx = x , yy = y ;
00732    frstpt_( &xx , &yy ) ;
00733 }
00734 
00735 /*---------------------------------------------------
00736   Draws next in the series of lines.
00737 -----------------------------------------------------*/
00738 void plotpak_vector( float x , float y )
00739 {
00740    real xx=x , yy=y ;
00741    vector_( &xx , &yy ) ;
00742 }
00743 
00744 /*-------------------
00745    Draws one line.
00746 ---------------------*/
00747 void plotpak_line( float x1 , float y1 , float x2 , float y2 )
00748 {
00749    real xx1 = x1 , yy1 = y1 , xx2 = x2 , yy2 = y2 ;
00750    line_(&xx1, &yy1, &xx2, &yy2);
00751 }
00752 
00753 /*--------------------------------------------------------------
00754  Converts (x1,y1) from user to memplot coordinates - 20 Nov 2001
00755 ----------------------------------------------------------------*/
00756 
00757 void plotpak_zzphys( float x1 , float y1 , float *x2 , float *y2 )
00758 {
00759    real xx1 = x1 , yy1 = y1 ;
00760    zzphys_( &xx1 , &yy1 ) ;
00761    if( x2 != NULL ) *x2 = xx1 ;
00762    if( y2 != NULL ) *y2 = yy1 ;
00763 }
00764 
00765 /*--------------------------------------------------------------
00766  Converts (x1,y1) from memplot to user coordinates - 20 Nov 2001
00767 ----------------------------------------------------------------*/
00768 
00769 void plotpak_unphys( float x1 , float y1 , float *x2 , float *y2 )
00770 {
00771    double rr ;
00772    if( x2 != NULL ){
00773       rr = (x1 - zzzplt_.betaxx) / zzzplt_.alphxx ;
00774       if( zzzplt_.ixcoor < 0 ) rr = pow(10.0,rr) ;
00775       *x2 = rr ;
00776    }
00777    if( y2 != NULL ){
00778       rr = (y1 - zzzplt_.betayy) / zzzplt_.alphyy ;
00779       if( zzzplt_.iycoor < 0 ) rr = pow(10.0,rr) ;
00780       *y2 = rr ;
00781    }
00782 }
00783 
00784 /*-----------------------------------------------------------------------------
00785   Establishes relationship between objective (memplot) coordinates
00786   (x range is from 0 to aspect, y from 0 to 1.0) and subjective (user)
00787   coordinates.  The range xs1..xs2 is mapped onto xo1..xo2, and similarly
00788   for y.  "code" establishes the form of the mapping:
00789     code = 1 => x linear, y linear
00790     code = 2 => x linear, y logarithmic
00791     code = 3 => x log   , y linear
00792     code = 4 => x log   , y log
00793 -------------------------------------------------------------------------------*/
00794 void plotpak_set( float xo1,float xo2 , float yo1,float yo2 ,
00795                   float xs1,float xs2 , float ys1,float ys2 , int code )
00796 {
00797    real xobj1=xo1, xobj2=xo2, yobj1=yo1, yobj2=yo2;
00798    real xsub1=xs1, xsub2=xs2, ysub1=ys1, ysub2=ys2 ;
00799    integer ltype = code ;
00800    set_(&xobj1, &xobj2, &yobj1, &yobj2, &xsub1, &xsub2, &ysub1, &ysub2, &ltype);
00801 }
00802 
00803 void plotpak_getset( float *xo1,float *xo2 , float *yo1,float *yo2 ,
00804                      float *xs1,float *xs2 , float *ys1,float *ys2  )
00805 {
00806    if( xo1 != NULL ) *xo1 = (float) zzzplt_.xbot ;
00807    if( xo2 != NULL ) *xo2 = (float) zzzplt_.xtop ;
00808    if( yo1 != NULL ) *yo1 = (float) zzzplt_.ybot ;
00809    if( yo1 != NULL ) *yo2 = (float) zzzplt_.ytop ;
00810 
00811    if( xs1 != NULL ) *xs1 = (float) zzzplt_.xmin ;
00812    if( xs2 != NULL ) *xs2 = (float) zzzplt_.xmax ;
00813    if( ys1 != NULL ) *ys1 = (float) zzzplt_.ymin ;
00814    if( ys1 != NULL ) *ys2 = (float) zzzplt_.ymax ;
00815 
00816    return ;
00817 }
00818 
00819 /*-----------------------------------------
00820   Set line type: 1 = solid (default)
00821                  2 = long dash
00822                  3 = short dash
00823                  4 = long - short - short
00824                  5 = very short
00825 -------------------------------------------*/
00826 void plotpak_setlin( int code )
00827 {
00828    integer ntype=code ;
00829    setlin_(&ntype);
00830 }
00831 
00832 /*------------------------------------------------
00833   Set the clipping window (objective coordinates)
00834 --------------------------------------------------*/
00835 void plotpak_setw( float xo1,float xo2 , float yo1,float yo2 )
00836 {
00837    real xobj1=xo1, xobj2=xo2, yobj1=yo1, yobj2=yo2;
00838    setw_( &xobj1, &xobj2, &yobj1, &yobj2 ) ;
00839 }
00840 
00841 /*------- Plotpak routines that I'm not documenting yet
00842           see the .f source code, or NCAR manual, if you care -------*/
00843 
00844 void plotpak_setdsh( int nd , float * xd )
00845 {
00846    integer nnd = nd ;
00847    setdsh_( &nnd , (real *) xd ) ;
00848 }
00849 
00850 void plotpak_setfrm( float xo1,float xo2 , float yo1,float yo2 )
00851 {
00852    real xobj1=xo1, xobj2=xo2, yobj1=yo1, yobj2=yo2;
00853    setfrm_( &xobj1, &xobj2, &yobj1, &yobj2 ) ;
00854 }
00855 
00856 void plotpak_phdot( float x1 , float y1 )
00857 {
00858    real xx1=x1 , yy1=y1 ;
00859    phdot_(&xx1,&yy1);
00860 }
00861 
00862 void plotpak_phline( float x1 , float y1 , float x2 , float y2 )
00863 {
00864    real xx1 = x1 , yy1 = y1 , xx2 = x2 , yy2 = y2 ;
00865    phline_(&xx1, &yy1, &xx2, &yy2);
00866 }
00867 
00868 void plotpak_point( float x1 , float y1 )
00869 {
00870    real xx1=x1 , yy1=y1 ;
00871    point_(&xx1,&yy1);
00872 }
00873 
00874 void plotpak_points( float *x , float *y , int n , int ipen )
00875 {
00876    integer nn=n , nipen=ipen , zero=0 ;
00877    points_( (real *)x , (real *)y , &nn , &zero , &nipen ) ;
00878 }
00879 
00880 void ppak_garbage_routine(void) ;
00881 void this_is_real_junk(void){ ppak_garbage_routine(); }