Changeset 27214 in project


Ignore:
Timestamp:
08/06/12 12:16:07 (9 years ago)
Author:
Ivan Raikov
Message:

neurolucida: beginnings of Neurolucida graph format transformation

Location:
release/4/neurolucida/trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • release/4/neurolucida/trunk/neurolucida.meta

    r27172 r27214  
    1818 ; A list of eggs neurolucida depends on.
    1919
    20  (needs sxml-transforms sxpath ssax getopt-long digraph format-graph)
     20 (needs sxml-transforms sxpath ssax getopt-long (digraph 1.16) graph-dfs format-graph)
    2121
    2222 (author "Ivan Raikov")
  • release/4/neurolucida/trunk/neurolucida.scm

    r27172 r27214  
    2121
    2222(require-extension files data-structures srfi-1 srfi-13  getopt-long sxml-transforms sxpath sxpath-lolevel ssax)
    23 (require-extension digraph format-graph  )
     23(require-extension digraph graph-dfs format-graph  )
    2424
    2525
     
    4040        (else str)))
    4141
     42
     43; Returns attr-list node for a given obj
     44;   or #f if it is absent
     45(define (sxml:attr-alist obj)
     46  (if (and (not (null? (cdr obj)))
     47            (pair? (cadr obj))
     48            (eq? '@ (caadr obj)))
     49         (cdr (cadr obj))
     50         #f))
     51
     52;; obtain  child named n of a node
     53(define (sxml:kidn name node)
     54  ((select-first-kid (lambda (x)  (eq? (car x) name))) node))
     55
     56;; obtain all children of a node named n
     57(define (sxml:kidsn name node)
     58  ((select-kids (lambda (x) (eq? (car x) name))) node))
    4259
    4360(define (mkdir dir)
     
    8299
    83100    (format
    84      "output format (swc, asc, vcg)"
     101     "output format (swc, asc, vcg, cyl)"
    85102     (single-char #\f)
    86103     (value       (required "FORMAT")
     
    163180  (let* ((sxml       `(*TOP* ,@(cdr x)))
    164181         (g          (make-digraph label #f))
    165          (add-node!  (g 'add-node!))
    166          (add-edge!  (g 'add-edge!)))
     182         (node-info      (g 'node-info))
     183         (node-info-set! (g 'node-info-set!))
     184         (add-node!      (g 'add-node!))
     185         (add-edge!      (g 'add-edge!))
     186         )
    167187
    168188    (add-node! 0 'soma)
     
    176196                                       ((sxpath `(nl:branch)) tree))))
    177197
    178         (let ((n (length points+spines)))
     198        ;; insert points in the dependency graph
     199        (let ((node-id1
     200               (fold (lambda (n i)
     201                       (case (car n)
     202                         ((nl:point)
     203                          (begin
     204                            (add-node! i n) (+ 1 i)))
     205                         ((nl:spine)
     206                          (begin
     207                            (node-info-set! (- i 1) (append (node-info (- i 1)) (list n))) i))
     208                         (else (error 'make-tree-graph "unknown node type" n))))
     209                     node-id points+spines)))
    179210         
    180           ;; insert points and spines in the dependency graph
    181           (let ((node-id1 (fold (lambda (n i) (add-node! i n) (+ 1 i)) node-id points+spines)))
    182            
    183             ;; insert edges
    184             (if (positive? n)
    185                 (let ((fin (- node-id1 1)))
    186                   (add-edge! (list parent node-id #f))
    187                   (let inner-recur ((i node-id))
    188                     (if (< i fin)
    189                         (let ((j (+ 1 i)))
    190                           (add-edge! (list i j  #f))
    191                           (inner-recur j))
    192                         ))
     211          ;; insert edges
     212          (if (pair? points+spines)
     213              (let ((fin (- node-id1 1)))
     214                (add-edge! (list parent node-id #f))
     215                (let inner-recur ((i node-id))
     216                  (if (< i fin)
     217                      (let ((j (+ 1 i)))
     218                        (add-edge! (list i j  #f))
     219                        (inner-recur j))
     220                      ))
     221                ))
     222         
     223          (if (pair? subtrees)
     224              (let ((branch-node (- node-id1 1)))
     225                (node-info-set! branch-node `(nl:branch ,(node-info branch-node)))
     226                (fold (lambda (subtree i)
     227                        (recur subtree branch-node (+ i 1))) node-id1 subtrees))
     228              node-id1)
     229          )))
     230    g))
     231
     232
     233
     234;;
     235;; NG (Neurolucida Graph) format
     236;;
     237;; ((Branch ... (Node ...) (Branch ...) ) ... )
     238;;
     239;; Where ... is one or more attributes:
     240;;
     241;;   origin x y z
     242;;   radius r
     243;;   length l
     244;;   spine-density-linear
     245;;   branch-order-Soma
     246;;   branch-order-Terminal
     247;;
     248
     249(define (make-ng g)
     250
     251  (define (sum lst) (fold + 0. lst))
     252
     253
     254  (define (pdist2 a b)
     255    (let ((diff2 (lambda (i) (let ((v (- (list-ref a i) (list-ref b i)))) (* v v)))))
     256      (sum (list-tabulate 3 diff2))))
     257
     258
     259  (define (get-node-origin x)
     260    (let ((x (if (equal? 'nl:point (car x)) x (sxml:kidn 'nl:point x))))
     261      (let ((attrs (sxml:attr-alist x)))
     262        (list (string->number (car (lookup-def 'x attrs)))
     263              (string->number (car (lookup-def 'y attrs)))
     264              (string->number (car (lookup-def 'z attrs))))
     265        )))
     266
     267
     268  (define (get-node-radius x)
     269    (let ((x (if (equal? 'nl:point (car x)) x (sxml:kidn 'nl:point x))))
     270      (let ((attrs (sxml:attr-alist x)))
     271        (/ (string->number (car (lookup-def 'd attrs))) 2.))))
     272
     273 
     274   
     275  ;;
     276  ;; Compute node order from the soma and from the terminals
     277  ;;
     278 
     279  (define (compute-node-orders g roots)
     280   
     281    (let ((ginv  (make-digraph (g 'name) (g 'graph-info))))
     282     
     283      ;; create an inverse of the graph
     284      (let ((add-node! (ginv 'add-node!))
     285            (add-edge! (ginv 'add-edge!)))
     286        ((g 'foreach-node) (lambda (i x) (add-node! i x)))
     287        ((g 'foreach-edge) (lambda (i ee)
     288                             (for-each
     289                              (lambda (e)
     290                                (add-edge! (list (cadr e) (car e) (caddr e))))
     291                              ee)))
     292        )
     293     
     294      (let (;; branch order from soma
     295            (node-depth (graph-dfs-depth g roots))
     296            ;; branch order from terminals
     297            (inverse-node-depth (graph-dfs-depth ginv ((ginv 'roots))))
     298            )
     299       
     300        `((node-depth . ,node-depth)
     301          (inverse-node-depth . ,inverse-node-depth)
     302          )
     303        )))
     304
     305 
     306  (let* ((roots ((g 'roots)))
     307         (out-edges (g 'out-edges))
     308         (node-info (g 'node-info))
     309         (node-orders (compute-node-orders g roots))
     310         (node-depths (lookup-def 'node-depth node-orders))
     311         (inverse-node-depths (lookup-def 'inverse-node-depth node-orders))
     312         )
     313
     314    (let recur ((node (car roots)) (ax '()))
     315
     316       (let* ((info   (node-info node))
     317              (spines (sxml:kidsn 'nl:spine info))
     318              )
     319
     320         (print "info = " info)
     321
     322         (if (equal? info 'soma)
     323
     324             (fold recur '() (map cadr (out-edges node)))
     325
     326             (case (car info)
     327               
     328               ((nl:point)
     329                (let* ((oes (out-edges node))
     330                       (next-node (and (pair? oes) (cadr (car oes))))
     331                       (dd (print "next-node = " next-node))
     332                       (dd (print "info next-node = " (and next-node (node-info next-node))))
     333                       (cylinder-origin (get-node-origin info))
     334                       (cylinder-radius (get-node-radius info))
     335                       (cylinder-length (and next-node
     336                                             (let ((node1-origin (get-node-origin (node-info next-node))))
     337                                                 (sqrt (pdist2 cylinder-origin node1-origin))
     338                                                 ))
     339                                             )
     340                       )
     341
     342                  (print "cylinder-length = " cylinder-length)
     343
     344                  (assert (<= (length oes) 1))
     345                  (and cylinder-length (assert (positive? cylinder-length)))
     346                 
     347                  (let ((ax1 (cons
     348                              (if next-node
     349                                  `(Node (origin . ,cylinder-origin)
     350                                         (radius . ,cylinder-radius)
     351                                         (length . ,cylinder-length)
     352                                         (branch-order-Soma . ,(s32vector-ref node-depths node))
     353                                         (branch-order-Terminal . ,(s32vector-ref inverse-node-depths node))
     354                                         ,@(if (null? spines) '() `((spine-density-linear . ,(/ (length spines) cylinder-length))))
     355                                         )
     356                                  `(Terminal (origin . ,cylinder-origin)
     357                                             (radius . ,cylinder-radius)
     358                                             (branch-order-Soma . ,(s32vector-ref node-depths node))
     359                                             (branch-order-Terminal . ,(s32vector-ref inverse-node-depths node))
     360                                   ))
     361                              ax)))
     362
     363                    (if next-node (recur next-node ax1) ax1)
     364
     365                    )
     366                ))
     367               
     368               ((nl:branch)
     369                (let* ((branch-origin (get-node-origin info))
     370                       (branch-radius (get-node-radius info))
     371                       )
     372                 
     373                  (cons
     374                   `(Branch (origin . ,branch-origin)
     375                            (radius . ,branch-radius)
     376                            (branch-order-Soma . ,(s32vector-ref node-depths node))
     377                            (branch-order-Terminal . ,(s32vector-ref inverse-node-depths node))
     378                            (children . ,(fold recur '() (map cadr (out-edges node))))
     379                            )
     380                   ax)
    193381                  ))
    194            
    195             (if (pair? subtrees)
    196                 (let ((branch-node (- node-id1 1))) ;; make the last node to be branching node
    197                   ((g 'node-info-set!) branch-node (cons 'branch ((g 'node-info) branch-node)))
    198                   (fold (lambda (subtree i) (recur subtree branch-node i)) node-id1 subtrees))
    199                 node-id1
    200                 ))
    201           ))
    202       )
    203     g))
     382             
     383             (else (recur (cdr nodes) ax)))
     384         
     385             ))
     386       ))
     387  )
     388       
     389     
     390 
    204391
    205392
     
    400587                     (map (lambda (x) (make-pathname ddir (car (string-split (cadar x) "#")) "vcg")) data-map)
    401588                     data-map))
     589
     590          ((ng)
     591           (for-each (lambda (f x)
     592                       (with-output-to-file f
     593                         (lambda ()
     594                           (pp  (make-ng (make-tree-graph 'neurolucida x)) (current-output-port)))))
     595                     (map (lambda (x) (make-pathname ddir (car (string-split (cadar x) "#")) "ng")) data-map)
     596                     data-map))
     597
    402598          ))
    403599      )))
Note: See TracChangeset for help on using the changeset viewer.