I think this gives the correct plot. It also uses viewmat, but just to plot the final array.
NB. Y is picture as an array with 1 at black points, x is values to test at NB. Returns r values of each theta value x. hough=:4 :0 ind=. ($ #: I.@,) y NB. Indices of points in y ind +/ .* 2 1 o./ x ) NB. Y is picture, x is resolution. Produces and sends to viewmat a density diagram NB. In hough space with resolution x drawhough=:4 :0 'w h'=.x val=. y hough~ +:o. (%~ 0.5+i.) w NB. Test r values for each pixel at each theta given by the centers of pixels 'min max'=. (<./@, , >./@,) val rel=.h * (max-min) %~ min -~ val NB. Rescale so that y values go from 0 to h pixels=.([: <:@(#/.~) (i.h) , ])"1&.|: <:@>. rel NB. Round y values; consolidate into picture viewmat pixels ) I have tested this on |. =/~ i. 100, but not on the actual picture. Note that it uses the top left corner as the natural coordinate center. Also, I made theta go from 0 to 2p1. Should these be the bounds, or just 0 to 1p1? Marshall -----Original Message----- From: [email protected] [mailto:[email protected]] On Behalf Of Aai Sent: Tuesday, September 14, 2010 10:16 AM To: Programming forum Subject: Re: [Jprogramming] Hough transform Hi David, I don't know exactly what you expect but here's an alternative using viewmat. require 'viewmat' Because you just increment coordinate values by 1 for every iteration it's sufficient to compose a triple X Y Z in which Z is the 'height' or 'density' of a point in hough space. hough2=: 4 : 0 radii=: (>:length>:$y)*_0.5+grid {. x angle=: 1r1p1*grid {: x b=: 0 < ,y inter=: b # , y coord=: ($y) #: I. b rhoug=: (coord,.1)+/ .*(1 0,0 _1,:_0.5 0.5*$y)+/ .*1 2 o./angle r=: (<:@# - I.&rhoug)radii 'X Y Z'=.|:({.,#)/.~,/r ,."1 i. {: x viewmat Z (;/ X,.Y) } 0 $~ 0 _1 { x ) try e.g. 250 hough2 data -- Met vriendelijke groet, =@@i ---------------------------------------------------------------------- For information about J forums see http://www.jsoftware.com/forums.htm ---------------------------------------------------------------------- For information about J forums see http://www.jsoftware.com/forums.htm
