class GSL::Contour2

Constants

VALID_CONNECTIONS

Edges: 1 __

0 | 4 \ /5  | 2
  | 7 / \ 6 |
    3 __

Public Class Methods

alloc(x, y, grid, function=nil) click to toggle source

Create a new Contour object. x and y are vectors of coordinates, and grid is a matrix of values on those coordinates. If a function is given, it will be used to evaluate gridpoints rather than the matrix (though a matrix must still be provided). Provide a function if your contours only cover a small part of the domain and your function is expensive to evaluate.

# File lib/gsl_extras.rb, line 560
def self.alloc(x, y, grid, function=nil)
        new(x, y, grid)
end
new(x, y, grid, function=nil) click to toggle source
# File lib/gsl_extras.rb, line 565
        def initialize(x, y, grid, function=nil)
                @function = function
                @evaluated = {}
                @x = x; @y=y; @grid=grid
#                       p @grid, @x, @y
                raise ArgumentError.new("Unmatching data sizes: #{x.size}, #{y.size}, #{grid.shape}") unless [x.size, y.size] == grid.shape
                @adaptive = false                     
        end

Public Instance Methods

contours(*values) click to toggle source

Create a series of contours at the given values. Returns a hash of {value => array_of_contours}. The array_of_contours is a an array of arrays, where each array is a list of [x, y] coordinates along the contour.

# File lib/gsl_extras.rb, line 582
def contours(*values)
        (values = (0..values[0]+1).to_a.map{|i| i.to_f * (@grid.max - @grid.min) / ( values[0]+1) + @grid.min}; values.pop; values.shift) if values.size==1 and values[0].kind_of? Integer
        cons = values.inject({}){|hash, val| hash[val] = []; hash}
        get_startpoints(cons)
        #@analysed = {}
        @found = {}
        #p cons
        starts = cons

        cons = values.inject({}){|hash, val| hash[val] = []; hash}
        temp_cons = values.inject({}){|hash, val| hash[val] = []; hash}
        starts.each do |val, arr|
                #p 'arr', arr
                arr.each do |start_con|
                #arr.map{|edges| edges[-1][1].slice(0..1)}.each do |starti, startj|
                        starti, startj = start_con[-1][1].slice(0..1)       
                        temp_cons[val] = [start_con]
                        p 'startj',  starti, startj
                        loop do
                                starti, startj = trace_contour(val, temp_cons, starti, startj)
                                break unless starti and startj
                        end
                        cons[val].push temp_cons[val][0]
                end
        end
        cons.keys.each{|val| cons[val] = connect_contours(cons[val]);
                cons[val].map!{|con| con.map{|point| point[0]}}}
        
        @last_contours = cons
        gk =  graphkit
        #gk.gp.style = 'data with linespoints'
        gk.data.each{|dk| dk.with = 'lp'}
        gk.gnuplot
        @last_contours
end
get_startpoints(cons) click to toggle source
# File lib/gsl_extras.rb, line 657
def get_startpoints(cons)      
        [0,((@x.size / 2).floor - 1)-1].each do |i|
                for j in 0...((@y.size / 2).floor - 1)
                        analyse_cell(i*2, j*2, cons)
                end
        end
        for i in 0...((@x.size / 2).floor - 1)
        [0,((@y.size / 2).floor - 1)-1].each do |j|
                        analyse_cell(i*2, j*2, cons)
                end
        end
end
graphkit(*args) click to toggle source

Create a GraphKit object of the contours.

# File lib/gsl_extras.rb, line 674
def graphkit(*args)
        if args.size == 0
                conts = @last_contours
        else
                conts = contours(*args)
        end
        graphs = conts.map do |val, cons|
                unless cons[0]
                        nil
                else
                        (cons.map do |con|
#                              p con
                                contour = con.transpose
                                kit = CodeRunner::GraphKit.autocreate({x: {data: contour[0]}, y: {data: contour[1], title: val.to_s}})
                                kit.data[0].with = "l"
                                kit
                        end).sum
                end
        end
        graphs.compact.reverse.sum
end
trace_contour(val, cons, starti, startj) click to toggle source
# File lib/gsl_extras.rb, line 618
def trace_contour(val, cons, starti, startj)
                        old_start = cons[val][0][0][1]; old_end = cons[val][0][-1][1]
                        #p 'old_start', old_start, 'old_end', old_end, 'starti', starti, 'startj', startj
                        for delti in -1..1
                                for deltj in -1..1
                                        celli, cellj = [starti + 2 * delti, startj + deltj * 2]
                                                #unless (@analysed[val] and @analysed[val][[celli, cellj]] ) or celli > (@x.size-3) or cellj >(@y.size - 3) or celli < 0 or cellj < 0
                                                unless  celli > (@x.size-3) or cellj >(@y.size - 3) or celli < 0 or cellj < 0

                                           #p 'analysing', celli, cellj
                                           analyse_cell(celli, cellj, cons, val) 
                                                end

                                end
                        end
                                cons[val] = connect_contours(cons[val])
                        #p cons[val]
                                new_contour =  cons[val].find{|cont| old_start ==  cont[0][1] or old_end == cont[-1][1]}
                                unless new_contour
                                        cons[val].map!{|con| con.reverse}
                                new_contour =  cons[val].find{|cont| old_start ==  cont[0][1] or old_end == cont[-1][1]}
                                end
                                raise "no new contour" unless new_contour
                                cons[val] = [new_contour]
                        #p cons[val]
                        #newtails  = cons[val].map{|con| con[-1]}
                        if cons[val][0][0][1] != old_start 
                                return cons[val][0][0][1].slice(0..1)
                        elsif cons[val][0][-1][1] != old_end
                                return cons[val][0][-1][1].slice(0..1)
                        else
                                #p cons[val]
                                
                        p 'old_start', old_start, 'old_end', old_end
                                #raise "Could not connect"
                                return nil, nil
                        end
end

Private Instance Methods

analyse_cell(i, j, cons, specific_val = nil) click to toggle source
# File lib/gsl_extras.rb, line 742
        def analyse_cell(i, j, cons, specific_val = nil)
                crossed_edges = get_crossed_edges(i, j, cons, specific_val)
                crossed_edges.each do |val, crossings|
                        if specific_val
                         next        unless val == specific_val
                         #@analysed[val] ||= {}
                         #@analysed[val][[i,j]] = true
                         @found[[i,j]] ||= {}
                        end
                        outer = crossings.keys.find_all{|edge| edge < 4}
                        inner = crossings.keys.find_all{|edge| edge > 4}
                        next if outer.size == 0 and inner.size == 0
                        VALID_CONNECTIONS.each do |connection|
                                path = crossings.values_at(*connection).compact.zip(connection.map{|edge| [i, j, edge]})
#                                       p path
                                next if path.size != connection.size
                                connection.each{|edge| crossings.delete(edge)}
                                if specific_val
                                 next if @found[[i,j]][path]
                                 @found[[i,j]][path] = true

                                end
                                cons[val].push path
                                
                        end

#                               p val, crossings, cons
                end
        end
connect_contours(contours) click to toggle source
# File lib/gsl_extras.rb, line 774
        def connect_contours(contours)
                return contours if contours.size == 1
                loop do
                        catch(:restart) do 
#                               old_contours = contours
#                               contours = []
                        joined = []
                        for i in 0...contours.size
                                break if i >= contours.size
                                for j in i+1...contours.size
                                        break if j >= contours.size
                                        coni = contours[i]
                                        conj = contours[j]
                                        if joins?(coni[-1], conj[0])
                                                contours[i] = coni + conj
                                                contours.delete_at(j)
#                                                       redo
                                                throw(:restart)
                                        elsif joins?(coni[0], conj[-1])
                                                contours[i] = conj + coni
                                                contours.delete_at(j)
#                                                       redo
                                                throw(:restart)
                                        elsif joins?(coni[0], conj[0])
                                                contours[i] = coni.reverse + conj
                                                contours.delete_at(j)
#                                                       redo
                                                throw(:restart)
                                        elsif joins?(coni[-1], conj[-1])
                                                contours[i] = conj + coni.reverse
                                                contours.delete_at(j)
#                                                       redo
                                                throw(:restart)
                                        end


                                        
                                end
                        end
                        contours.each{|con| con.uniq!}
                        return contours
                        end
                        
                end
        end
get_crossed_edges(i, j, cons, specific_val = nil) click to toggle source
# File lib/gsl_extras.rb, line 704
        def get_crossed_edges(i, j, cons, specific_val = nil)
                ce = {}
                edges = {0=>[[i, j], [i+2, j]],
                                1=>[[i, j], [i, j+2]],
                                2=>[[i, j+2], [i+2, j+2]],
                                3=>[[i+2, j], [i+2, j+2]],
                                4=>[[i, j], [i+1, j+1]],
                                5=>[[i, j+2], [i+1, j+1]],
                                6=>[[i+2, j+2], [i+1, j+1]],
                                7=>[[i+2, j], [i+1, j+1]]
                        }    
                cons.keys.each do |val|
                        next if specific_val and specific_val != val
                        ce[val] = {}
                        edges.each do |edge, (start, fin)|
                                if @function
                                        #p start, fin
                                        (@evaluated[start]= true, (@grid[*start] = @function.eval(@x[start[0]], @y[start[1]]))) unless @evaluated[start]
                                        (@evaluated[fin]= true, (@grid[*fin] = @function.eval(@x[fin[0]], @y[fin[1]])))  unless @evaluated[fin]
                                end
                                bounds = [@grid[*start], @grid[*fin]]
#                                       p edge, bounds if bounds.max < 4
                                if val <= bounds.max and val > bounds.min
                                        dx = @x[fin[0]] - @x[start[0]]
                                        dy = @y[fin[1]] - @y[start[1]]
                                        df = bounds[1] - bounds[0]
                                        xcross = @x[start[0]] + (val-bounds[0]) / df * dx
                                        ycross = @y[start[1]] + (val-bounds[0]) / df * dy
                                        ce[val][edge] = [xcross, ycross]
                                end
                                
                        end
                end
                ce
        end
joins?(path, contour) click to toggle source
# File lib/gsl_extras.rb, line 824
        def joins?(path, contour)      
#                       p @x.size
#                       dx = (@x.subvector(1, @x.size - 1) - @x.subvector(0, @x.size - 1)).sum / @x.size #The average dx between gridpoints
#                       dy = (@y.subvector(1, @y.size - 1) - @y.subvector(0, @y.size - 1)).sum / @y.size
#                       eps = 1.0
# #                     return true
#                       return ((path[0] - contour[0]).abs < dx.abs * eps and (path[1] - contour[1]).abs < dy.abs * eps)
#                       p dx; exit
                pi, pj, pedge = path[1]
                ci, cj, cedge = contour[1]
                joins = false
                joins = ((pi == ci and pj == cj + 2 and pedge == 0 and cedge == 2) or
                        (pi == ci and pj == cj - 2 and pedge == 2 and cedge == 0) or
                        (pi == ci + 2 and pj == cj  and pedge == 1 and cedge == 3) or
                        (pi == ci - 2 and pj == cj  and pedge == 3 and cedge == 1))
#                       unless joins
#                               p path[1], contour[1]
#                               STDIN.gets
#                       end
#                       if path[1] == [2,10,1]
#                               p path[1], contour[1], joins
#                               STDIN.gets
#                       end
                return joins
        end