Group Abstract Group Abstract

Message Boards Message Boards

0
|
2.3K Views
|
6 Replies
|
3 Total Likes
View groups...
Share
Share this post:

Fixing old code that doesn't work on V13.3

Posted 1 year ago

I have been given a file from an older mathmatica version and my task is to get it working in mathmatica 13.3. While using the //Echo, I noticed that the code below is not computing the values for the variables, but just passing off a place holder and giving error messages. I am struggling to understand why that is and how to fix it. This module should function as the main algorithm for a heap sorting and fast marching on a rasterized image. I would appreciate any assistance or advice

findDistances {bandheap,len=Length[ll],wid=Length[ll[[1]]],hsize,hindex=0,distancetable,statetable,frozen=-1.,far=300.,band=0.,outofbounds=1000.,dist,maxpixval=255.,minpixval=0.,midpixval,heaptable,j1,j2,nbrs,pt,x,y,x1,y1,x2,y2,next,prev,done}, (Make sure everything in here is intuallzed and called correctlly) (Initializations.)

midpixval=(maxpixval-minpixval); distancetable =Map[If[#>midpixval,maxpixval,minpixval]&,ll,{2}]; statetable=Map[If[#1>1,far,frozen]&,distancetable,{2}]; heaptable=Table[0,{len},{wid}]; hsize=len*wid; bandheap=Table[{0.,0.,0.},{hsize}];

(Initialize the band by placing all neighbors of frozen cells thereon.) Do[ If[statetable[[ii,jj]]=!=frozen,Continue[]]; nbrs={{ii,jj+1},{ii,jj-1},{ii-1,jj},{ii+1,jj}};

Do[ {x,y}=nbrs[[kk]]; If[! (0<x<=len&&0<y<=wid&&Echo[statetable[[x,y]]]==far),Continue[]]; hindex++; statetable[[x,y]]=band; dist=newDistance[x,y,statetable,distancetable]; distancetable[[x,y]]=dist; bandheap[[hindex]]={dist,N[x],N[y]}; j1=hindex;

(Place each neighbor to the frozen cell on the end of the band heap. Then move it upward to the appropriate location.) While[ (j2=Floor[j1/2])>=1&&bandheap[[j2,1]]>bandheap[[j1,1]], bandheap[[{j1,j2}]]=bandheap[[{j2,j1}]]; {x1,y1}=Round[Rest[bandheap[[j1]]]]; heaptable[[x1,y1]]=j1; j1=j2; ]; heaptable[[x,y]]=j1, {kk,Length[nbrs]} ], {ii,len},{jj,wid} ];

While[hindex>0, pt=bandheap[[1]]; (Get smallest element in heap and freeze it.) {x,y}=Round[Rest[pt]]; statetable[[x,y]]=frozen;

(Replace with last element in heap,then rearrange by moving it downward to restore the heap ordering property.) bandheap[[1]]=bandheap[[hindex]]; done=False; prev=1; {j1,j2}=2*prev+{0,1};

While[j1<hindex&&!done, If[j2<hindex, next=If[bandheap[[j1,1]]<=bandheap[[j2,1]],j1,j2], next=j1 ];

If[bandheap[[prev,1]]>bandheap[[next,1]],(Swap? If not, we are done.) bandheap[[{prev,next}]]=bandheap[[{next,prev}]]; {x1,y1}=Round[Rest[bandheap[[prev]]]]; heaptable[[x1,y1]]=prev; prev=next; {j1,j2}=2*prev+{0,1} ,(else) done=True,done=True ]; ]; {x1,y1}=Round[Rest[bandheap[[prev]]]]; heaptable[[x1,y1]]=prev;

(Find neighbors of the removed point that were already in the band, update their distances, and move them to new locations in the band heap as needed.) nbrs={{x,y+1},{x,y-1},{x-1,y},{x+1,y}}; Do[ {x2,y2}=nbrs[[kk]]; If[! (0<x2<=len&&0<y2<=wid&&statetable[[x2,y2]]==band),Continue[]]; dist=newDistance[x2,y2,statetable,distancetable]; distancetable[[x2,y2]]=dist; (Find position of neighbor in heap. Change its value to reflect new distance.) j1=heaptable[[x2,y2]]; bandheap[[j1]]={dist,N[x2],N[y2]};

(The new distance is no larger than the old, so it can only move upward (that is, we iteratively compare to parent and swap if it is smaller).) While[(j2=Floor[j1/2])>=1&&bandheap[[j2,1]]>bandheap[[j1,1]], bandheap[[{j1,j2}]]=bandheap[[{j2,j1}]]; {x1,y1}=Round[Rest[bandheap[[j1]]]]; heaptable[[x1,y1]]=j1; j1=j2; ]; heaptable[[x2,y2]]=j1 , {kk,Length[nbrs]} ];

hindex--;

(Find neighbors of the removed point that are not already in the band, and place them there.) Do[ {x2,y2}=nbrs[[kk]]; If[! (0<x2<=len&&0<y2<=wid&&statetable[[x2,y2]]==far),Continue[]]; hindex++; statetable[[x2,y2]]=band; dist=newDistance[x2,y2,statetable,distancetable]; distancetable[[x2,y2]]=dist; bandheap[[hindex]]={dist,N[x2],N[y2]}; j1=hindex;

(Place each neighbor on the band heap and move upward to the appropriate location.) While[ (j2=Floor[j1/2])>=1&&bandheap[[j2,1]]>bandheap[[j1,1]], bandheap[[{j1,j2}]]=bandheap[[{j2,j1}]]; {x1,y1}=Round[Rest[bandheap[[j1]]]]; heaptable[[x1,y1]]=j1; j1=j2; ]; heaptable[[x2,y2]]=j1, {kk,Length[nbrs]} ];

]; distancetable ]

POSTED BY: Aidan Vieson
6 Replies

I would recommend ditching the circa 2003 code and using one of the newer variants such as this one from Mathematica StackExchange. For added speed I would (if I were writing this today) likely use the built-in "PriorityQueue" data structure.

As for the fast marching part, to be honest I doubt the author remembers much about it. He might have a closer look if absolutely necessary. What is your intended use case? (I hope something more exciting than making contours around initials. Modeling some physical phenomenon perhaps.)

POSTED BY: Daniel Lichtblau
Posted 1 year ago

This is the file I have been working from. The Module I am having issues with is close to the bottom labeled "Main Algorithm" .

POSTED BY: Aidan Vieson
Posted 1 year ago
POSTED BY: Bill Nelson
Posted 1 year ago

The original document was made in 2003 but I cant seem to find any documents from then

POSTED BY: Aidan Vieson
Posted 1 year ago

I've been using the //Echo an the only this I see going wrong is with the module labeled Main algorithm. It seems to not be calculating the values for distance tables nor state tables. I am trying to input the a rasterized image, but I've also tried ImageData["image name"] and gotten the same results. is there a way to tell the program to solve here and not just carry place-holders?

POSTED BY: Aidan Vieson
Posted 1 year ago

If I just scrape-n-paste the code you show and I fix a few problems then it appears that you are missing a few characters or lines at the top of this because there is an unmatched ] at the bottom of this. There doesn't seem to be enough to scrape the code below, paste it into a fresh empty notebook, tap <shift><enter> to run this and see exactly the behavior that you are describing..

findDistances
    {bandheap, len=Length[ll], wid=Length[ll[[1]]], hsize, hindex=0,
    distancetable, statetable, frozen=-1., far=300., band=0.,
    outofbounds=1000., dist, maxpixval=255., minpixval=0., midpixval, heaptable,
    j1, j2, nbrs, pt, x, y, x1, y1, x2, y2, next, prev, done},
    (*Make sure everything in here is intuallzed and called correctlly*)
    (*Initializations.*)
midpixval=(maxpixval-minpixval);
distancetable =Map[If[#>midpixval,maxpixval,minpixval]&,ll,{2}];
statetable=Map[If[#1>1,far,frozen]&,distancetable,{2}];
heaptable=Table[0,{len},{wid}];
hsize=len*wid;
bandheap=Table[{0.,0.,0.},{hsize}];
(*Initialize the band by placing all neighbors of frozen cells thereon.*)
Do[
    If[statetable[[ii,jj]]=!=frozen,
        Continue[]
    ];
    nbrs={{ii,jj+1},{ii,jj-1},{ii-1,jj},{ii+1,jj}};
    Do[
        {x,y}=nbrs[[kk]];
        If[! (0<x<=len&&0<y<=wid&&Echo[statetable[[x,y]]]==far),
            Continue[]
        ];
        hindex++;
        statetable[[x,y]]=band;
        dist=newDistance[x,y,statetable,distancetable];
        distancetable[[x,y]]=dist;
        bandheap[[hindex]]={dist,N[x],N[y]};
        j1=hindex;
        (*Place each neighbor to the frozen cell on the end of the band heap. Then move it upward to the appropriate location.*)
        While[ (j2=Floor[j1/2])>=1&&bandheap[[j2,1]]>bandheap[[j1,1]],
            bandheap[[{j1,j2}]]=bandheap[[{j2,j1}]];
            {x1,y1}=Round[Rest[bandheap[[j1]]]];
            heaptable[[x1,y1]]=j1;
            j1=j2;
        ];
        heaptable[[x,y]]=j1,
        {kk,Length[nbrs]}
    ],
    {ii,len},{jj,wid}
];
While[hindex>0,
    pt=bandheap[[1]];
    (*Get smallest element in heap and freeze it.*)
    {x,y}=Round[Rest[pt]];
    statetable[[x,y]]=frozen;
    (*Replace with last element in heap,then rearrange by moving it downward to restore the heap ordering property.*)
    bandheap[[1]]=bandheap[[hindex]];
    done=False;
    prev=1;
    {j1,j2}=2*prev+{0,1};
    While[j1<hindex&&!done,
        If[j2<hindex,
            next=If[bandheap[[j1,1]]<=bandheap[[j2,1]],
                j1,
                j2
            ],
            next=j1
        ];
        If[bandheap[[prev,1]]>bandheap[[next,1]],
            (*Swap? If not, we are done.*)
            bandheap[[{prev,next}]]=bandheap[[{next,prev}]];
            {x1,y1}=Round[Rest[bandheap[[prev]]]];
            heaptable[[x1,y1]]=prev;
            prev=next;
            {j1,j2}=2*prev+{0,1} ,
            (*else*)
            done=True,
            done=True
        ];
    ];
    {x1,y1}=Round[Rest[bandheap[[prev]]]];
    heaptable[[x1,y1]]=prev;
    (*Find neighbors of the removed point that were already in the band, update their distances, and move them to new locations in the band heap as needed.*)
    nbrs={{x,y+1},{x,y-1},{x-1,y},{x+1,y}};
    Do[
        {x2,y2}=nbrs[[kk]];
        If[! (0<x2<=len&&0<y2<=wid&&statetable[[x2,y2]]==band),
            Continue[]
        ];
        dist=newDistance[x2,y2,statetable,distancetable];
        distancetable[[x2,y2]]=dist;
        (*Find position of neighbor in heap. Change its value to reflect new distance.*)
        j1=heaptable[[x2,y2]];
        bandheap[[j1]]={dist,N[x2],N[y2]};
        (*The new distance is no larger than the old, so it can only move upward (that is, we iteratively compare to parent and swap if it is smaller).*)
        While[(j2=Floor[j1/2])>=1&&bandheap[[j2,1]]>bandheap[[j1,1]],
            bandheap[[{j1,j2}]]=bandheap[[{j2,j1}]];
            {x1,y1}=Round[Rest[bandheap[[j1]]]];
            heaptable[[x1,y1]]=j1;
            j1=j2;
        ];
        heaptable[[x2,y2]]=j1 ,
        {kk,Length[nbrs]}
    ];
    hindex--;
    (*Find neighbors of the removed point that are not already in the band, and place them there.*)
    Do[
        {x2,y2}=nbrs[[kk]];
        If[! (0<x2<=len&&0<y2<=wid&&statetable[[x2,y2]]==far),
            Continue[]
        ];
        hindex++;
        statetable[[x2,y2]]=band;
        dist=newDistance[x2,y2,statetable,distancetable];
        distancetable[[x2,y2]]=dist;
        bandheap[[hindex]]={dist,N[x2],N[y2]};
        j1=hindex;
        (*Place each neighbor on the band heap and move upward to the appropriate location.*)
        While[ (j2=Floor[j1/2])>=1&&bandheap[[j2,1]]>bandheap[[j1,1]],
            bandheap[[{j1,j2}]]=bandheap[[{j2,j1}]];
            {x1,y1}=Round[Rest[bandheap[[j1]]]];
            heaptable[[x1,y1]]=j1;
            j1=j2;
        ];
        heaptable[[x2,y2]]=j1,
        {kk,Length[nbrs]}
    ];
];
distancetable ]
POSTED BY: Bill Nelson
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard