Search in sources :

Example 26 with Line

use of org.apache.commons.math3.geometry.euclidean.twod.Line in project chordatlas by twak.

the class GISGen method importMesh.

private void importMesh(int index) {
    LoopL<Point3d> polies = blocks.get(index);
    List<Vector2D> verts = polies.stream().flatMap(ll -> ll.streamAble()).map(x -> {
        Line3d l = new Line3d(x.get(), x.getNext().get());
        l.move(perp(l.dir(), EXPAND_MESH));
        return new Vector2D(l.start.x, l.start.z);
    }).collect(Collectors.toList());
    double tol = 0.0001;
    ConvexHull2D chull = null;
    while (tol < 10) {
        try {
            chull = new MonotoneChain(false, tol).generate(verts);
            tol = 1000;
        } catch (ConvergenceException e) {
            tol *= 10;
        }
    }
    if (chull == null) {
        System.out.println("unable to find hull");
        return;
    }
    Loop<Point3d> hull = new Loop<Point3d>((Arrays.stream(chull.getLineSegments()).map(x -> new Point3d(x.getStart().getX(), 0, x.getStart().getY())).collect(Collectors.toList())));
    File root = new File(Tweed.SCRATCH + "meshes" + File.separator);
    int i = 0;
    File l;
    while ((l = new File(root, "" + i)).exists()) i++;
    l.mkdirs();
    File croppedFile = new File(l, CROPPED_OBJ);
    boolean found = false;
    for (Gen gen : tweed.frame.gens(MiniGen.class)) {
        // minigen == optimised obj
        ((MiniGen) gen).clip(hull, croppedFile);
        found = true;
    }
    if (!found)
        for (Gen gen : tweed.frame.gens(MeshGen.class)) {
            // obj == just import whole obj
            ObjGen objg = (ObjGen) gen;
            try {
                Files.asByteSource(objg.getFile()).copyTo(Files.asByteSink(croppedFile));
                objg.setVisible(false);
                found = true;
            } catch (IOException e) {
                e.printStackTrace();
            }
        }
    if (found) {
        Graph2D g2 = new Graph2D();
        polies.stream().flatMap(ll -> ll.streamAble()).forEach(x -> g2.add(new Point2d(x.get().x, x.get().z), new Point2d(x.getNext().get().x, x.getNext().get().z)));
        g2.removeInnerEdges();
        // new Plot (true, g2 );
        UnionWalker uw = new UnionWalker();
        for (Point2d p : g2.map.keySet()) for (Line line : g2.map.get(p)) uw.addEdge(line.end, line.start);
        // new Plot (true, new ArrayList( uw.map.keySet()) );
        Loopz.writeXZObj(uw.findAll(), new File(l, "gis.obj"), true);
        Loopz.writeXZObj(Loopz.to2dLoop(polies, 1, null), new File(l, "gis_footprints.obj"), false);
        BlockGen bg = new BlockGen(l, tweed, polies);
        lastMesh.put(index, bg);
        tweed.frame.addGen(bg, true);
        tweed.frame.setSelected(bg);
    } else
        JOptionPane.showMessageDialog(tweed.frame(), "Failed to find mesh from minimesh or gml layers");
}
Also used : FactoryException(org.opengis.referencing.FactoryException) Arrays(java.util.Arrays) CRS(org.geotools.referencing.CRS) NoSuchAuthorityCodeException(org.opengis.referencing.NoSuchAuthorityCodeException) Matrix4d(javax.vecmath.Matrix4d) Closer(org.twak.viewTrace.Closer) ConvexHull2D(org.apache.commons.math3.geometry.euclidean.twod.hull.ConvexHull2D) Loop(org.twak.utils.collections.Loop) ObjRead(org.twak.utils.geom.ObjRead) Complete(org.twak.utils.Parallel.Complete) Map(java.util.Map) Point3d(javax.vecmath.Point3d) LoopL(org.twak.utils.collections.LoopL) SuperLoop(org.twak.utils.collections.SuperLoop) Parallel(org.twak.utils.Parallel) ListDownLayout(org.twak.utils.ui.ListDownLayout) Line(org.twak.utils.Line) Set(java.util.Set) ConvergenceException(org.apache.commons.math3.exception.ConvergenceException) DefaultGeocentricCRS(org.geotools.referencing.crs.DefaultGeocentricCRS) Collectors(java.util.stream.Collectors) SatUtils(org.twak.footprints.SatUtils) FacadeFinder(org.twak.viewTrace.FacadeFinder) List(java.util.List) Optional(java.util.Optional) GMLReader(org.twak.viewTrace.GMLReader) Line3d(org.twak.utils.geom.Line3d) JPanel(javax.swing.JPanel) GenHandlesSelect(org.twak.tweed.GenHandlesSelect) FacadeTool(org.twak.tweed.tools.FacadeTool) HashMap(java.util.HashMap) Graph2D(org.twak.utils.geom.Graph2D) Pair(org.twak.utils.Pair) Vector3d(javax.vecmath.Vector3d) Vector2D(org.apache.commons.math3.geometry.euclidean.twod.Vector2D) Tweed(org.twak.tweed.Tweed) ArrayList(java.util.ArrayList) TweedSettings(org.twak.tweed.TweedSettings) HashSet(java.util.HashSet) Files(com.google.common.io.Files) SelectTool(org.twak.tweed.tools.SelectTool) NoSuchElementException(java.util.NoSuchElementException) JComponent(javax.swing.JComponent) UnionWalker(org.twak.utils.geom.UnionWalker) Work(org.twak.utils.Parallel.Work) MonotoneChain(org.apache.commons.math3.geometry.euclidean.twod.hull.MonotoneChain) IOException(java.io.IOException) JOptionPane(javax.swing.JOptionPane) File(java.io.File) Loopz(org.twak.utils.collections.Loopz) Point2d(javax.vecmath.Point2d) ConvexHull2D(org.apache.commons.math3.geometry.euclidean.twod.hull.ConvexHull2D) Line3d(org.twak.utils.geom.Line3d) MonotoneChain(org.apache.commons.math3.geometry.euclidean.twod.hull.MonotoneChain) Point3d(javax.vecmath.Point3d) ConvergenceException(org.apache.commons.math3.exception.ConvergenceException) UnionWalker(org.twak.utils.geom.UnionWalker) Loop(org.twak.utils.collections.Loop) SuperLoop(org.twak.utils.collections.SuperLoop) IOException(java.io.IOException) Graph2D(org.twak.utils.geom.Graph2D) Line(org.twak.utils.Line) Vector2D(org.apache.commons.math3.geometry.euclidean.twod.Vector2D) Point2d(javax.vecmath.Point2d) File(java.io.File)

Example 27 with Line

use of org.apache.commons.math3.geometry.euclidean.twod.Line in project chordatlas by twak.

the class Prof method parameterize.

public Prof parameterize() {
    // find and subtract vertical lines
    Set<Line> lines = new HashSet<>();
    for (int i = 1; i < size(); i++) lines.add(new Line(get(i - 1), get(i)));
    double avgMinY = get(0).y;
    SliceParameters P = new SliceParameters(5);
    P.FL_REGRESS = false;
    P.FL_BINS = 20;
    // simple = 0.4
    double A = 0.4;
    // simple = 1;
    double B = 1;
    lines = new FindLines(lines, P) {

        protected double nextAngle(Set<Line> remaining, int iteration) {
            double delta = Math.PI / P.FL_BINS;
            // angle bin
            Bin<Line> aBin = new Bin(-Math.PI - delta, Math.PI + delta, P.FL_BINS * 2, true);
            for (Line l : remaining) {
                double len = l.length();
                double angle = l.aTan2();
                aBin.add(angle, len, l);
            }
            // return MUtils.PI2;
            if (iteration < 1 && aBin.getWeight(Mathz.PI2) >= 5)
                return Mathz.PI2;
            int aBinI = aBin.maxI();
            return aBin.val(aBinI);
        }

        protected double getTolNearLine(Point2d p) {
            return P.FL_NEAR_LINE * (p.y < avgMinY + 5 ? 5 : B);
        }

        protected double getTolNearLine2(Point2d p) {
            return P.FL_NEAR_LINE_2 * (p.y < avgMinY + 5 ? 10 : B);
        }
    }.result.all;
    clean = new Prof(this);
    clean.clear();
    if (lines.isEmpty()) {
        clean.add(new Point2d(0, 0));
        clean.add(new Point2d(0, 1));
        return clean;
    }
    List<Line> llines = new ArrayList(lines);
    llines.stream().filter(l -> l.start.y > l.end.y).forEach(l -> l.reverseLocal());
    Collections.sort(llines, new Comparator<Line>() {

        public int compare(Line o1, Line o2) {
            return Double.compare(o1.start.y + o1.end.y, o2.start.y + o2.end.y);
        }
    });
    for (int i = 0; i < llines.size(); i++) {
        Line l = llines.get(i);
        double angle = l.aTan2();
        if (angle < Mathz.PI2 + 0.1 && angle > Mathz.PI2 - 0.4)
            llines.set(i, FindLines.rotateToAngle(l, l.fromPPram(0.5), Mathz.PI2));
    }
    Line bottomLine = llines.get(0);
    llines.add(0, new Line(new Point2d(bottomLine.start.x, get(0).y), new Point2d(bottomLine.start.x, get(0).y)));
    double lastY = -Double.MAX_VALUE, lastX = Double.MAX_VALUE;
    for (Line l : llines) {
        boolean startAbove = l.start.y >= lastY && l.start.x <= lastX, endAbove = l.end.y >= lastY && l.end.x <= lastX;
        if (startAbove && endAbove) {
            clean.add(l.start);
            clean.add(l.end);
        } else if (!startAbove && endAbove) {
            if (l.start.y < lastY) {
                Point2d sec = new LinearForm(new Vector2d(1, 0)).findC(new Point2d(0, lastY)).intersect(new LinearForm(l));
                if (sec != null)
                    l.start = sec;
            }
            if (l.start.x > lastX) {
                Point2d sec = new LinearForm(new Vector2d(0, 1)).findC(new Point2d(lastX, 0)).intersect(new LinearForm(l));
                if (sec != null)
                    l.start = sec;
            }
            if (l.end.distanceSquared(l.start) < 100)
                clean.add(l.start);
            clean.add(l.end);
        } else {
            Vector2d dir = l.dir();
            if (Math.abs(dir.x) > Math.abs(dir.y)) {
                LinearForm x = new LinearForm(new Vector2d(1, 0)).findC(new Point2d(0, lastY));
                l.start = x.project(l.start);
                l.end = x.project(l.end);
                l.start.x = Math.min(l.start.x, lastX);
                l.end.x = Math.min(l.end.x, lastX);
            } else {
                LinearForm y = new LinearForm(new Vector2d(0, 1)).findC(new Point2d(lastX, 0));
                l.start = y.project(l.start);
                l.end = y.project(l.end);
                l.start.y = Math.max(l.start.y, lastY);
                l.end.y = Math.max(l.end.y, lastY);
            }
            clean.add(l.start);
            clean.add(l.end);
        }
        lastY = l.end.y;
        lastX = l.end.x;
    }
    clean.clean(0.2);
    return clean;
}
Also used : ConsecutivePairs(org.twak.utils.collections.ConsecutivePairs) Matrix4d(javax.vecmath.Matrix4d) ConsecutiveItPairs(org.twak.utils.collections.ConsecutiveItPairs) SliceParameters(org.twak.viewTrace.SliceParameters) Arrayz(org.twak.utils.collections.Arrayz) Node(com.jme3.scene.Node) ObjRead(org.twak.utils.geom.ObjRead) Map(java.util.Map) Cache(org.twak.utils.Cache) Material(com.jme3.material.Material) Point3d(javax.vecmath.Point3d) VertexBuffer(com.jme3.scene.VertexBuffer) IdentityHashMap(java.util.IdentityHashMap) Collection(java.util.Collection) Line(org.twak.utils.Line) FindLines(org.twak.viewTrace.FindLines) Set(java.util.Set) Vector2d(javax.vecmath.Vector2d) LinearForm(org.twak.utils.geom.LinearForm) Collectors(java.util.stream.Collectors) List(java.util.List) Rainbow(org.twak.utils.ui.Rainbow) Line3d(org.twak.utils.geom.Line3d) Mesh(com.jme3.scene.Mesh) Geometry(com.jme3.scene.Geometry) DBSCANClusterer(org.apache.commons.math3.ml.clustering.DBSCANClusterer) LinearForm3D(org.twak.utils.geom.LinearForm3D) Bin(org.twak.viewTrace.Bin) Pair(org.twak.utils.Pair) Vector3d(javax.vecmath.Vector3d) Clusterable(org.apache.commons.math3.ml.clustering.Clusterable) Tweed(org.twak.tweed.Tweed) ArrayList(java.util.ArrayList) TweedSettings(org.twak.tweed.TweedSettings) HashSet(java.util.HashSet) PanMouseAdaptor(org.twak.utils.PanMouseAdaptor) Graphics2D(java.awt.Graphics2D) ICanPaintU(org.twak.utils.PaintThing.ICanPaintU) Mathz(org.twak.utils.Mathz) PaintThing(org.twak.utils.PaintThing) Iterator(java.util.Iterator) Point2d(javax.vecmath.Point2d) SuperLine(org.twak.viewTrace.SuperLine) Cluster(org.apache.commons.math3.ml.clustering.Cluster) ColorRGBA(com.jme3.math.ColorRGBA) Comparator(java.util.Comparator) InAxDoubleArray(org.twak.utils.streams.InAxDoubleArray) Collections(java.util.Collections) ObjSlice(org.twak.viewTrace.ObjSlice) FindLines(org.twak.viewTrace.FindLines) Bin(org.twak.viewTrace.Bin) ArrayList(java.util.ArrayList) LinearForm(org.twak.utils.geom.LinearForm) Line(org.twak.utils.Line) SuperLine(org.twak.viewTrace.SuperLine) Vector2d(javax.vecmath.Vector2d) SliceParameters(org.twak.viewTrace.SliceParameters) Point2d(javax.vecmath.Point2d) HashSet(java.util.HashSet)

Example 28 with Line

use of org.apache.commons.math3.geometry.euclidean.twod.Line in project SeqMonk by s-andrews.

the class GeneSetIntensityDifferenceFilter method generateProbeList.

/* (non-Javadoc)
	 * @see uk.ac.babraham.SeqMonk.Filters.ProbeFilter#generateProbeList()
	 */
protected void generateProbeList() {
    try {
        applyMultipleTestingCorrection = optionsPanel.multipleTestingBox.isSelected();
        calculateCustomRegression = optionsPanel.calculateCustomRegressionBox.isSelected();
        if (calculateCustomRegression == false) {
            customRegressionValues = null;
        }
        if (calculateLinearRegression) {
            simpleRegression = new SimpleRegression();
        }
        Probe[] allProbes = startingList.getAllProbes();
        // We're not allowing multiple comparisons - this is a bit of a messy workaround so it's compatible with other methods.
        fromStore = fromStores[0];
        toStore = toStores[0];
        ArrayList<Probe> probeArrayList = new ArrayList<Probe>();
        int NaNcount = 0;
        // remove the invalid probes - the ones without a value
        for (int i = 0; i < allProbes.length; i++) {
            if ((Float.isNaN(fromStore.getValueForProbe(allProbes[i]))) || (Float.isNaN(toStore.getValueForProbe(allProbes[i])))) {
                NaNcount++;
            } else {
                probeArrayList.add(allProbes[i]);
            }
        }
        System.err.println("Found " + NaNcount + " probes that were invalid.");
        probes = probeArrayList.toArray(new Probe[0]);
        if (calculateCustomRegression == true) {
            customRegressionValues = new float[2][probes.length];
        }
        // We'll pull the number of probes to sample from the preferences if they've changed it
        Integer updatedProbesPerSet = optionsPanel.probesPerSet();
        if (updatedProbesPerSet != null)
            probesPerSet = updatedProbesPerSet;
        // we want a set of z-scores using the local distribution.
        probeZScoreLookupTable = new Hashtable<Probe, Double>();
        // Put something in the progress whilst we're ordering the probe values to make the comparison.
        progressUpdated("Generating background model", 0, 1);
        Comparator<Integer> comp = new AverageIntensityComparator(fromStore, toStore, probes);
        // We need to generate a set of probe indices that can be ordered by their average intensity
        Integer[] indices = new Integer[probes.length];
        for (int i = 0; i < probes.length; i++) {
            indices[i] = i;
            /* add the data to the linear regression object */
            if (calculateLinearRegression) {
                simpleRegression.addData((double) fromStore.getValueForProbe(probes[i]), (double) toStore.getValueForProbe(probes[i]));
            }
        }
        if (calculateLinearRegression) {
            System.out.println("intercept = " + simpleRegression.getIntercept() + ", slope = " + simpleRegression.getSlope());
        }
        Arrays.sort(indices, comp);
        /* This holds the indices to get the deduplicated probe from the original list of probes */
        deduplicatedIndices = new Integer[probes.length];
        /* The number of probes with different values */
        int dedupProbeCounter = 0;
        // the probes deduplicated by value
        ArrayList<Probe> deduplicatedProbes = new ArrayList<Probe>();
        // populate the first one so that we have something to compare to in the loop
        deduplicatedIndices[0] = 0;
        deduplicatedProbes.add(probes[indices[0]]);
        progressUpdated("Made 0 of 1 comparisons", 0, 1);
        for (int i = 1; i < indices.length; i++) {
            /* indices have been sorted, now we need to check whether adjacent pair values are identical */
            if ((fromStore.getValueForProbe(probes[indices[i]]) == fromStore.getValueForProbe(probes[indices[i - 1]])) && (toStore.getValueForProbe(probes[indices[i]]) == toStore.getValueForProbe(probes[indices[i - 1]]))) {
                /* If they are identical, do not add the probe to the deduplicatedProbes object, but have a reference for it so we can look up which deduplicated probe and 
					 * therefore which distribution slice to use for the duplicated probe. */
                deduplicatedIndices[i] = dedupProbeCounter;
            } else {
                deduplicatedProbes.add(probes[indices[i]]);
                dedupProbeCounter++;
                deduplicatedIndices[i] = dedupProbeCounter;
            }
        }
        Probe[] dedupProbes = deduplicatedProbes.toArray(new Probe[0]);
        // make sure we're not trying to use more probes than we've got in the analysis
        if (probesPerSet > dedupProbes.length) {
            probesPerSet = dedupProbes.length;
        }
        System.out.println("total number of probe values = " + probes.length);
        System.out.println("number of deduplicated probe values = " + dedupProbes.length);
        System.out.println("probesPerSet = " + probesPerSet);
        // I want this to contain all the differences, then from that I'm going to get the z-scores.
        double[] currentDiffSet = new double[probesPerSet];
        for (int i = 0; i < indices.length; i++) {
            if (cancel) {
                cancel = false;
                progressCancelled();
                return;
            }
            if (i % 1000 == 0) {
                int progress = (i * 100) / indices.length;
                // progress += 100*comparisonIndex;
                progressUpdated("Made 0 out of 1 comparisons", progress, 100);
            }
            /**
             * There are +1s in here because we skip over j when j == startingIndex.
             */
            // We need to make up the set of differences to represent this probe
            int startingIndex = deduplicatedIndices[i] - (probesPerSet / 2);
            if (startingIndex < 0)
                startingIndex = 0;
            if (startingIndex + (probesPerSet + 1) >= dedupProbes.length)
                startingIndex = dedupProbes.length - (probesPerSet + 1);
            try {
                for (int j = startingIndex; j < startingIndex + (probesPerSet + 1); j++) {
                    if (j == startingIndex) {
                        // Don't include the point being tested in the background model
                        continue;
                    }
                    double diff;
                    if (calculateLinearRegression == true) {
                        if (j > dedupProbes.length) {
                            System.err.println(" j is too big, it's " + j + " and dedupProbes.length = " + dedupProbes.length + ", starting index = " + startingIndex);
                        }
                        double x = fromStore.getValueForProbe(dedupProbes[j]);
                        double expectedY = (simpleRegression.getSlope() * x) + simpleRegression.getIntercept();
                        diff = toStore.getValueForProbe(dedupProbes[j]) - expectedY;
                    } else {
                        diff = toStore.getValueForProbe(dedupProbes[j]) - fromStore.getValueForProbe(dedupProbes[j]);
                    }
                    if (j < startingIndex) {
                        currentDiffSet[j - startingIndex] = diff;
                    } else {
                        currentDiffSet[(j - startingIndex) - 1] = diff;
                    }
                }
                if (calculateCustomRegression == true) {
                    // the average/ kind of centre line
                    float z = ((fromStore.getValueForProbe(probes[indices[i]]) + toStore.getValueForProbe(probes[indices[i]])) / 2);
                    customRegressionValues[0][indices[i]] = z - ((float) SimpleStats.mean(currentDiffSet) / 2);
                    customRegressionValues[1][indices[i]] = z + ((float) SimpleStats.mean(currentDiffSet) / 2);
                }
                double mean = 0;
                // Get the difference for this point
                double diff;
                if (calculateLinearRegression == true) {
                    double x = fromStore.getValueForProbe(probes[indices[i]]);
                    double expectedY = (simpleRegression.getSlope() * x) + simpleRegression.getIntercept();
                    diff = toStore.getValueForProbe(probes[indices[i]]) - expectedY;
                } else if (calculateCustomRegression == true) {
                    diff = toStore.getValueForProbe(probes[indices[i]]) - fromStore.getValueForProbe(probes[indices[i]]);
                    mean = SimpleStats.mean(currentDiffSet);
                } else {
                    diff = toStore.getValueForProbe(probes[indices[i]]) - fromStore.getValueForProbe(probes[indices[i]]);
                }
                double stdev;
                stdev = SimpleStats.stdev(currentDiffSet, mean);
                // if there are no reads in the probe for either of the datasets, should we set the zscore to 0??
                double zScore = (diff - mean) / stdev;
                // modified z score
                // median absolute deviation
                /*		double[] madArray = new double[currentDiffSet.length];
						double median = SimpleStats.median(currentDiffSet);					
						
						for(int d=0; d<currentDiffSet.length; d++){							
							madArray[d] = Math.abs(currentDiffSet[d] - median);							
						}
						
						double mad = SimpleStats.median(madArray);							
						zScore = (0.6745 * (diff - median))/mad;
					}
				*/
                probeZScoreLookupTable.put(probes[indices[i]], zScore);
            } catch (SeqMonkException sme) {
                progressExceptionReceived(sme);
                return;
            }
        }
        // make this an array list as we're kicking out the mapped gene sets that have zscores with variance of 0.
        ArrayList<MappedGeneSetTTestValue> pValueArrayList = new ArrayList<MappedGeneSetTTestValue>();
        MappedGeneSet[] mappedGeneSets = null;
        /* if we're using the gene set from a file, map the gene sets to the probes */
        if (optionsPanel.geneSetsFileRadioButton.isSelected()) {
            GeneSetCollectionParser geneSetCollectionParser = new GeneSetCollectionParser(minGenesInSet, maxGenesInSet);
            GeneSetCollection geneSetCollection = geneSetCollectionParser.parseGeneSetInformation(validGeneSetFilepath);
            MappedGeneSet[] allMappedGeneSets = geneSetCollection.getMappedGeneSets(probes);
            if (allMappedGeneSets == null) {
                JOptionPane.showMessageDialog(SeqMonkApplication.getInstance(), "No sets of genes could be matched to probes.\nCheck that the gmt file is for the right species. \nTo use gene sets from a file, probe names must contain the gene name. \nTry defining new probes over genes (e.g. using the RNA-seq quantitation pipeline) or use existing probes lists instead of a gene set file.", "No gene sets identified", JOptionPane.ERROR_MESSAGE);
                throw new SeqMonkException("No sets of genes could be matched to probes.\nTo use gene sets from a file, probe names must contain the gene name.\nTry defining new probes over genes or use existing probes lists instead of a gene set file.");
            } else {
                ArrayList<MappedGeneSet> mgsArrayList = new ArrayList<MappedGeneSet>();
                /* get rid of those that have fewer probes in the set than minGenesInSet. We shouldn't exceed maxGenesInSet unless probes have been made over something other than genes */
                for (int i = 0; i < allMappedGeneSets.length; i++) {
                    if (allMappedGeneSets[i].getProbes().length >= minGenesInSet) {
                        mgsArrayList.add(allMappedGeneSets[i]);
                    }
                }
                mappedGeneSets = mgsArrayList.toArray(new MappedGeneSet[0]);
            }
        } else /* or if we're using existing probelists, create mappedGeneSets from them */
        if (optionsPanel.probeListRadioButton.isSelected() && selectedProbeLists != null) {
            mappedGeneSets = new MappedGeneSet[selectedProbeLists.length];
            for (int i = 0; i < selectedProbeLists.length; i++) {
                mappedGeneSets[i] = new MappedGeneSet(selectedProbeLists[i]);
            }
        } else {
            throw new SeqMonkException("Haven't got any genesets to use, shouldn't have got here without having any selected.");
        }
        if (mappedGeneSets == null || mappedGeneSets.length == 0) {
            throw new SeqMonkException("Couldn't map gene sets to the probes, try again with a different probe set");
        } else {
            System.err.println("there are " + mappedGeneSets.length + " mappedGeneSets");
            System.err.println("size of zScore lookup table = " + probeZScoreLookupTable.size());
        }
        System.err.println("total number of mapped gene sets = " + mappedGeneSets.length);
        // deduplicate our mappedGeneSets
        if (optionsPanel.deduplicateGeneSetBox.isSelected()) {
            mappedGeneSets = deduplicateMappedGeneSets(mappedGeneSets);
        }
        System.err.println("deduplicated mapped gene sets = " + mappedGeneSets.length);
        /* 
			 * we need to go through the mapped gene set and get all the values for the matched probes 
			 * 
			 */
        for (int i = 0; i < mappedGeneSets.length; i++) {
            Probe[] geneSetProbes = mappedGeneSets[i].getProbes();
            // to contain all the z-scores for the gene set
            double[] geneSetZscores = new double[geneSetProbes.length];
            // Find the z-scores for each of the probes in the mappedGeneSet
            for (int gsp = 0; gsp < geneSetProbes.length; gsp++) {
                if (probeZScoreLookupTable.containsKey(geneSetProbes[gsp])) {
                    geneSetZscores[gsp] = probeZScoreLookupTable.get(geneSetProbes[gsp]);
                }
            }
            // this is just temporary to check the stats - there's some duplication with this.... is there??
            mappedGeneSets[i].zScores = geneSetZscores;
            if (geneSetZscores.length > 1) {
                // the number of probes in the mappedGeneSet should always be > 1 anyway as the mappedGeneSet shouldn't be created if there are < 2 matched probes.
                double pVal;
                if (optionsPanel.statisticalTestBox.getSelectedItem().toString().equals("t-test")) {
                    pVal = TTest.calculatePValue(geneSetZscores, 0);
                } else if (optionsPanel.statisticalTestBox.getSelectedItem().toString().equals("Kolmorogov-Smirnov test")) {
                    double[] allZscores = getAllZScores(probeZScoreLookupTable);
                    KolmogorovSmirnovTest ksTest = new KolmogorovSmirnovTest();
                    pVal = ksTest.kolmogorovSmirnovTest(geneSetZscores, allZscores);
                } else if (optionsPanel.statisticalTestBox.getSelectedItem().toString().equals("Kolmorogov-Smirnov non-directional test")) {
                    double[] allZscores = getAllZScores(probeZScoreLookupTable);
                    KolmogorovSmirnovTest ksTest = new KolmogorovSmirnovTest();
                    pVal = ksTest.kolmogorovSmirnovTest(convertToAbsoluteValues(geneSetZscores), convertToAbsoluteValues(allZscores));
                } else {
                    throw new IllegalArgumentException("Didn't recognise statistical test " + optionsPanel.statisticalTestBox.getSelectedItem().toString());
                }
                mappedGeneSets[i].meanZScore = SimpleStats.mean(geneSetZscores);
                // check the variance - we don't want variances of 0.
                double stdev = SimpleStats.stdev(geneSetZscores);
                if ((stdev * stdev) < 0.00000000000001) {
                    continue;
                } else // if all the differences between the datasets are 0 then just get rid of them
                if (Double.isNaN(pVal)) {
                    continue;
                } else {
                    pValueArrayList.add(new MappedGeneSetTTestValue(mappedGeneSets[i], pVal));
                }
            } else {
                System.err.println("Fell through the net somewhere, why does the set of zscores contain fewer than 2 values?");
                continue;
            }
        }
        MappedGeneSetTTestValue[] filterResultpValues = pValueArrayList.toArray(new MappedGeneSetTTestValue[0]);
        ArrayList<MappedGeneSetTTestValue> filteredPValueArrayList = new ArrayList<MappedGeneSetTTestValue>();
        // Now we've got all the p values they need to be corrected.
        if (applyMultipleTestingCorrection) {
            BenjHochFDR.calculateQValues(filterResultpValues);
        }
        System.err.println(filterResultpValues.length + " p-values calculated, multtest = " + applyMultipleTestingCorrection + ", pval limit = " + pValueLimit);
        for (int i = 0; i < filterResultpValues.length; i++) {
            double pOrQvalue;
            if (applyMultipleTestingCorrection) {
                pOrQvalue = filterResultpValues[i].q;
            } else {
                pOrQvalue = filterResultpValues[i].p;
            }
            // check whether it passes the p/q-value and z-score cut-offs
            if (optionsPanel.reportAllResults.isSelected()) {
                filteredPValueArrayList.add(filterResultpValues[i]);
            } else {
                if ((pOrQvalue < pValueLimit) && (Math.abs(filterResultpValues[i].mappedGeneSet.meanZScore) > zScoreThreshold)) {
                    filteredPValueArrayList.add(filterResultpValues[i]);
                }
            }
        }
        // convert the ArrayList to MappedGeneSetTTestValue
        filterResultpValues = filteredPValueArrayList.toArray(new MappedGeneSetTTestValue[0]);
        if (filterResultpValues.length == 0) {
            JOptionPane.showMessageDialog(SeqMonkApplication.getInstance(), "No sets of genes were identified using the selected parameters, \ntry changing the gene sets or relaxing the p-value/z-score thresholds.", "No gene sets identified", JOptionPane.INFORMATION_MESSAGE);
        } else {
            geneSetDisplay = new GeneSetDisplay(dataCollection, listDescription(), fromStore, toStore, probes, probeZScoreLookupTable, filterResultpValues, startingList, customRegressionValues, simpleRegression, // optionsPanel.bidirectionalRadioButton.isSelected());
            false);
            geneSetDisplay.addWindowListener(this);
        }
        // We don't want to save the probe list here, we're bringing up the intermediate display from which probe lists can be saved.
        progressCancelled();
    } catch (Exception e) {
        // TODO Auto-generated catch block
        e.printStackTrace();
    }
}
Also used : GeneSetCollection(uk.ac.babraham.SeqMonk.GeneSets.GeneSetCollection) ArrayList(java.util.ArrayList) Probe(uk.ac.babraham.SeqMonk.DataTypes.Probes.Probe) MappedGeneSetTTestValue(uk.ac.babraham.SeqMonk.Analysis.Statistics.MappedGeneSetTTestValue) KolmogorovSmirnovTest(org.apache.commons.math3.stat.inference.KolmogorovSmirnovTest) SeqMonkException(uk.ac.babraham.SeqMonk.SeqMonkException) GeneSetCollectionParser(uk.ac.babraham.SeqMonk.GeneSets.GeneSetCollectionParser) SeqMonkException(uk.ac.babraham.SeqMonk.SeqMonkException) SimpleRegression(org.apache.commons.math3.stat.regression.SimpleRegression) MappedGeneSet(uk.ac.babraham.SeqMonk.GeneSets.MappedGeneSet)

Example 29 with Line

use of org.apache.commons.math3.geometry.euclidean.twod.Line in project tutorials by eugenp.

the class GeometryUnitTest method whenLineIntersection_thenCorrect.

@Test
public void whenLineIntersection_thenCorrect() {
    Line l1 = new Line(new Vector2D(0, 0), new Vector2D(1, 1), 0);
    Line l2 = new Line(new Vector2D(0, 1), new Vector2D(1, 1.5), 0);
    Vector2D intersection = l1.intersection(l2);
    Assert.assertEquals(2, intersection.getX(), 1e-7);
    Assert.assertEquals(2, intersection.getY(), 1e-7);
}
Also used : Line(org.apache.commons.math3.geometry.euclidean.twod.Line) Vector2D(org.apache.commons.math3.geometry.euclidean.twod.Vector2D) Test(org.junit.Test)

Example 30 with Line

use of org.apache.commons.math3.geometry.euclidean.twod.Line in project imagej-ops by imagej.

the class DefaultMinimumFeret method calculate.

@Override
public Pair<RealLocalizable, RealLocalizable> calculate(Polygon2D input) {
    final List<? extends RealLocalizable> points = GeomUtils.vertices(function.calculate(input));
    double distance = Double.POSITIVE_INFINITY;
    RealLocalizable p0 = points.get(0);
    RealLocalizable p1 = points.get(0);
    double tmpDist = 0;
    RealLocalizable tmpP0 = p0;
    RealLocalizable tmpP1 = p1;
    for (int i = 0; i < points.size() - 2; i++) {
        final RealLocalizable lineStart = points.get(i);
        final RealLocalizable lineEnd = points.get(i + 1);
        tmpDist = 0;
        final Line l = new Line(new Vector2D(lineStart.getDoublePosition(0), lineStart.getDoublePosition(1)), new Vector2D(lineEnd.getDoublePosition(0), lineEnd.getDoublePosition(1)), 10e-12);
        for (int j = 0; j < points.size(); j++) {
            if (j != i && j != i + 1) {
                final RealLocalizable ttmpP0 = points.get(j);
                final double tmp = l.distance(new Vector2D(ttmpP0.getDoublePosition(0), ttmpP0.getDoublePosition(1)));
                if (tmp > tmpDist) {
                    tmpDist = tmp;
                    final Vector2D vp = (Vector2D) l.project(new Vector2D(ttmpP0.getDoublePosition(0), ttmpP0.getDoublePosition(1)));
                    tmpP0 = new RealPoint(vp.getX(), vp.getY());
                    tmpP1 = ttmpP0;
                }
            }
        }
        if (tmpDist < distance) {
            distance = tmpDist;
            p0 = tmpP0;
            p1 = tmpP1;
        }
    }
    final RealLocalizable lineStart = points.get(points.size() - 1);
    final RealLocalizable lineEnd = points.get(0);
    final Line l = new Line(new Vector2D(lineStart.getDoublePosition(0), lineStart.getDoublePosition(1)), new Vector2D(lineEnd.getDoublePosition(0), lineEnd.getDoublePosition(1)), 10e-12);
    tmpDist = 0;
    for (int j = 0; j < points.size(); j++) {
        if (j != points.size() - 1 && j != 0 + 1) {
            final RealLocalizable ttmpP0 = points.get(j);
            final double tmp = l.distance(new Vector2D(ttmpP0.getDoublePosition(0), ttmpP0.getDoublePosition(1)));
            if (tmp > tmpDist) {
                tmpDist = tmp;
                final Vector2D vp = (Vector2D) l.project(new Vector2D(ttmpP0.getDoublePosition(0), ttmpP0.getDoublePosition(1)));
                tmpP0 = new RealPoint(vp.getX(), vp.getY());
                tmpP1 = ttmpP0;
            }
        }
    }
    if (tmpDist < distance) {
        distance = tmpDist;
        p0 = tmpP0;
        p1 = tmpP1;
    }
    return new ValuePair<>(p0, p1);
}
Also used : RealLocalizable(net.imglib2.RealLocalizable) Line(org.apache.commons.math3.geometry.euclidean.twod.Line) Vector2D(org.apache.commons.math3.geometry.euclidean.twod.Vector2D) RealPoint(net.imglib2.RealPoint) ValuePair(net.imglib2.util.ValuePair) RealPoint(net.imglib2.RealPoint)

Aggregations

ArrayList (java.util.ArrayList)9 Line (org.twak.utils.Line)8 Plot2 (ij.gui.Plot2)7 List (java.util.List)7 Map (java.util.Map)7 Set (java.util.Set)7 Point2d (javax.vecmath.Point2d)7 HashSet (java.util.HashSet)6 Collectors (java.util.stream.Collectors)6 Pair (org.twak.utils.Pair)6 Iterator (java.util.Iterator)5 Point3d (javax.vecmath.Point3d)5 Vector2d (javax.vecmath.Vector2d)5 PointValuePair (org.apache.commons.math3.optim.PointValuePair)5 Tweed (org.twak.tweed.Tweed)5 TweedSettings (org.twak.tweed.TweedSettings)5 Material (com.jme3.material.Material)4 ColorRGBA (com.jme3.math.ColorRGBA)4 Geometry (com.jme3.scene.Geometry)4 Mesh (com.jme3.scene.Mesh)4