Compare View
Commits (7)
Showing
9 changed files
Show diff stats
python/classify.py
@@ -12,7 +12,7 @@ import sklearn.metrics | @@ -12,7 +12,7 @@ import sklearn.metrics | ||
12 | import scipy | 12 | import scipy |
13 | import scipy.misc | 13 | import scipy.misc |
14 | import envi | 14 | import envi |
15 | -import spectral | 15 | +import hyperspectral |
16 | import random | 16 | import random |
17 | import progressbar | 17 | import progressbar |
18 | import matplotlib.pyplot as plt | 18 | import matplotlib.pyplot as plt |
@@ -177,7 +177,7 @@ def envi_batch_predict(E, C, batch=10000): | @@ -177,7 +177,7 @@ def envi_batch_predict(E, C, batch=10000): | ||
177 | else: | 177 | else: |
178 | Tv = numpy.concatenate((Tv, C.predict(Fv.transpose()).transpose()), 0) | 178 | Tv = numpy.concatenate((Tv, C.predict(Fv.transpose()).transpose()), 0) |
179 | tempmask = E.batchmask() | 179 | tempmask = E.batchmask() |
180 | - Lv = spectral.unsift2(Tv, tempmask) | 180 | + Lv = hyperspectral.unsift2(Tv, tempmask) |
181 | Cv = label2class(Lv.squeeze(), background=0) | 181 | Cv = label2class(Lv.squeeze(), background=0) |
182 | RGB = class2color(Cv) | 182 | RGB = class2color(Cv) |
183 | plt.imshow(RGB) | 183 | plt.imshow(RGB) |
python/digitalstain.py
@@ -5,7 +5,7 @@ Created on Tue Jul 25 16:28:37 2017 | @@ -5,7 +5,7 @@ Created on Tue Jul 25 16:28:37 2017 | ||
5 | @author: david | 5 | @author: david |
6 | """ | 6 | """ |
7 | 7 | ||
8 | -import spectral | 8 | +import hyperspectral |
9 | import envi | 9 | import envi |
10 | import classify | 10 | import classify |
11 | import numpy | 11 | import numpy |
@@ -18,22 +18,36 @@ import glob | @@ -18,22 +18,36 @@ import glob | ||
18 | import matplotlib.pyplot as plt | 18 | import matplotlib.pyplot as plt |
19 | import random | 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 | mask = classify.random_mask(E.mask, N) | 28 | mask = classify.random_mask(E.mask, N) |
29 | + scipy.misc.imsave("random.bmp", mask) | ||
24 | 30 | ||
25 | Ft = E.loadmask(mask).transpose() | 31 | Ft = E.loadmask(mask).transpose() |
26 | 32 | ||
27 | stain = numpy.rollaxis(scipy.misc.imread(stainfile), 2) | 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 | CLASS = sklearn.neural_network.MLPRegressor(solver='lbfgs', alpha=1e-5, hidden_layer_sizes=(), random_state=1, verbose=True) | 37 | CLASS = sklearn.neural_network.MLPRegressor(solver='lbfgs', alpha=1e-5, hidden_layer_sizes=(), random_state=1, verbose=True) |
31 | CLASS.fit(Ft, Tt) | 38 | CLASS.fit(Ft, Tt) |
32 | 39 | ||
33 | if validate == False: | 40 | if validate == False: |
34 | return CLASS | 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 | Fv = E.loadbatch(batch_size) #load the first batch | 51 | Fv = E.loadbatch(batch_size) #load the first batch |
38 | n = 0 | 52 | n = 0 |
39 | while not Fv == []: #loop until an empty batch is returned | 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,7 +55,7 @@ def generate_stain(envifile, stainfile, N=5000, batch_size=10000, validate=True) | ||
41 | Tv = CLASS.predict(Fv.transpose()).transpose() | 55 | Tv = CLASS.predict(Fv.transpose()).transpose() |
42 | else: | 56 | else: |
43 | Tv = numpy.append(Tv, CLASS.predict(Fv.transpose()).transpose(), 1) #append the predicted labels from this batch to those of previous batches | 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 | RGB = numpy.rollaxis(COLORS, 0, 3).astype(numpy.ubyte) | 59 | RGB = numpy.rollaxis(COLORS, 0, 3).astype(numpy.ubyte) |
46 | plt.imshow(RGB) #display it | 60 | plt.imshow(RGB) #display it |
47 | plt.pause(0.05) | 61 | plt.pause(0.05) |
python/envi.py
@@ -100,11 +100,11 @@ class envi_header: | @@ -100,11 +100,11 @@ class envi_header: | ||
100 | #t = map(str.strip, t) #strip all of the tokens in the token list | 100 | #t = map(str.strip, t) #strip all of the tokens in the token list |
101 | 101 | ||
102 | #handle the simple conditions | 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 | self.samples = int(l[li].split()[-1]) | 108 | self.samples = int(l[li].split()[-1]) |
109 | elif l[li].startswith("lines"): | 109 | elif l[li].startswith("lines"): |
110 | self.lines = int(l[li].split()[-1]) | 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,6 +473,146 @@ public: | ||
473 | return E[e]; //return the specified edge (casting it to a fiber) | 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 | //load a network from an OBJ file | 616 | //load a network from an OBJ file |
477 | void load_obj(std::string filename){ | 617 | void load_obj(std::string filename){ |
478 | 618 | ||
@@ -715,7 +855,7 @@ public: | @@ -715,7 +855,7 @@ public: | ||
715 | edge new_edge(C3); // new edge | 855 | edge new_edge(C3); // new edge |
716 | 856 | ||
717 | //create an edge from the given centerline | 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 | //get the first and last vertex IDs for the line | 860 | //get the first and last vertex IDs for the line |
721 | i[0] = S.E[l].front(); | 861 | i[0] = S.E[l].front(); |
@@ -1113,7 +1253,7 @@ public: | @@ -1113,7 +1253,7 @@ public: | ||
1113 | unsigned int id = 0; // split value | 1253 | unsigned int id = 0; // split value |
1114 | for(unsigned e = 0; e < E.size(); e++){ // for every edge | 1254 | for(unsigned e = 0; e < E.size(); e++){ // for every edge |
1115 | for(unsigned p = 0; p < E[e].size() - 1; p++){ // for every point in each edge | 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 | if (t <= 20) | 1257 | if (t <= 20) |
1118 | threshold_fac = E[e].size(); | 1258 | threshold_fac = E[e].size(); |
1119 | else | 1259 | else |
stim/grids/image_stack.h
@@ -206,9 +206,12 @@ public: | @@ -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 | void read(std::string file, unsigned int X, unsigned int Y, unsigned int Z, unsigned int C = 1, unsigned int header = 0){ | 211 | void read(std::string file, unsigned int X, unsigned int Y, unsigned int Z, unsigned int C = 1, unsigned int header = 0){ |
210 | read(file, stim::vec<unsigned long>(C, X, Y, Z), header); | 212 | read(file, stim::vec<unsigned long>(C, X, Y, Z), header); |
211 | } | 213 | } |
214 | + */ | ||
212 | 215 | ||
213 | T* data(){ | 216 | T* data(){ |
214 | return ptr; | 217 | return ptr; |
stim/image/image.h
@@ -653,7 +653,7 @@ public: | @@ -653,7 +653,7 @@ public: | ||
653 | 653 | ||
654 | //crop regions given by an array of 1D index values | 654 | //crop regions given by an array of 1D index values |
655 | std::vector< image<T> > crop_idx(size_t w, size_t h, std::vector<size_t> idx) { | 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 | for (size_t i = 0; i < idx.size(); i++) { //for each specified index point | 657 | for (size_t i = 0; i < idx.size(); i++) { //for each specified index point |
658 | size_t y = idx[i] / X(); //calculate the y coordinate from the 1D index (center of ROI) | 658 | size_t y = idx[i] / X(); //calculate the y coordinate from the 1D index (center of ROI) |
659 | size_t x = idx[i] - y * X(); //calculate the x coordinate (center of ROI) | 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,13 +56,12 @@ public: | ||
56 | init(); | 56 | init(); |
57 | } | 57 | } |
58 | 58 | ||
59 | -/* | ||
60 | cylinder(stim::centerline<T> p, std::vector<T> s) | 59 | cylinder(stim::centerline<T> p, std::vector<T> s) |
61 | { | 60 | { |
62 | - d = p; | 61 | + //was d = s; |
62 | + p = s; | ||
63 | init(); | 63 | init(); |
64 | } | 64 | } |
65 | -*/ | ||
66 | 65 | ||
67 | //cylinder(centerline<T>c, T r) : centerline(c) { | 66 | //cylinder(centerline<T>c, T r) : centerline(c) { |
68 | // init(); | 67 | // init(); |
stim/visualization/gl_network.h
1 | #ifndef STIM_GL_NETWORK | 1 | #ifndef STIM_GL_NETWORK |
2 | #define STIM_GL_NETWORK | 2 | #define STIM_GL_NETWORK |
3 | 3 | ||
4 | +#include <GL/glut.h> | ||
4 | #include <stim/biomodels/network.h> | 5 | #include <stim/biomodels/network.h> |
5 | #include <stim/visualization/aaboundingbox.h> | 6 | #include <stim/visualization/aaboundingbox.h> |
6 | 7 |