Iterative extraction of xy coordinates of a given pixel value from image

Hello all,

I have TIFF images with 1 pixel thick lines on them which have unique value. I am trying to extract the xy coordinates of each value and fit a line to them and store the results of the fit. I need to do this iteratively, for each unique value and for hundreds of images. I have written a procedure that loads the TIFF and then calls another function (which is shown below).

My problem: it works for the first iteration, but then fails on subsequent iterations. I can't figure out what is wrong. The debugger shows on the second loop the conditional wave assignment results in all NaN waves. I can't figure out why. The value is present in the image. I've tried rounding to ensure that 2 is an integer and not 2.000001 or something.

Can anybody help? Please feel free to tell me if there is a simpler/faster way to extract the coordinates!

Function Extractor(m0)
    Wave m0
    ImageStats m0
    NVAR /Z nZ = fileIndex  // global variable
    Variable lastMT = V_max
    String wName
    DoWindow/K allPlot
    Display /N=allPlot
    Variable i
    for(i = 1; i < (lastMT + 1); i += 1) // MT, 1-based
        Duplicate/O m0, tempXw
        Duplicate/O m0, tempYw
        tempXw = (m0 == i) ? p : NaN
        tempYw = (m0 == i) ? q : NaN
        Redimension/N=(V_npnts) tempXw,tempYw
        WaveTransform zapnans tempXw
        WaveTransform zapnans tempYw
        Print i, ":", numpnts(tempXw) // for debugging
        if(numpnts(tempXw) > 2)
            CurveFit/Q/NTHR=0 line tempYw /X=tempXw /D
            WAVE /Z fit_tempYw // calling inside loop in case this is the problem
            wName = "vec_" + num2str(nZ) + "_" + num2str(i)
            Make/O/N=(2,2) $wName
            Wave m1 = $wName
            m1[0][0] = leftx(fit_tempYw)
            m1[1][0] = rightx(fit_tempYw)
            m1[0][1] = fit_tempYw[0]
            m1[1][1] = fit_tempYw[1]
            AppendToGraph/W=allPlot m1[][1] vs m1[][0] // for debugging
        KillWaves tempXw,tempYw,fit_tempYw // added this but it doesn't help
//  KillWaves m0,tempXw,tempYw,tempfitw

Edit: if I remove the if-endif code the first part of the loop (conditional asssignments) runs OK.
Edit2: if I comment out the CurveFit, function runs OK.

Problem solved! I put the curve fit in a separate function and this solved the problem. Picked up that tip here: no idea why it works.

I'll leave this here in case it helps anyone or if anyone has a better idea to pull the co-ordinates out. On a 768 * 768 image with 16 values it takes less than a second, which hopefully will be OK when I process all the images.
Presumably m0 is a 2D matrix. Wave assignments with more than a single dimension should explicitly use [] notation on the RHS, like:
        tempXw = (m0[p][q] == i) ? p : NaN
        tempYw = (m0[p][q] == i) ? q : NaN

But I would think that would cause problems with the first iteration as well.

John Weeks
WaveMetrics, Inc.