Examples with Voxel3D mcib3d.geom.Voxel3D used on opensource projects

Example 21 with Voxel3D

use of mcib3d.geom.Voxel3D in project mcib3d-core by mcib3d.

the class ManualSpot method segmentSpots.

// End of variables declaration//GEN-END:variables
private void segmentSpots() {
    Roi roi;
    RoiManager manager = RoiManager.getInstance();
    Roi[] rois = manager.getRoisAsArray();
    if ((manager == null) || (rois.length == 0)) {
        IJ.showMessage("Put rois in the RoiManager");
        return;
    }
    if (objs == null) {
        objs = new Object3D[rois.length];
    }
    String text = "";
    for (int r = 0; r < rois.length; r++) {
        roi = rois[r];
        Rectangle rect = roi.getBounds();
        int x = rect.x;
        int y = rect.y;
        int z = manager.getSliceNumber(manager.getName(r));
        text = text.concat("\nx=" + x + " " + y + " " + z);
        jTextFieldLog.setText(text);
        Segment3DSpots segSpots = new Segment3DSpots(spot3DImage, null);
        segSpots.setLabelImage(labelImage);
        segSpots.show = true;
        // find local maximum
        Voxel3D v = segSpots.getLocalMaximum(x, y, z, 2, 2, 2);
        if (v != null) {
            text = text.concat("\nMax local " + v.getX() + " " + v.getY() + " " + v.getZ() + " " + v.getValue());
            jTextFieldLog.setText(text);
            x = v.getRoundX();
            y = v.getRoundY();
            z = v.getRoundZ();
        }
        // set parameter
        if (threshold_method == 0) {
            segSpots.setMethodLocal(Segment3DSpots.LOCAL_CONSTANT);
        } else if (threshold_method == 1) {
            segSpots.setMethodLocal(Segment3DSpots.LOCAL_MEAN);
        } else {
            segSpots.setMethodLocal(Segment3DSpots.LOCAL_GAUSS);
        }
        if (spot_method == 0) {
            segSpots.setMethodSeg(Segment3DSpots.SEG_CLASSICAL);
        } else if (spot_method == 0) {
            segSpots.setMethodSeg(Segment3DSpots.SEG_MAX);
        } else {
            segSpots.setMethodSeg(Segment3DSpots.SEG_BLOCK);
        }
        segSpots.setLocalThreshold(threshold_constant);
        segSpots.setGaussMaxr(rad_gauss);
        segSpots.setGaussPc(sd_gauss / 10.0);
        segSpots.setRadiusLocalMean(r0, r1, r2, weight);
        segSpots.setWatershed(false);
        // reset labelimage
        if (objs[r] != null) {
            objs[r].draw(labelImage, 0);
        // IJ.log("deleting old spot " + r);
        }
        // redraw population (in case spot has erased objects)
        popOrig.draw(labelImage);
        segSpots.setLabelImage(labelImage);
        ArrayList<Voxel3D> spot = segSpots.segmentOneSpot(x, y, z, col);
        // IJ.log("spot " + r + " at " + x + " " + y + " " + z + " : " + spot);
        if (spot != null) {
            // IJ.log("spot found, size " + spot.size());
            if (spot.size() > 0) {
                objs[r] = new Object3DVoxels(spot);
                objs[r].draw(labelImage, col);
            }
            if (plusres == null) {
                plusres = labelImage.getImagePlus();
                plusres.show();
                IJ.run(plusres, "3-3-2 RGB", "");
            }
            IJ.setMinAndMax(plusres, 0, col);
            plusres.setSlice((int) objs[r].getCenterZ() + 1);
            plusres.updateAndDraw();
        } else {
            jTextFieldLog.setText(jTextFieldLog.getText() + "\nno spot found !");
        }
    }
    labelImage.updateDisplay();
}

Example 22 with Voxel3D

use of mcib3d.geom.Voxel3D in project mcib3d-core by mcib3d.

the class MaximaFinder method computePeaks.

private void computePeaks() {
    imagePeaks = new ImageShort("peaks", img.sizeX, img.sizeY, img.sizeZ);
    if (verbose) {
        IJ.log("Finding all peaks");
    }
    ArrayList<Voxel3DComparable> maximaTmp = FastFilters3D.getListMaxima(img, radXY, radXY, radZ, nbCpus, false);
    Collections.sort(maximaTmp);
    for (Voxel3DComparable V : maximaTmp) {
        imagePeaks.setPixel(V, (float) V.getValue());
    }
    if (verbose) {
        IJ.log(maximaTmp.size() + " peaks found");
    }
    if (verbose) {
        IJ.log("Removing peaks below noise");
    }
    maxima = new ArrayList<Voxel3D>();
    int c = 1;
    int nb = maximaTmp.size();
    Date start = new Date();
    Date temp;
    for (Voxel3DComparable V : maximaTmp) {
        if (img.getPixel(V) > 0) {
            if (V.getValue() > noiseTolerance) {
                maxima.add(V);
                if (verbose) {
                    temp = new Date();
                    if ((temp.getTime() - start.getTime()) > 100) {
                        IJ.showStatus("Processing peak " + c + "/" + nb + " " + V);
                        start = new Date();
                    }
                    c++;
                }
                Flood3D.flood3DNoise26(img, new IntCoord3D(V.getRoundX(), V.getRoundY(), V.getRoundZ()), (int) (Math.max(1, V.getValue() - noiseTolerance)), 0);
            }
        }
    }
    if (verbose) {
        IJ.log(maxima.size() + " peaks found");
    }
    if (verbose) {
        IJ.log("Creating final peaks");
    }
    imagePeaks.fill(0);
    for (Voxel3D V : maxima) {
        imagePeaks.setPixel(V, (float) V.getValue());
    }
    if (verbose) {
        IJ.log("MaximaFinder3D finished.");
    }
}

Example 23 with Voxel3D

use of mcib3d.geom.Voxel3D in project mcib3d-core by mcib3d.

the class RegionGrowing method mergeLabels.

public int mergeLabels(int[] oldLabels, ArrayList<Voxel3D> borders) {
    // get borders between these labels
    // ArrayList<Voxel3D> borders = this.getBordersMerge(oldLabels, outsideBorders);
    watershedImage.resetStats(null);
    int newLabel = (int) watershedImage.getMax() + 1;
    for (Voxel3D V : borders) {
        watershedImage.setPixel(V, newLabel);
    }
    // update image
    watershedImage.replacePixelsValue(oldLabels, newLabel);
    // update associations
    if (computeAssociation) {
        ArrayList<String> toRemove = new ArrayList();
        ArrayList<String> toAdd = new ArrayList();
        for (String S : associations) {
            boolean reprep = false;
            String[] asso = S.split("_");
            ArrayUtil assoTab = new ArrayUtil(asso.length);
            for (int i = 0; i < assoTab.getSize(); i++) {
                int val = Integer.parseInt(asso[i]);
                boolean rep = false;
                for (int old : oldLabels) {
                    if (val == old) {
                        rep = true;
                        reprep = true;
                        break;
                    }
                }
                if (rep) {
                    assoTab.putValue(i, newLabel);
                } else {
                    assoTab.putValue(i, val);
                }
            }
            if (reprep) {
                assoTab = assoTab.distinctValues();
                // array from max to min with unique values
                assoTab.reverse();
                if (assoTab.getSize() == 1) {
                    toRemove.add(S);
                } else {
                    String assoRep = "" + (int) assoTab.getValue(0);
                    for (int i = 1; i < assoTab.getSize(); i++) {
                        assoRep = assoRep.concat("_" + (int) assoTab.getValue(i));
                    }
                    toRemove.add(S);
                    toAdd.add(assoRep);
                }
            }
        }
        if (!toRemove.isEmpty()) {
            associations.removeAll(toRemove);
        }
        if (!toAdd.isEmpty()) {
            for (String A : toAdd) {
                // test if association exists already
                boolean ok = true;
                for (String S : associations) {
                    if ((S.compareTo(A)) == 0) {
                        ok = false;
                        break;
                    }
                }
                if (ok) {
                    associations.add(A);
                }
            }
        }
        // With allRegions
        assoRegions.replaceRegion(oldLabels, newLabel);
    }
    // ??
    if (computeUpdatedLabels) {
    }
    // update volumes
    if (computeVolumes) {
        double newVol = 0;
        for (int i : oldLabels) {
            newVol += volumeLabels.get(i);
            // volumeLabels.get(i)[1] = 0;
            volumeLabels.set(i, 0.0);
        }
        newVol += borders.size();
        volumeLabels.add(newVol);
    }
    return newLabel;
}

Example 24 with Voxel3D

use of mcib3d.geom.Voxel3D in project mcib3d-core by mcib3d.

the class Watershed3D method createNeigList.

private void createNeigList() {
    voxels = new LinkedList<Voxel3DComparable>();
    int sx = rawImage.sizeX;
    int sy = rawImage.sizeY;
    int sz = rawImage.sizeZ;
    // watershed images
    watershedImage = new ImageShort("watershed", sx, sy, sz);
    labelQueueImage = new ImageShort("labelQ", sx, sy, sz);
    // okseeds = false;
    float pix;
    float se;
    // compute the seeds image
    // threshold, // TODO 32-bits seeds ?
    ImageHandler seedsLabel = seedsImage.duplicate();
    seedsLabel.thresholdCut(seedsThreshold, false, true);
    if ((labelSeeds)) {
        log.log("Labelling ");
        ImageLabeller labeller = new ImageLabeller();
        seedsLabel = labeller.getLabels(seedsLabel);
    }
    // since seeds Label starts at 1 and watershed at 2, replace values
    int max = (int) seedsLabel.getMax();
    if (seedsLabel.hasOneValue(QUEUE)) {
        seedsLabel.replacePixelsValue(QUEUE, max + 1);
        seedsValue.put(max + 1, QUEUE);
    }
    if (seedsLabel.hasOneValue(DAM)) {
        seedsLabel.replacePixelsValue(DAM, max + 2);
        seedsValue.put(max + 2, DAM);
    }
    if (!seedsValue.isEmpty())
        seedsLabel.resetStats(null);
    for (int z = sz - 1; z >= 0; z--) {
        IJ.showStatus("Processing watershed " + (z + 1));
        for (int y = sy - 1; y >= 0; y--) {
            for (int x = sx - 1; x >= 0; x--) {
                pix = rawImage.getPixel(x, y, z);
                se = seedsLabel.getPixel(x, y, z);
                if (pix > rawThreshold) {
                    if (se > 0) {
                        watershedImage.setPixel(x, y, z, se);
                        // okseeds = true;
                        ArrayList<Voxel3D> list = watershedImage.getNeighborhood3x3x3ListNoCenter(x, y, z);
                        Collections.shuffle(list);
                        for (Voxel3D N : list) {
                            int vx = (int) N.getX();
                            int vy = (int) N.getY();
                            int vz = (int) N.getZ();
                            float raw = rawImage.getPixel(vx, vy, vz);
                            if ((raw > rawThreshold) && (seedsLabel.getPixel(vx, vy, vz) == 0) && (watershedImage.getPixel(vx, vy, vz) != QUEUE)) {
                                voxels.add(new Voxel3DComparable(vx, vy, vz, raw, se));
                                watershedImage.setPixel(vx, vy, vz, QUEUE);
                                labelQueueImage.setPixel(vx, vy, vz, se);
                            }
                        }
                    }
                }
            }
        }
    }
    log.log("Watershed...");
}

Example 25 with Voxel3D

use of mcib3d.geom.Voxel3D in project mcib3d-core by mcib3d.

the class Segment3DSpots method splitSpotProjection.

/**
 * NOT USED, see splitSpotWatershed
 * @param ori
 * @param val
 * @param f
 * @param new1
 * @param new2
 * @param dist
 * @param rad
 * @param dim
 * @return
 */
public static ArrayList<Voxel3D>[] splitSpotProjection(IntImage3D ori, int val, Object3D f, int new1, int new2, double dist, float rad, int dim) {
    // only for voxels object
    if (!(f instanceof Object3DVoxels)) {
        return null;
    }
    boolean debug = false;
    int p1 = f.getZmin();
    int p2 = f.getZmax();
    IntImage3D proj;
    if (dim == 2) {
        proj = ori.projectionZ(-1, p1, p2);
    } else {
        proj = ori;
    }
    // proj = (IntImage3D) proj.medianFilter(1, 1, 0);
    float radz = (dim == 2) ? 0 : rad;
    IntImage3D maxi = proj.createLocalMaximaImage(rad, rad, radz, false);
    // IntImage3D seg = f.getSegImage();
    if (debug) {
        System.out.println("Separe2D " + val);
    }
    // IntImage3D seg2D = seg.projectionZ(val, p1, p2);
    double cx1;
    double cy1;
    double cz1;
    double nb1;
    double cx2;
    double cy2;
    double cz2;
    double nb2;
    double dist1;
    double dist2;
    int nb;
    if (debug) {
        System.out.println("separe spots 2D : " + val);
    }
    // compter nb de max locaux
    /*
         * for (int j = 2; j < maxi.getSizey() - 2; j++) { for (int i = 2; i <
         * maxi.getSizex() - 2; i++) { if ((maxi.getPixel(i, j, 0) > 0) &&
         * (seg2D.getPixel(i, j, 0) > 0)) { maxi.putPixel(i, j, 0, 255); nb++; }
         * else { } } }
         *
         */
    // with ArrayList
    ArrayList<Voxel3D> list = f.getVoxels();
    ArrayList<Voxel3D> maxlo = new ArrayList<Voxel3D>();
    Iterator it = list.iterator();
    int xx, yy, zz;
    Voxel3D vox;
    Voxel3D vox1;
    while (it.hasNext()) {
        vox = (Voxel3D) it.next();
        xx = vox.getRoundX();
        yy = vox.getRoundY();
        zz = vox.getRoundZ();
        if (dim == 2) {
            zz = 0;
        }
        if (maxi.getPix(xx, yy, zz) > 0) {
            // maxi.putPixel(xx, yy, 0, 255);
            vox1 = new Voxel3D(vox);
            vox1.setZ(0);
            maxlo.add(vox1);
        // nb++;
        // } else {
        // maxi.putPixel(xx, yy, 0, 0);
        // }
        }
    }
    nb = maxlo.size();
    IJ.log("max loco " + nb);
    if (debug) {
        new ImagePlus("max loco-" + val, maxi.getStack()).show();
        // new ImagePlus("seg2D-" + val, seg2D.getStack()).show();
        // new ImagePlus("seg-" + val, seg.getStack()).show();
        new ImagePlus("proj-" + val, proj.getStack()).show();
    }
    // si 1 seul pixel on retourne false
    if (maxlo.size() < 2) {
        return null;
    }
    /*
         * Voxel3D tablo[] = new Voxel3D[nb]; nb = 0; // ranger max locaux dans
         * le tablo for (int j = 0; j < maxi.getSizey(); j++) { for (int i = 0;
         * i < maxi.getSizex(); i++) { if ((maxi.getPixel(i, j, 0) > 0)) {
         * tablo[nb] = new Voxel3D(i, j, 0, 0); nb++; } } } /*
         */
    // chercher deux max locaux les plus eloignes
    int i1 = 0;
    int i2 = 0;
    double dmax = 0.0;
    double d;
    for (int i = 0; i < nb; i++) {
        for (int j = i + 1; j < nb; j++) {
            d = maxlo.get(i).distance(maxlo.get(j));
            if (d > dmax) {
                dmax = d;
                i1 = i;
                i2 = j;
            }
        }
    }
    // Ranger les max locaux en deux classes et calculer barycentre des deux classes
    cx1 = 0;
    cy1 = 0;
    cz1 = 0;
    cx2 = 0;
    cy2 = 0;
    cz2 = 0;
    double d1;
    double d2;
    /*
         * if (debug) { System.out.println("max locaux eloignes:");
         * System.out.println("centre1: " + tablo[i1].getX() + " " +
         * tablo[i1].getY() + " " + tablo[i1].getZ());
         * System.out.println("centre2: " + tablo[i2].getX() + " " +
         * tablo[i2].getY() + " " + tablo[i2].getZ());
         * System.out.println("distance : " + tablo[i1].distance(tablo[i2])); }
         *
         */
    Voxel3D PP1 = new Voxel3D(maxlo.get(i1).getX(), maxlo.get(i1).getY(), maxlo.get(i1).getZ(), 0.0);
    Voxel3D PP2 = new Voxel3D(maxlo.get(i2).getX(), maxlo.get(i2).getY(), maxlo.get(i2).getZ(), 0.0);
    Voxel3D P1 = new Voxel3D(cx1, cy1, cz1, 0);
    Voxel3D P2 = new Voxel3D(cx2, cy2, cz2, 0);
    int nbite = 0;
    while (((P1.distance(PP1) > 1) || (P2.distance(PP2) > 1)) && (nbite < 100)) {
        nbite++;
        cx1 = 0;
        cy1 = 0;
        cx2 = 0;
        cy2 = 0;
        cz1 = 0;
        cz2 = 0;
        nb1 = 0;
        nb2 = 0;
        P1.setX(PP1.getX());
        P1.setY(PP1.getY());
        P1.setZ(PP1.getZ());
        P2.setX(PP2.getX());
        P2.setY(PP2.getY());
        P2.setZ(PP2.getZ());
        for (int i = 0; i < nb; i++) {
            d1 = P1.distance(maxlo.get(i));
            d2 = P2.distance(maxlo.get(i));
            if (d1 < d2) {
                cx1 += maxlo.get(i).getX();
                cy1 += maxlo.get(i).getY();
                cz1 += maxlo.get(i).getZ();
                nb1++;
            } else {
                cx2 += maxlo.get(i).getX();
                cy2 += maxlo.get(i).getY();
                cz2 += maxlo.get(i).getZ();
                nb2++;
            }
        }
        cx1 /= nb1;
        cy1 /= nb1;
        cx2 /= nb2;
        cy2 /= nb2;
        cz1 /= nb1;
        cz2 /= nb2;
        PP1.setX(cx1);
        PP1.setY(cy1);
        PP2.setX(cx2);
        PP2.setY(cy2);
        PP1.setZ(cz1);
        PP2.setZ(cz2);
    }
    if (debug) {
        System.out.println("max locaux centres:");
        System.out.println("centre1: " + cx1 + " " + cy1 + " " + cz1);
        System.out.println("centre2: " + cx2 + " " + cy2 + " " + cz2);
        System.out.println("distance: " + PP1.distance(PP2));
    }
    // prendre z milieu
    if (dim == 2) {
        double zmin1 = f.getZmax();
        double zmax1 = f.getZmin();
        double zmin2 = f.getZmax();
        double zmax2 = f.getZmin();
        double di, z;
        // With ArrayList
        it = list.iterator();
        while (it.hasNext()) {
            vox = (Voxel3D) it.next();
            z = vox.getZ();
            // PP1
            PP1.setZ(z);
            di = PP1.distanceSquare(vox);
            if (di < 1) {
                if (z < zmin1) {
                    zmin1 = z;
                }
                if (z > zmax1) {
                    zmax1 = z;
                }
            }
            // PP2
            PP2.setZ(z);
            di = PP2.distanceSquare(vox);
            if (di < 1) {
                if (z < zmin2) {
                    zmin2 = z;
                }
                if (z > zmax2) {
                    zmax2 = z;
                }
            }
        }
        cz1 = 0.5 * (zmin1 + zmax1);
        cz2 = 0.5 * (zmin2 + zmax2);
        PP1.setZ(cz1);
        PP2.setZ(cz2);
    }
    /*
         * // z pour spot 2 k = f.getZmin(); while ((seg.getPixel((int) cx2,
         * (int) cy2, k) != val) && (k < sizez)) { k++; } if (k >= sizez) { k1 =
         * f.getZmin(); k2 = f.getZmax(); } else { k1 = k; while
         * ((seg.getPixel((int) cx2, (int) cy2, k) == val) && (k < sizez)) {
         * k++; } if (k >= sizez) { k2 = f.getZmax(); } else { k2 = k; } } cz2 =
         * (k1 + k2) / 2.0;
         *
         */
    if (debug) {
        System.out.println("max locaux centres en Z:");
        System.out.println("centre1: " + cx1 + " " + cy1 + " " + cz1);
        System.out.println("centre2: " + cx2 + " " + cy2 + " " + cz2);
        System.out.println("distance: " + PP1.distance(PP2));
    }
    // si distance trop petite, on retourne false
    if (PP1.distance(PP2) < dist) {
        return null;
    }
    // hybrid watershed test
    double db1 = f.distPixelBorder(PP1);
    double db2 = f.distPixelBorder(PP2);
    // System.out.println("db " + db1 + " " + db2);
    double db = db1 - db2;
    boolean oneIsLarger = true;
    // object 2 larger than object 1
    if (db2 > db1) {
        db = db2 - db1;
        oneIsLarger = false;
    }
    // //////////////////////////////////////////////////////
    // /// separate the two objects   ///////////////////////
    // //////////////////////////////////////////////////////
    // ArrayUtil tab;
    // IntImage3D seg1 = new IntImage3D(sizex, sizey, sizez);
    // IntImage3D seg2 = new IntImage3D(sizex, sizey, sizez);
    ArrayList<Voxel3D> ob1 = new ArrayList<Voxel3D>();
    ArrayList<Voxel3D> ob2 = new ArrayList<Voxel3D>();
    // IntImage3D tmp;
    // with ArrayList
    it = list.iterator();
    while (it.hasNext()) {
        vox = (Voxel3D) it.next();
        dist1 = PP1.distance(vox);
        dist2 = PP2.distance(vox);
        if (oneIsLarger) {
            dist1 -= db;
        } else {
            dist2 -= db;
        }
        if (dist1 < dist2) {
            vox1 = new Voxel3D(vox);
            vox1.setValue(new1);
            ob1.add(new Voxel3D(vox1));
        } else {
            vox1 = new Voxel3D(vox);
            vox1.setValue(new2);
            ob2.add(new Voxel3D(vox1));
        }
    }
    /*
         * for (k = 0; k < seg.getSizez(); k++) { for (int j = 0; j <
         * seg.getSizey(); j++) { for (int i = 0; i < seg.getSizex(); i++) { if
         * (seg.getPixel(i, j, k) == val) { dist1 = (i - cx1) * (i - cx1) + (j -
         * cy1) * (j - cy1) + (k - cz1) * (k - cz1); dist2 = (i - cx2) * (i -
         * cx2) + (j - cy2) * (j - cy2) + (k - cz2) * (k - cz2); if
         * (oneIsLarger) { dist1 -= db; } else { dist2 -= db; }
         *
         * }
         * }
         * }
         * }
         */
    ArrayList[] tab = new ArrayList[2];
    tab[0] = ob1;
    tab[1] = ob2;
    return tab;
}