Compare View

switch
from
...
to
 
Commits (7)
python/classify.py
... ... @@ -12,7 +12,7 @@ import sklearn.metrics
12 12 import scipy
13 13 import scipy.misc
14 14 import envi
15   -import spectral
  15 +import hyperspectral
16 16 import random
17 17 import progressbar
18 18 import matplotlib.pyplot as plt
... ... @@ -177,7 +177,7 @@ def envi_batch_predict(E, C, batch=10000):
177 177 else:
178 178 Tv = numpy.concatenate((Tv, C.predict(Fv.transpose()).transpose()), 0)
179 179 tempmask = E.batchmask()
180   - Lv = spectral.unsift2(Tv, tempmask)
  180 + Lv = hyperspectral.unsift2(Tv, tempmask)
181 181 Cv = label2class(Lv.squeeze(), background=0)
182 182 RGB = class2color(Cv)
183 183 plt.imshow(RGB)
... ...
python/digitalstain.py
... ... @@ -5,7 +5,7 @@ Created on Tue Jul 25 16:28:37 2017
5 5 @author: david
6 6 """
7 7  
8   -import spectral
  8 +import hyperspectral
9 9 import envi
10 10 import classify
11 11 import numpy
... ... @@ -18,22 +18,36 @@ import glob
18 18 import matplotlib.pyplot as plt
19 19 import random
20 20  
21   -def generate_stain(envifile, stainfile, N=5000, batch_size=10000, validate=True):
22   - E = envi.envi(envifile)
  21 +def generate_stain(envifile, stainfile, maskfile="", trainmask="", N=5000, batch_size=10000, validate=True):
  22 + if trainmask == "":
  23 + E = envi.envi(envifile)
  24 + else:
  25 + mask = scipy.misc.imread(trainmask, flatten=True)
  26 + E = envi.envi(envifile, mask=mask)
  27 +
23 28 mask = classify.random_mask(E.mask, N)
  29 + scipy.misc.imsave("random.bmp", mask)
24 30  
25 31 Ft = E.loadmask(mask).transpose()
26 32  
27 33 stain = numpy.rollaxis(scipy.misc.imread(stainfile), 2)
28   - Tt = spectral.sift2(stain, mask).transpose()
  34 + Tt = hyperspectral.sift2(stain, mask).transpose()
29 35  
  36 + print("Training MLPRegressor...")
30 37 CLASS = sklearn.neural_network.MLPRegressor(solver='lbfgs', alpha=1e-5, hidden_layer_sizes=(), random_state=1, verbose=True)
31 38 CLASS.fit(Ft, Tt)
32 39  
33 40 if validate == False:
34 41 return CLASS
35 42  
36   - plt.ion()
  43 + print("Validating Stain...")
  44 + plt.ion()
  45 + if not maskfile == "":
  46 + E.close() #close the ENVI file
  47 + mask = scipy.misc.imread(maskfile, flatten=True)
  48 + print(numpy.count_nonzero(mask))
  49 + E = envi.envi(envifile, mask=mask)
  50 +
37 51 Fv = E.loadbatch(batch_size) #load the first batch
38 52 n = 0
39 53 while not Fv == []: #loop until an empty batch is returned
... ... @@ -41,7 +55,7 @@ def generate_stain(envifile, stainfile, N=5000, batch_size=10000, validate=True)
41 55 Tv = CLASS.predict(Fv.transpose()).transpose()
42 56 else:
43 57 Tv = numpy.append(Tv, CLASS.predict(Fv.transpose()).transpose(), 1) #append the predicted labels from this batch to those of previous batches
44   - COLORS = spectral.unsift2(Tv, E.batchmask()) #convert the matrix of class labels to a 2D array
  58 + COLORS = hyperspectral.unsift2(Tv, E.batchmask()) #convert the matrix of class labels to a 2D array
45 59 RGB = numpy.rollaxis(COLORS, 0, 3).astype(numpy.ubyte)
46 60 plt.imshow(RGB) #display it
47 61 plt.pause(0.05)
... ...
python/envi.py
... ... @@ -100,11 +100,11 @@ class envi_header:
100 100 #t = map(str.strip, t) #strip all of the tokens in the token list
101 101  
102 102 #handle the simple conditions
103   - if l[li].startswith("file type"):
104   - if not l[li].strip().endswith("ENVI Standard"):
105   - print("ERROR: unsupported ENVI file format: " + l[li].strip())
106   - return
107   - elif l[li].startswith("samples"):
  103 + #if l[li].startswith("file type"):
  104 + # if not l[li].strip().endswith("ENVI Standard"):
  105 + # print("ERROR: unsupported ENVI file format: " + l[li].strip())
  106 + # return
  107 + if l[li].startswith("samples"):
108 108 self.samples = int(l[li].split()[-1])
109 109 elif l[li].startswith("lines"):
110 110 self.lines = int(l[li].split()[-1])
... ...
python/stim_spectral.py renamed to python/hyperspectral.py
stim/biomodels/network.h
... ... @@ -473,6 +473,146 @@ public:
473 473 return E[e]; //return the specified edge (casting it to a fiber)
474 474 }
475 475  
  476 + /// subdivide current network
  477 + void subdivision() {
  478 +
  479 + std::vector<unsigned> ori_index; // original index
  480 + std::vector<unsigned> new_index; // new index
  481 + std::vector<edge> nE; // new edge
  482 + std::vector<vertex> nV; // new vector
  483 + unsigned id = 0;
  484 +
  485 + for (unsigned i = 0; i < num_edge; i++) {
  486 + if (E[i].size() == 2) { // if current edge can't be subdivided
  487 + stim::centerline<T> line(2);
  488 + for (unsigned k = 0; k < 2; k++)
  489 + line[k] = E[i][k];
  490 + line.update();
  491 +
  492 + edge new_edge(line);
  493 +
  494 + vertex new_vertex = new_edge[0];
  495 + id = E[i].v[0];
  496 + auto position = std::find(ori_index.begin(), ori_index.end(), id);
  497 + if (position == ori_index.end()) { // new vertex
  498 + ori_index.push_back(id);
  499 + new_index.push_back(nV.size());
  500 +
  501 + new_vertex.e[0].push_back(nE.size());
  502 + new_edge.v[0] = nV.size();
  503 + nV.push_back(new_vertex); // push back vertex as a new vertex
  504 + }
  505 + else { // existing vertex
  506 + int k = std::distance(ori_index.begin(), position);
  507 + new_edge.v[0] = new_index[k];
  508 + nV[new_index[k]].e[0].push_back(nE.size());
  509 + }
  510 +
  511 + new_vertex = new_edge[1];
  512 + id = E[i].v[1];
  513 + position = std::find(ori_index.begin(), ori_index.end(), id);
  514 + if (position == ori_index.end()) { // new vertex
  515 + ori_index.push_back(id);
  516 + new_index.push_back(nV.size());
  517 +
  518 + new_vertex.e[1].push_back(nE.size());
  519 + new_edge.v[1] = nV.size();
  520 + nV.push_back(new_vertex); // push back vertex as a new vertex
  521 + }
  522 + else { // existing vertex
  523 + int k = std::distance(ori_index.begin(), position);
  524 + new_edge.v[1] = new_index[k];
  525 + nV[new_index[k]].e[1].push_back(nE.size());
  526 + }
  527 +
  528 + nE.push_back(new_edge);
  529 +
  530 + nE[nE.size() - 1].cylinder<T>::set_r(0, E[i].cylinder<T>::r(0));
  531 + nE[nE.size() - 1].cylinder<T>::set_r(1, E[i].cylinder<T>::r(1));
  532 + }
  533 + else { // subdivide current edge
  534 + for (unsigned j = 0; j < E[i].size() - 1; j++) {
  535 + stim::centerline<T> line(2);
  536 + for (unsigned k = 0; k < 2; k++)
  537 + line[k] = E[i][j + k];
  538 + line.update();
  539 +
  540 + edge new_edge(line);
  541 +
  542 + if (j == 0) { // edge contains original starting point
  543 + vertex new_vertex = new_edge[0];
  544 + id = E[i].v[0];
  545 + auto position = std::find(ori_index.begin(), ori_index.end(), id);
  546 + if (position == ori_index.end()) { // new vertex
  547 + ori_index.push_back(id);
  548 + new_index.push_back(nV.size());
  549 +
  550 + new_vertex.e[0].push_back(nE.size());
  551 + new_edge.v[0] = nV.size();
  552 + nV.push_back(new_vertex); // push back vertex as a new vertex
  553 + }
  554 + else { // existing vertex
  555 + int k = std::distance(ori_index.begin(), position);
  556 + new_edge.v[0] = new_index[k];
  557 + nV[new_index[k]].e[0].push_back(nE.size());
  558 + }
  559 +
  560 + new_vertex = new_edge[1];
  561 + new_vertex.e[1].push_back(nE.size());
  562 + new_edge.v[1] = nV.size();
  563 + nV.push_back(new_vertex); // push back internal point as a new vertex
  564 +
  565 + nE.push_back(new_edge);
  566 + }
  567 +
  568 + else if (j == E[i].size() - 2) { // edge contains original ending point
  569 +
  570 + vertex new_vertex = new_edge[1];
  571 + nV[nV.size() - 1].e[0].push_back(nE.size());
  572 + new_edge.v[0] = nV.size() - 1;
  573 +
  574 + id = E[i].v[1];
  575 + auto position = std::find(ori_index.begin(), ori_index.end(), id);
  576 + if (position == ori_index.end()) { // new vertex
  577 + ori_index.push_back(id);
  578 + new_index.push_back(nV.size());
  579 +
  580 + new_vertex.e[1].push_back(nE.size());
  581 + new_edge.v[1] = nV.size();
  582 + nV.push_back(new_vertex); // push back vertex as a new vertex
  583 + }
  584 + else { // existing vertex
  585 + int k = std::distance(ori_index.begin(), position);
  586 + new_edge.v[1] = new_index[k];
  587 + nV[new_index[k]].e[1].push_back(nE.size());
  588 + }
  589 +
  590 + nE.push_back(new_edge);
  591 + }
  592 +
  593 + else {
  594 + vertex new_vertex = new_edge[1];
  595 +
  596 + nV[nV.size() - 1].e[0].push_back(nE.size());
  597 + new_vertex.e[1].push_back(nE.size());
  598 + new_edge.v[0] = nV.size() - 1;
  599 + new_edge.v[1] = nV.size();
  600 + nV.push_back(new_vertex);
  601 +
  602 + nE.push_back(new_edge);
  603 + }
  604 +
  605 + // get radii
  606 + nE[nE.size() - 1].cylinder<T>::set_r(0, E[i].cylinder<T>::r(j));
  607 + nE[nE.size() - 1].cylinder<T>::set_r(1, E[i].cylinder<T>::r(j + 1));
  608 + }
  609 + }
  610 + }
  611 +
  612 + (*this).E = nE;
  613 + (*this).V = nV;
  614 + }
  615 +
476 616 //load a network from an OBJ file
477 617 void load_obj(std::string filename){
478 618  
... ... @@ -715,7 +855,7 @@ public:
715 855 edge new_edge(C3); // new edge
716 856  
717 857 //create an edge from the given centerline
718   - unsigned int I = new_edge.size(); //calculate the number of points on the centerline
  858 + unsigned int I = (unsigned)new_edge.size(); //calculate the number of points on the centerline
719 859  
720 860 //get the first and last vertex IDs for the line
721 861 i[0] = S.E[l].front();
... ... @@ -1113,7 +1253,7 @@ public:
1113 1253 unsigned int id = 0; // split value
1114 1254 for(unsigned e = 0; e < E.size(); e++){ // for every edge
1115 1255 for(unsigned p = 0; p < E[e].size() - 1; p++){ // for every point in each edge
1116   - int t = E[e].length() / sigma * 2;
  1256 + int t = (int)(E[e].length() / sigma * 2);
1117 1257 if (t <= 20)
1118 1258 threshold_fac = E[e].size();
1119 1259 else
... ...
stim/grids/image_stack.h
... ... @@ -206,9 +206,12 @@ public:
206 206 }
207 207  
208 208  
  209 + /* This was causing compiler errors. I don't think this function call exists anywhere:
  210 +
209 211 void read(std::string file, unsigned int X, unsigned int Y, unsigned int Z, unsigned int C = 1, unsigned int header = 0){
210 212 read(file, stim::vec<unsigned long>(C, X, Y, Z), header);
211 213 }
  214 + */
212 215  
213 216 T* data(){
214 217 return ptr;
... ...
stim/image/image.h
... ... @@ -653,7 +653,7 @@ public:
653 653  
654 654 //crop regions given by an array of 1D index values
655 655 std::vector< image<T> > crop_idx(size_t w, size_t h, std::vector<size_t> idx) {
656   - std::vector<image<T>> result(idx.size()); //create an array of image files to return
  656 + std::vector< image<T> > result(idx.size()); //create an array of image files to return
657 657 for (size_t i = 0; i < idx.size(); i++) { //for each specified index point
658 658 size_t y = idx[i] / X(); //calculate the y coordinate from the 1D index (center of ROI)
659 659 size_t x = idx[i] - y * X(); //calculate the x coordinate (center of ROI)
... ...
stim/visualization/cylinder.h
... ... @@ -56,13 +56,12 @@ public:
56 56 init();
57 57 }
58 58  
59   -/*
60 59 cylinder(stim::centerline<T> p, std::vector<T> s)
61 60 {
62   - d = p;
  61 + //was d = s;
  62 + p = s;
63 63 init();
64 64 }
65   -*/
66 65  
67 66 //cylinder(centerline<T>c, T r) : centerline(c) {
68 67 // init();
... ...
stim/visualization/gl_network.h
1 1 #ifndef STIM_GL_NETWORK
2 2 #define STIM_GL_NETWORK
3 3  
  4 +#include <GL/glut.h>
4 5 #include <stim/biomodels/network.h>
5 6 #include <stim/visualization/aaboundingbox.h>
6 7  
... ...