src/xmaxima_def.lisp

   1 ;; xmaxima.lisp: routines for Maxima's interface to xmaxima
   2 ;; Copyright (C) 2007-2021 J. Villate
   3 ;; Time-stamp: "2024-03-25 08:58:31 villate"
   4 ;;
   5 ;; This program is free software; you can redistribute it and/or
   6 ;; modify it under the terms of the GNU General Public License
   7 ;; as published by the Free Software Foundation; either version 2
   8 ;; of the License, or (at your option) any later version.
   9 ;;
  10 ;; This program is distributed in the hope that it will be useful,
  11 ;; but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 ;; MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  13 ;; GNU General Public License for more details.
  14 ;;
  15 ;; You should have received a copy of the GNU General Public License
  16 ;; along with this program; if not, write to the Free Software
  17 ;; Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
  18 ;; MA  02110-1301, USA
  19
  20 (in-package :maxima)
  21
  22 ;; Given a list of valid colors (see rgb-color function) and an object c
  23 ;; that can be a real number or a string, produces an rgb color
  24 ;; specification for c; when c is real, its nearest integer is assigned
  25 ;; to one of the numbers in the list, using modulo length of the list.
  26 (defun xmaxima-color (colors c)
  27   (unless (listp colors) (setq colors (list colors)))
  28   (when (realp c)
  29     (unless (integerp c) (setq c (round c)))
  30     (setq c (nth (mod (1- c) (length colors)) colors)))
  31   (rgb-color c))
  32
  33 ;; style is a list starting with a symbol from the list: points, lines,
  34 ;; linespoints or dots,
  35 ;; The meaning of the numbers that follow the symbol are:
  36 ;;
  37 ;;   lines, linewidth, color
  38 ;;   points, radius, color
  39 ;;   linespoints, linewidth, radius, color
  40 ;;   dots, color
  41 ;;
  42 ;; linewidth and radius are measured in the same units and can be
  43 ;; floating-point numbers.
  44 ;;
  45 ;; type must be an integer
  46 ;; color can be an integer, used as index to get one of the colors defined
  47 ;; by the color option, or a 6-digit hexadecimal number #rrggbb
  48
  49 (defun xmaxima-curve-style (style colors i)
  50   (unless (listp style) (setq style (list style)))
  51   (unless (listp colors) (setq colors (list colors)))
  52   (with-output-to-string
  53     (st)
  54     (case (first style)
  55       ($dots
  56        (format st "\{ nolines 1 \} \{ plotpoints 1 \} \{ pointsize 0.7 \}")
  57        (if (second style)
  58          (format st " \{ color ~a \}" (xmaxima-color colors (second style)))
  59          (format st " \{ color ~a \}" (xmaxima-color colors i))))
  60       ($lines
  61        (format st "\{ nolines 0 \} \{ plotpoints 0 \}")
  62        (if (realp (second style))
  63          (format st " \{ linewidth ~f \}" (second style)))
  64        (if (third style)
  65          (format st " \{ color ~a \}" (xmaxima-color colors (third style)))
  66          (format st " \{ color ~a \}" (xmaxima-color colors i))))
  67       ($points
  68        (format st "\{ nolines 1 \} \{ plotpoints 1 \}")
  69        (if (realp (second style))
  70          (format st " \{ pointsize ~f \}" (second style))
  71          (format st " \{ pointsize 3 \}"))
  72        (if (third style)
  73          (format st " \{ color ~a \}" (xmaxima-color colors (third style)))
  74          (format st " \{ color ~a \}" (xmaxima-color colors i))))
  75       ($linespoints
  76        (format st "\{ nolines 0 \} \{ plotpoints 1 \}")
  77        (if (realp (second style))
  78          (format st " \{ linewidth ~f \}" (second style)))
  79        (if (realp (third style))
  80          (format st " \{ pointsize ~f \}" (third style))
  81          (format st " \{ pointsize 3 \}"))
  82        (if (fourth style)
  83          (format st " \{ color ~a \}" (xmaxima-color colors (fourth style)))
  84          (format st " \{ color ~a \}" (xmaxima-color colors i))))
  85       (t
  86        (format st "\{ nolines 0 \} \{ plotpoints 0 \} \{ color ~a \}"
  87                (xmaxima-color colors i))))))
  88
  89 (defun xmaxima-palette (palette)
  90 ;; palette should be a list starting with one of the symbols: hue,
  91 ;; saturation, value, gray or gradient.
  92 ;;
  93 ;; If the symbol is gray, it should be followed by two floating point
  94 ;; numbers that indicate the initial gray level and the interval of
  95 ;; gray values.
  96 ;;
  97 ;; If the symbol is one of hue, saturation or value, it must be followed
  98 ;; by three numbers that specify the hue, saturation and value for the
  99 ;; initial color, and a fourth number that gives the range of values for
 100 ;; the increment of hue, saturation or value.
 101 ;; The values for the initial hue, saturation, value and grayness should
 102 ;; be within 0 and 1, while the range can be higher or even negative.
 103 ;;
 104 ;; If the symbol is gradient, it must be followed by either a list of valid
 105 ;; colors or by a list of lists with two elements, a number and a valid color.
 106
 107   (unless (listp palette) (setq palette (list palette)))
 108   (let (hue sat val gray range fun)
 109     (case (first palette)
 110       ($gray
 111        (case (length (rest palette))
 112          (2 (setq gray (second palette)) (setq range (third palette)))
 113          (t (merror
 114              (intl:gettext
 115               "palette: gray must be followed by two numbers."))))
 116        (when (or (< gray 0) (> gray 1))
 117          (setq gray (- gray (floor gray))))
 118        (setq fun (format nil "{value ~,,,,,,'eg} {colorrange ~,,,,,,'eg}"
 119                          gray range)))
 120       (($hue $saturation $value)
 121        (case (length (rest palette))
 122          (4 (setq hue (second palette))
 123             (setq sat (third palette))
 124             (setq val (fourth palette))
 125             (setq range (fifth palette)))
 126          (t (merror
 127              (intl:gettext
 128               "palette: ~M must be followed by four numbers.")
 129               (first palette))))
 130        (when (or (< hue 0) (> hue 1)) (setq hue (- hue (floor hue))))
 131        (when (or (< sat 0) (> sat 1)) (setq sat (- sat (floor sat))))
 132        (when (or (< val 0) (> val 1)) (setq val (- val (floor val))))
 133        (setq fun
 134              (format nil " {hue ~,,,,,,'eg} {saturation ~,,,,,,'eg} {value ~,,,,,,'eg} {colorrange ~,,,,,,'eg}"
 135                      hue sat val range))))
 136     (with-output-to-string (st)
 137       (case (first palette)
 138         ($hue (format st "~&~a {colorscheme hue}" fun))
 139         ($saturation (format st "~&~a {colorscheme saturation}" fun))
 140         ($value (format st "~&~a {colorscheme value}" fun))
 141         ($gray (format st "~&~a {colorscheme gray}" fun))
 142         ($gradient
 143          (let* ((colors (rest palette)) (n (length colors)) (map nil))
 144            ;; map is constructed as (n1 c1 n2 c2 ... nj cj) where ni is a
 145            ;; decreasing sequence of numbers (n1=1, nj=0) and ci are colors
 146            (cond
 147              ;; Maxima list of numbers and colors (((mlist) ni ci) ...)
 148              ((listp (first colors))
 149               (setq colors (sort colors #'< :key #'cadr))
 150               (dotimes (i n)
 151                 (setq map (cons (rgb-color (third (nth i colors))) ;; color
 152                                 (cons
 153                                  (/ (- (second (nth i colors))   ;; ni minus
 154                                        (second (first colors)))  ;; smallest ni
 155                                     (- (second (nth (- n 1) colors));; biggest
 156                                        (second (first colors)))) ;; - smallest
 157                                  map)))))
 158              ;; list of only colors
 159              (t (dotimes (i n)
 160                   (setq map (cons (rgb-color (nth i colors))  ;; color i
 161                                   (cons (/ i (1- n)) map))))))    ;; number i
 162
 163            ;; prints map with the format:  nj, "cj", ...,n1, "c1"
 164            (setq fun (format nil "~{{ ~,,,,,,'eg ~s}~^ ~}" (reverse map)))
 165            (format st "~&{colorscheme gradient} ")
 166            ;; writes: {gradlist {{nj "cj"} ...{n1 "c1"}}}
 167            (format st "{gradlist {~a}}" fun)))
 168         (t
 169          (merror
 170           (intl:gettext
 171            "palette: wrong keyword ~M. Must be hue, saturation, value, gray or gradient.")
 172           (first palette)))))))
 173
 174 (defun xmaxima-palettes (palette n)
 175   (unless (integerp n) (setq n (round n)))
 176   (if (find 'mlist palette :key #'car) (setq palette (list palette)))
 177   (xmaxima-palette (rest (nth (mod (- n 1) (length palette)) palette))))
 178
 179 (defmethod plot-preamble ((plot xmaxima-plot) plot-options)
 180   (let (outfile zmin zmax)
 181     (setf
 182      (slot-value plot 'data)
 183      (concatenate
 184       'string
 185       (slot-value plot 'data)
 186       (with-output-to-string (dest)
 187         (cond ($show_openplot
 188                (format dest "~a -data {~%" (getf plot-options '$type)))
 189               (t (format dest "{~a " (getf plot-options '$type))))
 190         (when (string= (getf plot-options '$type) "plot3d")
 191           (let ((palette (getf plot-options '$palette))
 192                 (meshcolor (if (member '$mesh_lines_color plot-options)
 193                                (getf plot-options '$mesh_lines_color)
 194                                '$black))
 195                 (elev (getf plot-options '$elevation))
 196                 (azim (getf plot-options '$azimuth)))
 197             (if (find 'mlist palette :key #'car) (setq palette (list palette)))
 198             (if palette
 199                 (progn
 200                   (if meshcolor
 201                       (format dest " {mesh_lines ~a}" (rgb-color meshcolor))
 202                       (format dest " {mesh_lines 0}")))
 203                 (format dest " {colorscheme 0}~%"))
 204             (when elev (format dest " {el ~d}" elev))
 205             (when azim (format dest " {az ~d}" azim))
 206             (format dest "~%")))
 207         (when (getf plot-options '$ps_file)
 208           (setq outfile (plot-file-path (getf plot-options '$ps_file) t))
 209           (format dest " {psfile ~s}" outfile))
 210         (when (member '$legend plot-options)
 211           (unless (getf plot-options '$legend)
 212             (format dest " {nolegend 1}")))
 213         (when (member '$box plot-options)
 214           (unless (getf plot-options '$box)
 215             (format dest " {nobox 1}")))
 216         (if (getf plot-options '$axes)
 217             (case (getf plot-options '$axes)
 218               ($x (format dest " {axes {x} }"))
 219               ($y (format dest " {axes {y} }"))
 220               (t (format dest " {axes {xy} }")))
 221             (format dest " {axes 0}"))
 222         (when (getf plot-options '$x)
 223           (format dest " {xrange ~{~,,,,,,'eg~^ ~}}" (getf plot-options '$x)))
 224         (when (getf plot-options '$y)
 225           (format dest " {yrange ~{~,,,,,,'eg~^ ~}}" (getf plot-options '$y)))
 226         (when (getf plot-options '$z)
 227           (setq zmin (first (getf plot-options '$z)))
 228           (setq zmax (second (getf plot-options '$z)))
 229           (format dest " {zcenter ~,,,,,,'eg }" (/ (+ zmax zmin) 2.0))
 230           (format dest " {zradius ~,,,,,,'eg }" (/ (- zmax zmin) 2.0)))
 231         (when (getf plot-options '$xlabel)
 232           (format dest " {xaxislabel ~s}" (getf plot-options '$xlabel)))
 233         (when (getf plot-options '$ylabel)
 234           (format dest " {yaxislabel ~s}" (getf plot-options '$ylabel)))
 235         (when (getf plot-options '$z)
 236           (format $pstream " {zcenter ~,,,,,,'eg }"
 237                   (/ (apply #'+ (getf plot-options '$z)) 2))
 238           (format $pstream " {zradius ~,,,,,,'eg }~%"
 239                   (/ (apply #'- (getf plot-options '$z)) -2)))
 240         (format dest "~%"))))
 241     ;;returns a list with the name of the file to be created, or nil
 242     (if (null outfile) nil (list outfile))))
 243
 244 (defmethod plot2d-command ((plot xmaxima-plot) fun options range)
 245   (let (points-lists)
 246     (setq points-lists
 247           (mapcar #'(lambda (f) (cdr (draw2d f range options))) (cdr fun)))
 248     (when (= (count-if #'(lambda (x) x) points-lists) 0)
 249       (merror (intl:gettext "plot2d: nothing to plot.~%")))
 250     (let ((legends-new) (legends (getf options '$legend)))
 251       (unless (null legends)
 252         (dotimes (i (length legends))
 253           (unless (null (cdr (nth i points-lists)))
 254             (push (nth i legends) legends-new)))
 255         (setf (getf options '$legend) (reverse legends-new))))
 256     (setf
 257      (slot-value plot 'data)
 258      (concatenate
 259       'string
 260       (slot-value plot 'data)
 261       (with-output-to-string (st)
 262         (unless (or (getf options '$logy)
 263                     (and (getf options '$y) (listp (getf options '$y))))
 264           (let (x y ymin ymax (xmin +most-negative-flonum+)
 265                   (xmax +most-positive-flonum+))
 266             (when (getf options '$x)
 267               (setq xmin (first (getf options '$x)))
 268               (setq xmax (second (getf options '$x))))
 269             (dolist (points-list points-lists)
 270               (dotimes (i (/ (length points-list) 2))
 271                 (setq x (nth (* i 2) points-list))
 272                 (setq y (nth (1+ (* i 2)) points-list))
 273                 (when (and (numberp x) (>= x xmin) (<= x xmax))
 274                   (when (numberp y)
 275                     (if (numberp ymin)
 276                         (if (numberp ymax)
 277                             (progn
 278                               (when (< y ymin) (setq ymin y))
 279                               (when (> y ymax) (setq ymax y)))
 280                           (if (< y ymin)
 281                               (setq ymax ymin ymin y)
 282                             (setq ymax y)))
 283                       (if (numberp ymax)
 284                           (if (> y ymax)
 285                               (setq ymin ymax ymax y)
 286                             (setq ymin y))
 287                         (setq ymin y)))))))
 288             (when (and (numberp ymin) (numberp ymax) (< ymin ymax))
 289               (psetq ymin (- (* 1.05 ymin) (* 0.05 ymax))
 290                      ymax (- (* 1.05 ymax) (* 0.05 ymin)))
 291               (format st " {yrange ~,,,,,,'eg ~,,,,,,'eg}~%" ymin ymax))))
 292         (let ((legend (getf options '$legend))
 293               (colors (getf options '$color))
 294               (styles (getf options '$style)) (i 0) j style plot-name)
 295           (unless (listp legend) (setq legend (list legend)))
 296           (unless (listp colors) (setq colors (list colors)))
 297           (unless (listp styles) (setq styles (list styles)))
 298           (loop for v in (cdr fun) for points-list in points-lists do
 299                (when points-list
 300                  ;; case "contour" with several plots in one list
 301                  (when ($listp (car points-list))
 302                    (setq j 0)
 303                    (dolist (level (cdar points-list))
 304                      (if styles
 305                          (setq style (nth (mod i (length styles)) styles))
 306                          (setq style nil))
 307                      (when ($listp style) (setq style (cdr style)))
 308                      (incf i)
 309                      (incf j)
 310                      (format st " {label ~s} " (ensure-string level))
 311                      (format st (xmaxima-curve-style style colors i))
 312                      (format st "~%{xversusy~%")
 313                      (let ((lis (cdr (nth j points-list))))
 314                        (loop while lis do
 315                             (loop while (and lis (not (eq (car lis) 'moveto)))
 316                                collecting (car lis) into xx
 317                                collecting (cadr lis) into yy
 318                                do (setq lis (cddr lis))
 319                                finally
 320                                ;; only output if at least two points for line
 321                                  (when (cdr xx)
 322                                    (tcl-output-list st xx)
 323                                    (tcl-output-list st yy)))
 324                           ;; remove the moveto
 325                             (setq lis (cddr lis))))
 326                      (format st "}"))
 327                    (return))
 328                  ;; other cases with only one plot per list
 329                  (if styles
 330                      (setq style (nth (mod i (length styles)) styles))
 331                      (setq style nil))
 332                  (when ($listp style) (setq style (cdr style)))
 333                  (incf i)
 334                  ;; label the expression according to the legend,
 335                  ;; unless it is "false" or there is only one expression
 336                  (if (member '$legend options)
 337                      (setq plot-name
 338                            (if (first legend)
 339                                (ensure-string
 340                                 (nth (mod (- i 1) (length legend)) legend)) nil))
 341                      (if (= 2 (length fun))
 342                          (setq plot-name nil)
 343                          (progn
 344                            (setq
 345                             plot-name
 346                             (with-output-to-string (pn)
 347                               (cond ((atom v) (format pn "~a" ($sconcat v)))
 348                                     ((eq (second v) '$parametric)
 349                                      (format pn "~a, ~a"
 350                                              ($sconcat (third v))
 351                                              ($sconcat (fourth v))))
 352                                     ((eq (second v) '$discrete)
 353                                      (format pn "discrete~a" i))
 354                                     (t (format pn "~a" ($sconcat v))))))
 355                            (when (> (length plot-name) 50)
 356                              (setq plot-name (format nil "fun~a" i))))))
 357                  (if plot-name
 358                      (format st " {label ~s} " plot-name)
 359                      (format st " {nolegend 1} "))
 360                  (format st (xmaxima-curve-style style colors i))
 361                  (format st "~%{xversusy~%")
 362                  (let ((lis points-list))
 363                    (loop while lis
 364                       do
 365                         (loop while (and lis (not (eq (car lis) 'moveto)))
 366                            collecting (car lis) into xx
 367                            collecting (cadr lis) into yy
 368                            do (setq lis (cddr lis))
 369                            finally
 370                            ;; only output if at least two points for line
 371                              (when (cdr xx)
 372                                (tcl-output-list st xx)
 373                                (tcl-output-list st yy)))
 374                       ;; remove the moveto
 375                         (setq lis (cddr lis))))
 376                  (format st "}")))
 377           (format st "} ")))))))
 378
 379 (defmethod plot3d-command ((plot xmaxima-plot) functions options titles)
 380   (let ((i 0) fun xrange yrange lvars trans)
 381     (setf
 382      (slot-value plot 'data)
 383      (concatenate
 384       'string
 385       (slot-value plot 'data)
 386       (with-output-to-string ($pstream)
 387         ;; generate the mesh points for each surface in the functions stack
 388         (dolist (f functions)
 389           (setq i (+ 1 i))
 390           (setq fun (first f))
 391           (setq xrange (second f))
 392           (setq yrange (third f))
 393           (if ($listp fun)
 394               (progn
 395                 (setq trans
 396                       ($make_transform `((mlist) ,(second xrange)
 397                                          ,(second yrange) $z)
 398                                        (second fun) (third fun) (fourth fun)))
 399                 (setq fun '$zero_fun))
 400               (let*
 401                   ((x0 (third xrange))
 402                    (x1 (fourth xrange))
 403                    (y0 (third yrange))
 404                    (y1 (fourth yrange))
 405                    (xmid (+ x0 (/ (- x1 x0) 2)))
 406                    (ymid (+ y0 (/ (- y1 y0) 2))))
 407                 (setq lvars `((mlist) ,(second xrange) ,(second yrange)))
 408                 (setq fun (coerce-float-fun fun lvars))
 409                 ;; Evaluate FUN at the middle point of the range.
 410                 ;; Looking at a single point is somewhat unreliable.
 411                 ;; Call FUN with numerical arguments (symbolic arguments may
 412                 ;; fail due to trouble computing real/imaginary parts for
 413                 ;; complicated expressions, or it may be a numerical function)
 414                 (when (cdr ($listofvars (mfuncall fun xmid ymid)))
 415                   (mtell
 416                    (intl:gettext
 417                     "plot3d: expected <expr. of v1 and v2>, [v1,min,max], [v2,min,max]~%"))
 418                   (mtell
 419                    (intl:gettext
 420                     "plot3d: keep going and hope for the best.~%")))))
 421           (let* ((pl
 422                   (draw3d
 423                    fun (third xrange) (fourth xrange) (third yrange)
 424                    (fourth yrange) (first (getf options '$grid))
 425                    (second (getf options '$grid))))
 426                  (ar (polygon-pts pl))
 427                  (colors (getf options '$color))
 428                  (palettes (getf options '$palette)))
 429             (declare (type (cl:array t) ar))
 430             (when trans (mfuncall trans ar))
 431             (when (getf options '$transform_xy)
 432                 (mfuncall (getf options '$transform_xy) ar))
 433             (if palettes
 434                    (format $pstream " ~a~%" (xmaxima-palettes palettes i))
 435                    (format $pstream " {mesh_lines ~a}" (xmaxima-color colors i)))
 436               (output-points-tcl $pstream pl (first (getf options '$grid)))))
 437         (format $pstream "}~%"))))))
 438
 439 (defmethod plot-shipout ((plot xmaxima-plot) options &optional output-file)
 440   (let ((file (plot-file-path (format nil "~a.xmaxima" (random-name 16)))))
 441     (cond ($show_openplot
 442            (with-open-file (fl
 443                             #+sbcl (sb-ext:native-namestring file)
 444                             #-sbcl file
 445                             :direction :output :if-exists :supersede)
 446              (princ (slot-value plot 'data) fl))
 447            ($system (concatenate 'string *maxima-prefix*
 448                                  (if (string= *autoconf-windows* "true")
 449                                      "\\bin\\" "/bin/")
 450                                  $xmaxima_plot_command)
 451                     #-(or (and sbcl win32) (and sbcl win64) (and ccl windows))
 452                     (format nil " ~s &" file)
 453                     #+(or (and sbcl win32) (and sbcl win64) (and ccl windows))
 454                     file))
 455           (t (princ (slot-value plot 'data)) ""))
 456     (cons '(mlist) (cons file output-file))))