1 // Released under the OSGPL license, as part of the OpenSceneGraph distribution. 2 // 3 // specification : http://www.trackvis.org/docs/?subsect=fileformat 4 // 5 #include <osg/Geode> 6 #include <osg/io_utils> 7 8 #include <osgDB/FileNameUtils> 9 #include <osgDB/FileUtils> 10 #include <osgDB/ReadFile> 11 #include <osgDB/Registry> 12 13 struct TrkHeader 14 { 15 char id_string[6]; 16 short dim[3]; 17 float voxel_size[3]; 18 float origin[3]; 19 short n_scalars; 20 char scalar_name[10][20]; 21 short n_properties; 22 char property_name[10][20]; 23 float vox_to_ras[4][4]; 24 char reserved[444]; 25 char voxel_order[4]; 26 char pad2[4]; 27 float image_orientation_patient[6]; 28 char pad1[2]; 29 unsigned char invert_x; 30 unsigned char invert_y; 31 unsigned char invert_z; 32 unsigned char swap_xy; 33 unsigned char swap_yz; 34 unsigned char swap_zx; 35 int n_count; 36 int version; 37 int hdr_size; 38 }; 39 40 const char vert_shader_str[] = 41 "void main(void)\n" 42 "{\n" 43 " vec4 eye = gl_ModelViewMatrixInverse * vec4(0.0,0.0,0.0,1.0);\n" 44 " vec3 rayVector = normalize(gl_Vertex.xyz-eye.xyz);\n" 45 "\n" 46 " vec3 dv = gl_Normal;\n" 47 " float d = dot(rayVector, dv);\n" 48 " float d2 = abs(d);//*d;\n" 49 " const float base=1.5;\n" 50 " float l = (base-d2)/base;\n" 51 " float half_l = l*0.5;\n" 52 "\n" 53 " // gl_FrontColor = vec4( (dv.x+1.0)*half_l, (dv.y+1.0)*half_l, (dv.z+1.0)*half_l, 1.0);\n" 54 " // gl_FrontColor = vec4( abs(dv.x)*half_l, abs(dv.y)*half_l, abs(dv.z)*half_l, 1.0);\n" 55 " gl_FrontColor = vec4( abs(dv.x), abs(dv.y), abs(dv.z), 1.0);\n" 56 " gl_Position = gl_ModelViewProjectionMatrix * gl_Vertex;\n" 57 " gl_ClipVertex = gl_ModelViewMatrix * gl_Vertex;\n" 58 "}\n"; 59 60 61 struct AssignDirectionColour 62 { AssignDirectionColourAssignDirectionColour63 AssignDirectionColour() {} 64 assignAssignDirectionColour65 void assign(osg::Geometry* geometry) 66 { 67 osg::Vec3Array* vertices = dynamic_cast<osg::Vec3Array*>(geometry->getVertexArray()); 68 if (!vertices) return; 69 70 // allocate colours 71 osg::ref_ptr<osg::Vec4Array> colours = dynamic_cast<osg::Vec4Array*>(geometry->getColorArray()); 72 if (!colours) 73 { 74 colours = new osg::Vec4Array; 75 geometry->setColorArray(colours.get()); 76 } 77 colours->setBinding(osg::Array::BIND_PER_VERTEX); 78 colours->resize(vertices->size(), osg::Vec4(0.0,0.0,0.0,0.0)); 79 #if 1 80 // allocate normals 81 osg::ref_ptr<osg::Vec3Array> normals = dynamic_cast<osg::Vec3Array*>(geometry->getNormalArray()); 82 if (!normals) 83 { 84 normals = new osg::Vec3Array; 85 geometry->setNormalArray(normals.get()); 86 } 87 normals->setBinding(osg::Array::BIND_PER_VERTEX); 88 normals->resize(vertices->size(), osg::Vec3(0.0,0.0,0.0)); 89 #endif 90 91 typedef std::vector<float> Divisors; 92 Divisors divisors; 93 divisors.resize(vertices->size(), 0.0f); 94 95 for(unsigned int i=0; i<geometry->getNumPrimitiveSets(); ++i) 96 { 97 osg::PrimitiveSet* pset = geometry->getPrimitiveSet(i); 98 if (pset->getMode()==GL_LINES) 99 { 100 for(unsigned int pi=0; pi<pset->getNumIndices()-1; pi+=2) 101 { 102 unsigned int vi_0 = pset->index(pi); 103 unsigned int vi_1 = pset->index(pi+1); 104 osg::Vec3& v0 = (*vertices)[vi_0]; 105 osg::Vec3& v1 = (*vertices)[vi_1]; 106 107 osg::Vec3 dv(v1-v0); 108 dv.normalize(); 109 110 // assign colours 111 osg::Vec4& c0 = (*colours)[vi_0]; 112 osg::Vec4& c1 = (*colours)[vi_1]; 113 114 osg::Vec4 colour( (dv.x()+1.0f/2.0f), (dv.y()+1.0f/2.0f), (dv.z()+1.0f/2.0f), 1.0f ); 115 116 if (divisors[vi_0]==0.0f) c0 = colour; 117 else c0 += colour; 118 119 if (divisors[vi_1]==0.0f) c1 = colour; 120 else c1 += colour; 121 122 // assign normals 123 osg::Vec3& n0 = (*normals)[vi_0]; 124 osg::Vec3& n1 = (*normals)[vi_1]; 125 126 #if 1 127 if (divisors[vi_0]==0.0f) n0 = dv; 128 else n0 += dv; 129 130 if (divisors[vi_1]==0.0f) n1 = dv; 131 else n1 += dv; 132 #else 133 float len_xy = sqrtf(dv.x()*dv.x() + dv.y()*dv.y()); 134 osg::Vec3 normal(0.0f, 0.0f, 0.0f); 135 if (len_xy!=0.0) 136 { 137 normal.set((-dv.x()/len_xy)*dv.z(), (-dv.y()/len_xy)*dv.z(), len_xy); 138 } 139 else 140 { 141 normal.set(1.0f, 0.0f, 0.0f); 142 } 143 144 if (divisors[vi_0]==0.0f) n0 = normal; 145 else n0 += normal; 146 147 if (divisors[vi_1]==0.0f) n1 = normal; 148 else n1 += normal; 149 #endif 150 151 divisors[vi_0] += 1.0f; 152 divisors[vi_1] += 1.0f; 153 } 154 } 155 } 156 157 for(unsigned int vi=0; vi<vertices->size(); ++vi) 158 { 159 if (divisors[vi]>0.0f) 160 { 161 (*colours)[vi] /= divisors[vi]; 162 (*normals)[vi].normalize(); 163 164 } 165 } 166 167 std::string vertexShaderFile("track.vert"); 168 169 osg::ref_ptr<osg::StateSet> stateset = geometry->getOrCreateStateSet(); 170 osg::ref_ptr<osg::Program> program = new osg::Program; 171 172 osg::ref_ptr<osg::Shader> vertexShader = osgDB::readRefShaderFile(osg::Shader::VERTEX, vertexShaderFile); 173 if (!vertexShader) 174 { 175 vertexShader = new osg::Shader(osg::Shader::VERTEX, vert_shader_str); 176 } 177 178 program->addShader(vertexShader.get()); 179 180 stateset->setAttribute(program.get()); 181 } 182 }; 183 184 185 class ReaderWriterTRK : public osgDB::ReaderWriter 186 { 187 public: 188 ReaderWriterTRK()189 ReaderWriterTRK() 190 { 191 supportsExtension("trk","Track file format"); 192 193 OSG_NOTICE<<"sizeof(TrkHeader)="<<sizeof(TrkHeader)<<std::endl; 194 } 195 className() const196 virtual const char* className() const { return "Track Reader/Writer"; } 197 198 readObject(std::istream & fin,const osgDB::ReaderWriter::Options * options=NULL) const199 virtual ReadResult readObject(std::istream& fin,const osgDB::ReaderWriter::Options* options =NULL) const 200 { 201 return readNode(fin, options); 202 } 203 readObject(const std::string & file,const osgDB::ReaderWriter::Options * options=NULL) const204 virtual ReadResult readObject(const std::string& file, const osgDB::ReaderWriter::Options* options =NULL) const 205 { 206 return readNode(file, options); 207 } 208 readNode(std::istream & fin,const osgDB::ReaderWriter::Options * =NULL) const209 virtual ReadResult readNode(std::istream& fin,const osgDB::ReaderWriter::Options* =NULL) const 210 { 211 TrkHeader header; 212 fin.read(reinterpret_cast<char*>(&header), sizeof(TrkHeader)); 213 214 if (fin.fail()) return ReadResult::ERROR_IN_READING_FILE; 215 216 OSG_NOTICE<<"Read header successfully ["<<header.id_string<<"]"<<std::endl; 217 bool requiresByteSwap = header.hdr_size!=1000; 218 if (requiresByteSwap) 219 { 220 OSG_NOTICE<<"Requires byteswap, but not implemented yet."<<std::endl; 221 return ReadResult::FILE_NOT_HANDLED; 222 } 223 else 224 { 225 OSG_NOTICE<<"No byteswap required"<<std::endl; 226 } 227 228 OSG_NOTICE<<"version = "<<header.version<<std::endl; 229 OSG_NOTICE<<"dim = "<<header.dim[0]<<", "<<header.dim[1]<<", "<<header.dim[2]<<std::endl; 230 OSG_NOTICE<<"voxel_size = "<<header.voxel_size[0]<<", "<<header.voxel_size[1]<<", "<<header.voxel_size[2]<<std::endl; 231 OSG_NOTICE<<"origin = "<<header.origin[0]<<", "<<header.origin[1]<<", "<<header.origin[2]<<std::endl; 232 OSG_NOTICE<<"n_scalars = "<<header.n_scalars<<std::endl; 233 OSG_NOTICE<<"n_properties = "<<header.n_properties<<std::endl; 234 235 osg::Matrixf matrix(reinterpret_cast<const float*>(header.vox_to_ras)); 236 OSG_NOTICE<<"vox_to_ras="<<matrix<<std::endl; 237 238 OSG_NOTICE<<"invert_x="<<static_cast<int>(header.invert_x)<<std::endl; 239 OSG_NOTICE<<"invert_y="<<static_cast<int>(header.invert_y)<<std::endl; 240 OSG_NOTICE<<"invert_z="<<static_cast<int>(header.invert_z)<<std::endl; 241 OSG_NOTICE<<"swap_xy="<<static_cast<int>(header.swap_xy)<<std::endl; 242 OSG_NOTICE<<"swap_yz="<<static_cast<int>(header.swap_yz)<<std::endl; 243 OSG_NOTICE<<"swap_zx="<<static_cast<int>(header.swap_zx)<<std::endl; 244 OSG_NOTICE<<"n_count="<<header.n_count<<std::endl; 245 OSG_NOTICE<<"hdr_size="<<header.hdr_size<<std::endl; 246 247 if (header.n_count==0) 248 { 249 OSG_NOTICE<<"No tracks stored."<<std::endl; 250 return ReadResult::FILE_NOT_HANDLED; 251 } 252 253 int n_s = header.n_scalars; 254 int n_p = header.n_properties; 255 256 osg::ref_ptr<osg::Geode> geode = new osg::Geode; 257 258 osg::ref_ptr<osg::Geometry> geometry = new osg::Geometry; 259 geode->addDrawable(geometry.get()); 260 #if 0 261 osg::ref_ptr<osg::StateSet> stateset = new osg::StateSet; 262 stateset->setMode(GL_LIGHTING, osg::StateAttribute::OFF); 263 geometry->setStateSet(stateset.get()); 264 #endif 265 osg::ref_ptr<osg::Vec4Array> colours = new osg::Vec4Array; 266 geometry->setColorArray(colours.get(), osg::Array::BIND_OVERALL); 267 colours->push_back(osg::Vec4(1.0f,1.0f,1.0f,1.0f)); 268 269 osg::ref_ptr<osg::Vec3Array> vertices = new osg::Vec3Array; 270 geometry->setVertexArray(vertices.get()); 271 272 osg::ref_ptr<osg::DrawElementsUInt> lines = new osg::DrawElementsUInt(GL_LINES); 273 geometry->addPrimitiveSet(lines.get()); 274 275 for(int i=0; i<header.n_count; ++i) 276 { 277 // read track 278 // OSG_NOTICE<<"reading track #"<<i<<std::endl; 279 280 int n_points=0; 281 fin.read(reinterpret_cast<char*>(&n_points),4); 282 if (fin.fail()) break; 283 284 // OSG_NOTICE<<" n_points="<<n_points<<std::endl; 285 286 int n_floats_per_vertex = (3+n_s); 287 int n_point_floats = n_floats_per_vertex*n_points; 288 float* point_data = new float[n_point_floats]; 289 fin.read(reinterpret_cast<char*>(point_data), n_point_floats*4); 290 if (fin.fail()) break; 291 292 if (n_p!=0) 293 { 294 float* property_data = new float[n_p]; 295 fin.read(reinterpret_cast<char*>(&property_data), n_p*4); 296 if (fin.fail()) break; 297 298 delete [] property_data; 299 300 } 301 302 if (n_points>0) 303 { 304 // record the index of the first vertex of the new track 305 int vi = vertices->size(); 306 307 // add the vertices of the track 308 for(int pi=0; pi<n_points; ++pi) 309 { 310 float* vptr = &point_data[pi*n_floats_per_vertex]; 311 vertices->push_back(osg::Vec3(vptr[0],vptr[1],vptr[2])); 312 } 313 314 // add the line segmenets for track 315 for(int pi=0; pi<n_points-1; ++pi, ++vi) 316 { 317 lines->push_back(vi); 318 lines->push_back(vi+1); 319 } 320 } 321 322 delete [] point_data; 323 } 324 325 if (fin.fail()) 326 { 327 OSG_NOTICE<<"Error on reading track file."<<std::endl; 328 return ReadResult::ERROR_IN_READING_FILE; 329 } 330 331 AssignDirectionColour colourer; 332 colourer.assign(geometry.get()); 333 334 return geode.release(); 335 } 336 readNode(const std::string & file,const osgDB::ReaderWriter::Options * options) const337 virtual ReadResult readNode(const std::string& file, const osgDB::ReaderWriter::Options* options) const 338 { 339 std::string ext = osgDB::getLowerCaseFileExtension(file); 340 if (!acceptsExtension(ext)) return ReadResult::FILE_NOT_HANDLED; 341 342 std::string fileName = osgDB::findDataFile( file, options ); 343 if (fileName.empty()) return ReadResult::FILE_NOT_FOUND; 344 345 osgDB::ifstream fin(fileName.c_str(), std::ios::in | std::ios::binary); 346 if(!fin) return ReadResult::FILE_NOT_HANDLED; 347 348 return readNode(fin, options); 349 } 350 }; 351 352 // now register with Registry to instantiate the above 353 // reader/writer. 354 REGISTER_OSGPLUGIN(trk, ReaderWriterTRK) 355