diff --git a/Modules/DiffusionImaging/MiniApps/NetworkStatistics.cpp b/Modules/DiffusionImaging/MiniApps/NetworkStatistics.cpp index 187e3de41e..db1b4a7ef9 100644 --- a/Modules/DiffusionImaging/MiniApps/NetworkStatistics.cpp +++ b/Modules/DiffusionImaging/MiniApps/NetworkStatistics.cpp @@ -1,337 +1,509 @@ /*=================================================================== The Medical Imaging Interaction Toolkit (MITK) Copyright (c) German Cancer Research Center, Division of Medical and Biological Informatics. All rights reserved. This software is distributed WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See LICENSE.txt or http://www.mitk.org for details. ===================================================================*/ -#include "MiniAppManager.h" +#include // std includes #include #include #include #include #include +#include + +// boost includes +#include // ITK includes #include // CTK includes #include "ctkCommandLineParser.h" // MITK includes #include #include +#include #include int NetworkStatistics(int argc, char* argv[]) { ctkCommandLineParser parser; parser.setArgumentPrefix("--", "-"); parser.addArgument("inputNetwork", "i", ctkCommandLineParser::String, "input connectomics network (.cnf)", us::Any(), false); parser.addArgument("outputFile", "o", ctkCommandLineParser::String, "name of output file", us::Any(), false); parser.addArgument("noGlobalStatistics", "g", ctkCommandLineParser::Bool, "Do not calculate global statistics"); parser.addArgument("createConnectivityMatriximage", "I", ctkCommandLineParser::Bool, "Write connectivity matrix image"); parser.addArgument("binaryConnectivity", "b", ctkCommandLineParser::Bool, "Whether to create a binary connectivity matrix"); parser.addArgument("rescaleConnectivity", "r", ctkCommandLineParser::Bool, "Whether to rescale the connectivity matrix"); parser.addArgument("localStatistics", "L", ctkCommandLineParser::StringList, "Provide a list of node labels for local statistics", us::Any()); - + parser.addArgument("regionList", "R", ctkCommandLineParser::StringList, "A space separated list of regions. Each region has the format\n regionname;label1;label2;...;labelN", us::Any()); + parser.addArgument("granularity", "gr", ctkCommandLineParser::Int, "How finely to test the density range and how many thresholds to consider"); + parser.addArgument("startDensity", "d", ctkCommandLineParser::Bool, "Largest density for the range"); + parser.addArgument("thresholdStepSize", "t", ctkCommandLineParser::Int, "Distance of two adjacent thresholds"); map parsedArgs = parser.parseArguments(argc, argv); if (parsedArgs.size()==0) return EXIT_FAILURE; //default values bool noGlobalStatistics( false ); bool binaryConnectivity( false ); bool rescaleConnectivity( false ); bool createConnectivityMatriximage( false ); + unsigned int granularity( 1 ); + double startDensity( 1.0 ); + unsigned int thresholdStepSize( 3 ); // parse command line arguments std::string networkName = us::any_cast(parsedArgs["inputNetwork"]); std::string outName = us::any_cast(parsedArgs["outputFile"]); ctkCommandLineParser::StringContainerType localLabels; if(parsedArgs.count("localStatistics")) { localLabels = us::any_cast(parsedArgs["localStatistics"]); } + ctkCommandLineParser::StringContainerType unparsedRegions; + std::map< std::string, std::vector > parsedRegions; + std::map< std::string, std::vector >::iterator parsedRegionsIterator; + + if(parsedArgs.count("regionList")) + { + unparsedRegions = us::any_cast(parsedArgs["regionList"]); + + for(unsigned int index(0); index < unparsedRegions.size(); index++ ) + { + std::vector< std::string > tempRegionVector; + + boost::split(tempRegionVector, unparsedRegions.at(index), boost::is_any_of(";")); + + std::vector< std::string >::const_iterator begin = tempRegionVector.begin(); + std::vector< std::string >::const_iterator last = tempRegionVector.begin() + tempRegionVector.size(); + std::vector< std::string > insertRegionVector(begin + 1, last); + + if( parsedRegions.count( tempRegionVector.at(0) ) == 0 ) + { + parsedRegions.insert( std::pair< std::string, std::vector >( tempRegionVector.at(0), insertRegionVector) ); + } + else + { + MITK_ERROR << "Region already exists. Skipping second occurrence."; + } + } + } + if (parsedArgs.count("noGlobalStatistics")) noGlobalStatistics = us::any_cast(parsedArgs["noGlobalStatistics"]); if (parsedArgs.count("binaryConnectivity")) binaryConnectivity = us::any_cast(parsedArgs["binaryConnectivity"]); if (parsedArgs.count("rescaleConnectivity")) rescaleConnectivity = us::any_cast(parsedArgs["rescaleConnectivity"]); if (parsedArgs.count("createConnectivityMatriximage")) createConnectivityMatriximage = us::any_cast(parsedArgs["createConnectivityMatriximage"]); + if (parsedArgs.count("granularity")) + granularity = us::any_cast(parsedArgs["granularity"]); + if (parsedArgs.count("startDensity")) + startDensity = us::any_cast(parsedArgs["startDensity"]); + if (parsedArgs.count("thresholdStepSize")) + thresholdStepSize = us::any_cast(parsedArgs["thresholdStepSize"]); try { const std::string s1="", s2=""; // load network std::vector networkFile = mitk::BaseDataIO::LoadBaseDataFromFile( networkName, s1, s2, false ); if( networkFile.empty() ) { std::string errorMessage = "File at " + networkName + " could not be read. Aborting."; MITK_ERROR << errorMessage; return EXIT_FAILURE; } mitk::BaseData* networkBaseData = networkFile.at(0); mitk::ConnectomicsNetwork* network = dynamic_cast( networkBaseData ); if( !network ) { std::string errorMessage = "Read file at " + networkName + " could not be recognized as network. Aborting."; MITK_ERROR << errorMessage; return EXIT_FAILURE; } - mitk::ConnectomicsStatisticsCalculator::Pointer statisticsCalculator = mitk::ConnectomicsStatisticsCalculator::New(); - - statisticsCalculator->SetNetwork( network ); - statisticsCalculator->Update(); - - // global statistics - if( !noGlobalStatistics ) + // streams + std::stringstream globalHeaderStream; + globalHeaderStream << "NumberOfVertices " + << "NumberOfEdges " + << "AverageDegree " + << "ConnectionDensity " + << "NumberOfConnectedComponents " + << "AverageComponentSize " + << "LargestComponentSize " + << "RatioOfNodesInLargestComponent " + << "HopPlotExponent " + << "EffectiveHopDiameter " + << "AverageClusteringCoefficientsC " + << "AverageClusteringCoefficientsD " + << "AverageClusteringCoefficientsE " + << "AverageVertexBetweennessCentrality " + << "AverageEdgeBetweennessCentrality " + << "NumberOfIsolatedPoints " + << "RatioOfIsolatedPoints " + << "NumberOfEndPoints " + << "RatioOfEndPoints " + << "Diameter " + << "Diameter90 " + << "Radius " + << "Radius90 " + << "AverageEccentricity " + << "AverageEccentricity90 " + << "AveragePathLength " + << "NumberOfCentralPoints " + << "RatioOfCentralPoints " + << "SpectralRadius " + << "SecondLargestEigenValue " + << "AdjacencyTrace " + << "AdjacencyEnergy " + << "LaplacianTrace " + << "LaplacianEnergy " + << "LaplacianSpectralGap " + << "NormalizedLaplacianTrace " + << "NormalizedLaplacianEnergy " + << "NormalizedLaplacianNumberOf2s " + << "NormalizedLaplacianNumberOf1s " + << "NormalizedLaplacianNumberOf0s " + << "NormalizedLaplacianLowerSlope " + << "NormalizedLaplacianUpperSlope" + << "SmallWorldness" + << std::endl; + + std::stringstream localHeaderStream; + std::stringstream regionalHeaderStream; + + std::stringstream globalDataStream; + std::stringstream localDataStream; + std::stringstream regionalDataStream; + + std::string globalOutName = outName + "_global.txt"; + std::string localOutName = outName + "_local.txt"; + std::string regionalOutName = outName + "_regional.txt"; + + // iterate over all three possible methods + for(unsigned int method( 1 ); method < 3; method++) { - std::string globalOutName = outName + "_global.txt"; - std::stringstream headerStream; - headerStream << "NumberOfVertices " - << "NumberOfEdges " - << "AverageDegree " - << "ConnectionDensity " - << "NumberOfConnectedComponents " - << "AverageComponentSize " - << "LargestComponentSize " - << "RatioOfNodesInLargestComponent " - << "HopPlotExponent " - << "EffectiveHopDiameter " - << "AverageClusteringCoefficientsC " - << "AverageClusteringCoefficientsD " - << "AverageClusteringCoefficientsE " - << "AverageVertexBetweennessCentrality " - << "AverageEdgeBetweennessCentrality " - << "NumberOfIsolatedPoints " - << "RatioOfIsolatedPoints " - << "NumberOfEndPoints " - << "RatioOfEndPoints " - << "Diameter " - << "Diameter90 " - << "Radius " - << "Radius90 " - << "AverageEccentricity " - << "AverageEccentricity90 " - << "AveragePathLength " - << "NumberOfCentralPoints " - << "RatioOfCentralPoints " - << "SpectralRadius " - << "SecondLargestEigenValue " - << "AdjacencyTrace " - << "AdjacencyEnergy " - << "LaplacianTrace " - << "LaplacianEnergy " - << "LaplacianSpectralGap " - << "NormalizedLaplacianTrace " - << "NormalizedLaplacianEnergy " - << "NormalizedLaplacianNumberOf2s " - << "NormalizedLaplacianNumberOf1s " - << "NormalizedLaplacianNumberOf0s " - << "NormalizedLaplacianLowerSlope " - << "NormalizedLaplacianUpperSlope" - << "SmallWorldness" - << std::endl; - - std::stringstream dataStream; - dataStream << statisticsCalculator->GetNumberOfVertices() << " " - << statisticsCalculator->GetNumberOfEdges() << " " - << statisticsCalculator->GetAverageDegree() << " " - << statisticsCalculator->GetConnectionDensity() << " " - << statisticsCalculator->GetNumberOfConnectedComponents() << " " - << statisticsCalculator->GetAverageComponentSize() << " " - << statisticsCalculator->GetLargestComponentSize() << " " - << statisticsCalculator->GetRatioOfNodesInLargestComponent() << " " - << statisticsCalculator->GetHopPlotExponent() << " " - << statisticsCalculator->GetEffectiveHopDiameter() << " " - << statisticsCalculator->GetAverageClusteringCoefficientsC() << " " - << statisticsCalculator->GetAverageClusteringCoefficientsD() << " " - << statisticsCalculator->GetAverageClusteringCoefficientsE() << " " - << statisticsCalculator->GetAverageVertexBetweennessCentrality() << " " - << statisticsCalculator->GetAverageEdgeBetweennessCentrality() << " " - << statisticsCalculator->GetNumberOfIsolatedPoints() << " " - << statisticsCalculator->GetRatioOfIsolatedPoints() << " " - << statisticsCalculator->GetNumberOfEndPoints() << " " - << statisticsCalculator->GetRatioOfEndPoints() << " " - << statisticsCalculator->GetDiameter() << " " - << statisticsCalculator->GetDiameter90() << " " - << statisticsCalculator->GetRadius() << " " - << statisticsCalculator->GetRadius90() << " " - << statisticsCalculator->GetAverageEccentricity() << " " - << statisticsCalculator->GetAverageEccentricity90() << " " - << statisticsCalculator->GetAveragePathLength() << " " - << statisticsCalculator->GetNumberOfCentralPoints() << " " - << statisticsCalculator->GetRatioOfCentralPoints() << " " - << statisticsCalculator->GetSpectralRadius() << " " - << statisticsCalculator->GetSecondLargestEigenValue() << " " - << statisticsCalculator->GetAdjacencyTrace() << " " - << statisticsCalculator->GetAdjacencyEnergy() << " " - << statisticsCalculator->GetLaplacianTrace() << " " - << statisticsCalculator->GetLaplacianEnergy() << " " - << statisticsCalculator->GetLaplacianSpectralGap() << " " - << statisticsCalculator->GetNormalizedLaplacianTrace() << " " - << statisticsCalculator->GetNormalizedLaplacianEnergy() << " " - << statisticsCalculator->GetNormalizedLaplacianNumberOf2s() << " " - << statisticsCalculator->GetNormalizedLaplacianNumberOf1s() << " " - << statisticsCalculator->GetNormalizedLaplacianNumberOf0s() << " " - << statisticsCalculator->GetNormalizedLaplacianLowerSlope() << " " - << statisticsCalculator->GetNormalizedLaplacianUpperSlope() << " " - << statisticsCalculator->GetSmallWorldness() - << std::endl; - - std::ofstream outFile( globalOutName.c_str(), ios::out ); - - if( ! outFile.is_open() ) - { - std::string errorMessage = "Could not open " + globalOutName + " for writing."; - MITK_ERROR << errorMessage; - return EXIT_FAILURE; - } - - outFile << headerStream.str() << dataStream.str(); - outFile.close(); - } // end global statistics + // 0 - Random removal threshold + // 1 - Largest density below threshold + // 2 - Threshold based - //create connectivity matrix png - if( createConnectivityMatriximage ) - { - std::string connectivity_png_postfix = "_connectivity"; - if( binaryConnectivity ) - { - connectivity_png_postfix += "_binary"; - } - else if( rescaleConnectivity ) + // iterate over possible targets + for( unsigned int step( 0 ); step < granularity; step++ ) { - connectivity_png_postfix += "_rescaled"; - } - connectivity_png_postfix += ".png"; - - /* File format - * A png file depicting the binary connectivity matrix - */ - itk::ConnectomicsNetworkToConnectivityMatrixImageFilter::Pointer filter = itk::ConnectomicsNetworkToConnectivityMatrixImageFilter::New(); - filter->SetInputNetwork( network ); - filter->SetBinaryConnectivity( binaryConnectivity ); - filter->SetRescaleConnectivity( rescaleConnectivity ); - filter->Update(); - - typedef itk::ConnectomicsNetworkToConnectivityMatrixImageFilter::OutputImageType connectivityMatrixImageType; - - itk::ImageFileWriter< connectivityMatrixImageType >::Pointer connectivityWriter = itk::ImageFileWriter< connectivityMatrixImageType >::New(); - - connectivityWriter->SetInput( filter->GetOutput() ); - connectivityWriter->SetFileName( outName + connectivity_png_postfix); - connectivityWriter->Update(); + double targetValue( 0.0 ); + bool newStep( true ); - MITK_INFO << "Connectivity matrix image written."; - } // end create connectivity matrix png + switch ( method ) + { + case mitk::ConnectomicsNetworkThresholder::RandomRemovalOfWeakest : + case mitk::ConnectomicsNetworkThresholder::LargestLowerThanDensity : + targetValue = startDensity * (1 - static_cast( step ) / ( granularity + 0.5 ) ); + break; + case mitk::ConnectomicsNetworkThresholder::ThresholdBased : + targetValue = static_cast( thresholdStepSize * step ); + break; + default: + MITK_ERROR << "Invalid thresholding method called, aborting."; + return EXIT_FAILURE; + break; + } - // calculate local averages - if( localLabels.size() > 0 ) - { - std::string localOutName = outName + "_local.txt"; + mitk::ConnectomicsNetworkThresholder::Pointer thresholder = mitk::ConnectomicsNetworkThresholder::New(); + thresholder->SetNetwork( network ); + thresholder->SetTargetThreshold( targetValue ); + thresholder->SetTargetDensity( targetValue ); + // TEST + thresholder->SetThresholdingScheme( static_cast(1) ); + // ENDTEST + //thresholder->SetThresholdingScheme( static_cast(method) ); + mitk::ConnectomicsNetwork::Pointer thresholdedNetwork = thresholder->GetThresholdedNetwork(); + + mitk::ConnectomicsStatisticsCalculator::Pointer statisticsCalculator = mitk::ConnectomicsStatisticsCalculator::New(); + statisticsCalculator->SetNetwork( thresholdedNetwork ); + statisticsCalculator->Update(); + + // global statistics + if( !noGlobalStatistics ) + { + globalDataStream << statisticsCalculator->GetNumberOfVertices() << " " + << statisticsCalculator->GetNumberOfEdges() << " " + << statisticsCalculator->GetAverageDegree() << " " + << statisticsCalculator->GetConnectionDensity() << " " + << statisticsCalculator->GetNumberOfConnectedComponents() << " " + << statisticsCalculator->GetAverageComponentSize() << " " + << statisticsCalculator->GetLargestComponentSize() << " " + << statisticsCalculator->GetRatioOfNodesInLargestComponent() << " " + << statisticsCalculator->GetHopPlotExponent() << " " + << statisticsCalculator->GetEffectiveHopDiameter() << " " + << statisticsCalculator->GetAverageClusteringCoefficientsC() << " " + << statisticsCalculator->GetAverageClusteringCoefficientsD() << " " + << statisticsCalculator->GetAverageClusteringCoefficientsE() << " " + << statisticsCalculator->GetAverageVertexBetweennessCentrality() << " " + << statisticsCalculator->GetAverageEdgeBetweennessCentrality() << " " + << statisticsCalculator->GetNumberOfIsolatedPoints() << " " + << statisticsCalculator->GetRatioOfIsolatedPoints() << " " + << statisticsCalculator->GetNumberOfEndPoints() << " " + << statisticsCalculator->GetRatioOfEndPoints() << " " + << statisticsCalculator->GetDiameter() << " " + << statisticsCalculator->GetDiameter90() << " " + << statisticsCalculator->GetRadius() << " " + << statisticsCalculator->GetRadius90() << " " + << statisticsCalculator->GetAverageEccentricity() << " " + << statisticsCalculator->GetAverageEccentricity90() << " " + << statisticsCalculator->GetAveragePathLength() << " " + << statisticsCalculator->GetNumberOfCentralPoints() << " " + << statisticsCalculator->GetRatioOfCentralPoints() << " " + << statisticsCalculator->GetSpectralRadius() << " " + << statisticsCalculator->GetSecondLargestEigenValue() << " " + << statisticsCalculator->GetAdjacencyTrace() << " " + << statisticsCalculator->GetAdjacencyEnergy() << " " + << statisticsCalculator->GetLaplacianTrace() << " " + << statisticsCalculator->GetLaplacianEnergy() << " " + << statisticsCalculator->GetLaplacianSpectralGap() << " " + << statisticsCalculator->GetNormalizedLaplacianTrace() << " " + << statisticsCalculator->GetNormalizedLaplacianEnergy() << " " + << statisticsCalculator->GetNormalizedLaplacianNumberOf2s() << " " + << statisticsCalculator->GetNormalizedLaplacianNumberOf1s() << " " + << statisticsCalculator->GetNormalizedLaplacianNumberOf0s() << " " + << statisticsCalculator->GetNormalizedLaplacianLowerSlope() << " " + << statisticsCalculator->GetNormalizedLaplacianUpperSlope() << " " + << statisticsCalculator->GetSmallWorldness() + << std::endl; + + } // end global statistics + + //create connectivity matrix png + if( createConnectivityMatriximage ) + { + std::string connectivity_png_postfix = "_connectivity"; + if( binaryConnectivity ) + { + connectivity_png_postfix += "_binary"; + } + else if( rescaleConnectivity ) + { + connectivity_png_postfix += "_rescaled"; + } + connectivity_png_postfix += ".png"; + + /* File format + * A png file depicting the binary connectivity matrix + */ + itk::ConnectomicsNetworkToConnectivityMatrixImageFilter::Pointer filter = itk::ConnectomicsNetworkToConnectivityMatrixImageFilter::New(); + filter->SetInputNetwork( network ); + filter->SetBinaryConnectivity( binaryConnectivity ); + filter->SetRescaleConnectivity( rescaleConnectivity ); + filter->Update(); + + typedef itk::ConnectomicsNetworkToConnectivityMatrixImageFilter::OutputImageType connectivityMatrixImageType; + + itk::ImageFileWriter< connectivityMatrixImageType >::Pointer connectivityWriter = itk::ImageFileWriter< connectivityMatrixImageType >::New(); + + connectivityWriter->SetInput( filter->GetOutput() ); + connectivityWriter->SetFileName( outName + connectivity_png_postfix); + connectivityWriter->Update(); + + MITK_INFO << "Connectivity matrix image written."; + } // end create connectivity matrix png + + /* + * We can either calculate local indices for specific nodes, or specific regions + */ + + // Create LabelToIndex translation + std::map< std::string, int > labelToIdMap; + std::vector< mitk::ConnectomicsNetwork::NetworkNode > nodeVector = thresholdedNetwork->GetVectorOfAllNodes(); + for(int loop(0); loop < nodeVector.size(); loop++) + { + labelToIdMap.insert( std::pair< std::string, int>(nodeVector.at(loop).label, nodeVector.at(loop).id) ); + } - // Create LabelToIndex translation - std::map< std::string, int > labelToIdMap; - std::vector< mitk::ConnectomicsNetwork::NetworkNode > nodeVector = network->GetVectorOfAllNodes(); - for(unsigned int loop(0); loop < nodeVector.size(); loop++) - { - labelToIdMap.insert( std::pair< std::string, int>(nodeVector.at(loop).label, nodeVector.at(loop).id) ); - } + std::vector< int > degreeVector = thresholdedNetwork->GetDegreeOfNodes(); + std::vector< double > ccVector = thresholdedNetwork->GetLocalClusteringCoefficients( ); + std::vector< double > bcVector = thresholdedNetwork->GetNodeBetweennessVector( ); - std::vector< int > degreeVector = network->GetDegreeOfNodes(); - std::vector< double > ccVector = network->GetLocalClusteringCoefficients( ); - std::vector< double > bcVector = network->GetNodeBetweennessVector( ); - double sumDegree( 0 ); - double sumCC( 0 ); - double sumBC( 0 ); - double count( 0 ); + // calculate local indices - for(unsigned int loop(0); loop < localLabels.size(); loop++ ) - { - if( network->CheckForLabel(localLabels.at( loop )) ) { - sumDegree = sumDegree + degreeVector.at( labelToIdMap.find( localLabels.at( loop ) )->second ); - sumCC = sumCC + ccVector.at( labelToIdMap.find( localLabels.at( loop ) )->second ); - sumBC = sumBC + bcVector.at( labelToIdMap.find( localLabels.at( loop ) )->second ); - count = count + 1; + // only add to header for the first step of the first method + if(step == 0 && method == 0) + { + localHeaderStream << "Th_method " << "Th_target " << "density"; + } + + double density = statisticsCalculator->GetConnectionDensity(); + + localDataStream << "\n" + << method << " " + << targetValue << " " + << density; + + for(unsigned int loop(0); loop < localLabels.size(); loop++ ) + { + if( network->CheckForLabel(localLabels.at( loop )) ) + { + localHeaderStream << localLabels.at( loop ) << "_Degree " + << localLabels.at( loop ) << "_CC " + << localLabels.at( loop ) << "_BC "; + + localDataStream << degreeVector.at( labelToIdMap.find( localLabels.at( loop ) )->second ) << " " + << ccVector.at( labelToIdMap.find( localLabels.at( loop ) )->second ) << " " + << bcVector.at( labelToIdMap.find( localLabels.at( loop ) )->second ) << " "; + } + else + { + MITK_ERROR << "Illegal label. Label: \"" << localLabels.at( loop ) << "\" not found."; + } + } } - else + + // calculate regional indices + { - MITK_ERROR << "Illegal label. Label: \"" << localLabels.at( loop ) << "\" not found."; + // only add to header for the first step of the first method + if(step == 0 && method == 0) + { + regionalHeaderStream << "Th_method " << "Th_target " << "density"; + } + + double density = statisticsCalculator->GetConnectionDensity(); + + regionalDataStream << "\n" + << method << " " + << targetValue << " " + << density; + + for( parsedRegionsIterator = parsedRegions.begin(); parsedRegionsIterator != parsedRegions.end(); parsedRegionsIterator++ ) + { + std::vector regionLabelsVector = parsedRegionsIterator->second; + std::string regionName = parsedRegionsIterator->first; + + double sumDegree( 0 ); + double sumCC( 0 ); + double sumBC( 0 ); + double count( 0 ); + + for( int loop(0); loop < regionLabelsVector.size(); loop++ ) + { + if( thresholdedNetwork->CheckForLabel(regionLabelsVector.at( loop )) ) + { + sumDegree = sumDegree + degreeVector.at( labelToIdMap.find( regionLabelsVector.at( loop ) )->second ); + sumCC = sumCC + ccVector.at( labelToIdMap.find( regionLabelsVector.at( loop ) )->second ); + sumBC = sumBC + bcVector.at( labelToIdMap.find( regionLabelsVector.at( loop ) )->second ); + count = count + 1; + } + else + { + MITK_ERROR << "Illegal label. Label: \"" << regionLabelsVector.at( loop ) << "\" not found."; + } + } + + // only add to header for the first step of the first method + if(step == 0 && method == 0) + { + regionalHeaderStream << " " << regionName << "_LocalAverageDegree " + << regionName << "_LocalAverageCC " + << regionName << "_LocalAverageBC " + << regionName << "_NumberOfNodes"; + } + + regionalDataStream << " " << sumDegree / count << " " + << sumCC / count << " " + << sumBC / count << " " + << count; + } } } - std::stringstream headerStream; - headerStream << "LocalAverageDegree " - << "LocalAverageCC " - << "LocalAverageBC " - << "NumberOfNodes" - << std::endl; - - std::stringstream dataStream; - dataStream << sumDegree / count << " " - << sumCC / count << " " - << sumBC / count << " " - << count - << std::endl; + }// end calculate local averages - std::ofstream outFile( localOutName.c_str(), ios::out ); + if( !noGlobalStatistics ) + { + MITK_INFO << "Writing to " << globalOutName; + std::ofstream glocalOutFile( globalOutName.c_str(), ios::out ); + if( ! glocalOutFile.is_open() ) + { + std::string errorMessage = "Could not open " + globalOutName + " for writing."; + MITK_ERROR << errorMessage; + return EXIT_FAILURE; + } + glocalOutFile << globalHeaderStream.str() << globalDataStream.str(); + glocalOutFile.close(); + } - if( ! outFile.is_open() ) + if( localLabels.size() > 0 ) + { + MITK_INFO << "Writing to " << localOutName; + std::ofstream localOutFile( localOutName.c_str(), ios::out ); + if( ! localOutFile.is_open() ) { std::string errorMessage = "Could not open " + localOutName + " for writing."; MITK_ERROR << errorMessage; return EXIT_FAILURE; } + localOutFile << localHeaderStream.str() << localDataStream.str(); + localOutFile.close(); + } - outFile << headerStream.str() << dataStream.str(); - outFile.close(); - }// end calculate local averages + if( parsedRegions.size() > 0 ) + { + MITK_INFO << "Writing to " << regionalOutName; + std::ofstream regionalOutFile( regionalOutName.c_str(), ios::out ); + if( ! regionalOutFile.is_open() ) + { + std::string errorMessage = "Could not open " + regionalOutName + " for writing."; + MITK_ERROR << errorMessage; + return EXIT_FAILURE; + } + regionalOutFile << regionalHeaderStream.str() << regionalDataStream.str(); + regionalOutFile.close(); + } return EXIT_SUCCESS; } catch (itk::ExceptionObject e) { MITK_INFO << e; return EXIT_FAILURE; } catch (std::exception e) { MITK_INFO << e.what(); return EXIT_FAILURE; } catch (...) { MITK_INFO << "ERROR!?!"; return EXIT_FAILURE; } MITK_INFO << "DONE"; return EXIT_SUCCESS; } RegisterDiffusionMiniApp(NetworkStatistics);