Changeset 25615


Ignore:
Timestamp:
09/29/20 18:21:19 (4 years ago)
Author:
Eric.Larour
Message:

CHG: diverse

Location:
issm/branches/trunk-larour-SLPS2020/src/c
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • issm/branches/trunk-larour-SLPS2020/src/c/main/issm_dakota.cpp

    r25605 r25615  
    44
    55#include "./issm.h"
    6 #include <sys/stat.h>
    76
    87/*Dakota includes: */
    … …  
    1413#include "DakotaInterface.hpp"
    1514#endif
    16 
    17 /*prototypes:*/
    18 bool dirstructure(int argc,char** argv);
    19 int issm_dakota_statistics(int argc,char** argv);
    2015
    2116int main(int argc,char **argv){ /*{{{*/
    … …  
    4641
    4742        /*Create directory structure for model outputs:*/
    48         statistics=dirstructure(argc,argv);
     43        statistics=DakotaDirStructure(argc,argv);
    4944
    5045        /* Parse input and construct Dakota LibraryEnvironment, performing input data checks*/
    … …  
    9388
    9489        /* Run statistics if requested:*/
    95         if(statistics)issm_dakota_statistics(argc,argv);
     90        if(statistics)DakotaStatistics(argc,argv);
    9691
    9792        /*free allocations:*/
    … …  
    108103
    109104} /*}}}*/
    110 bool dirstructure(int argc,char** argv){ /*{{{*/
    111 
    112         char* input_file;
    113         FILE* fid;
    114         IoModel* iomodel=NULL;
    115         int check;
    116 
    117         //qmu statistics
    118         bool statistics    = false;
    119         int  numdirectories = 0;
    120 
    121         /*First things first, set the communicator as a global variable: */
    122         IssmComm::SetComm(MPI_COMM_WORLD);
    123 
    124         /*Barrier:*/
    125         ISSM_MPI_Barrier(IssmComm::GetComm());
    126         _printf0_("Preparing directory structure for model outputs:" << "\n");
    127 
    128         //open model input file for reading
    129         input_file=xNew<char>((strlen(argv[2])+strlen(argv[3])+strlen(".bin")+2));
    130         sprintf(input_file,"%s/%s%s",argv[2],argv[3],".bin");
    131         fid=fopen(input_file,"rb");
    132         if (fid==NULL) Cerr << "dirstructure error message: could not open model " << input_file << " to retrieve qmu statistics parameters" << std::endl;
    133 
    134         //initialize IoModel, but light version, we just need it to fetch one constant:
    135         iomodel=new IoModel();
    136         iomodel->fid=fid;
    137         iomodel->FetchConstants();
    138 
    139         //early return if statistics not requested:
    140         iomodel->FindConstant(&statistics,"md.qmu.statistics");
    141         if(!statistics){
    142                 delete iomodel;
    143                 fclose(fid);
    144                 return false; //important return value!
    145         }
    146 
    147         iomodel->FindConstant(&numdirectories,"md.qmu.statistics.ndirectories");
    148 
    149         /*Ok, we have everything we need to create the directory structure:*/
    150         if(IssmComm::GetRank()==0){
    151                 for (int i=0;i<numdirectories;i++){
    152                         char directory[1000];
    153                         sprintf(directory,"./%i",i+1);
    154 
    155                         check = mkdir(directory,ACCESSPERMS);
    156                         if (check) _error_("dirstructure error message: could not create directory " << directory << "\n");
    157                 }
    158         }
    159 
    160         //close model file:
    161         fclose(fid);
    162 
    163         //return value:
    164         return true; //statistics computation on!
    165 } /*}}}*/
    166 int issm_dakota_statistics(int argc,char** argv){ /*{{{*/
    167 
    168         char* input_file;
    169         FILE* fid;
    170         IoModel* iomodel=NULL;
    171         ISSM_MPI_Comm statcomm;
    172         int my_rank;
    173 
    174         //qmu statistics
    175         bool statistics    = false;
    176         int  numstatistics = 0;
    177         int  numdirectories = 0;
    178         int  nfilesperdirectory = 0;
    179         char string[1000];
    180         char* name = NULL;
    181         char** fields = NULL;
    182         int    nfields;
    183         int*   steps=NULL;
    184         int    nsteps;
    185         int    nbins;
    186         int*   indices=NULL;
    187         int    nindices;
    188         int    nsamples;
    189         int    dummy;
    190         char*  directory=NULL;
    191         char*  model=NULL;
    192         Results* results=NULL;
    193         Parameters* parameters=NULL;
    194         int color;
    195 
    196         /*First things first, set the communicator as a global variable: */
    197         IssmComm::SetComm(MPI_COMM_WORLD);
    198         my_rank=IssmComm::GetRank();
    199 
    200         /*Barrier:*/
    201         ISSM_MPI_Barrier(IssmComm::GetComm());
    202         _printf0_("Dakota Statistic Computation" << "\n");
    203 
    204         //open model input file for reading
    205         input_file=xNew<char>((strlen(argv[2])+strlen(argv[3])+strlen(".bin")+2));
    206         sprintf(input_file,"%s/%s%s",argv[2],argv[3],".bin");
    207         fid=fopen(input_file,"rb");
    208         if (fid==NULL) Cerr << "issm_dakota_statistics error message: could not open model " << input_file << " to retrieve qmu statistics parameters" << std::endl;
    209 
    210         //initialize IoModel, but light version, we'll need it to fetch constants:
    211         iomodel=new IoModel();
    212         iomodel->fid=fid;
    213         iomodel->FetchConstants();
    214 
    215         //early return if statistics not requested:
    216         iomodel->FindConstant(&statistics,"md.qmu.statistics");
    217         if(!statistics){
    218                 delete iomodel;
    219                 fclose(fid);
    220                 return 0;
    221         }
    222 
    223         //create parameters datasets with al the qmu statistics settings we need:
    224         if(statistics){
    225 
    226                 /*Initialize parameters and results:*/
    227                 results   = new Results();
    228                 parameters=new Parameters();
    229                
    230                 //solution type:
    231                 parameters->AddObject(new IntParam(SolutionTypeEnum,StatisticsSolutionEnum));
    232        
    233                 //root  directory
    234                 directory=xNew<char>(strlen(argv[2])+1);
    235                 xMemCpy<char>(directory,argv[2],strlen(argv[2])+1);
    236                 parameters->AddObject(new StringParam(DirectoryNameEnum,directory));
    237 
    238                 //model  name
    239                 model=xNew<char>(strlen(argv[3])+1);
    240                 xMemCpy<char>(model,argv[3],strlen(argv[3])+1);
    241                 parameters->AddObject(new StringParam(InputFileNameEnum,model));
    242 
    243                 //nsamples
    244                 iomodel->FindConstant(&nsamples,"md.qmu.method.params.samples");
    245                 parameters->AddObject(new IntParam(QmuNsampleEnum,nsamples));
    246 
    247                 //ndirectories
    248                 iomodel->FindConstant(&numdirectories,"md.qmu.statistics.ndirectories");
    249                 parameters->AddObject(new IntParam(QmuNdirectoriesEnum,numdirectories));
    250 
    251                 //nfiles per directory
    252                 iomodel->FindConstant(&nfilesperdirectory,"md.qmu.statistics.nfiles_per_directory");
    253                 parameters->AddObject(new IntParam(QmuNfilesPerDirectoryEnum,nfilesperdirectory));
    254 
    255                 //At this point, we don't want to go forward any longer, we want to create an MPI
    256                 //communicator on which to carry out the computations:
    257                 if ((my_rank+1)*nfilesperdirectory>nsamples)color=MPI_UNDEFINED;
    258                 else color=0;
    259                 ISSM_MPI_Comm_split(ISSM_MPI_COMM_WORLD,color, my_rank, &statcomm);
    260 
    261 
    262                 iomodel->FindConstant(&numstatistics,"md.qmu.statistics.numstatistics");
    263                 for (int i=1;i<=numstatistics;i++){
    264 
    265                         char* directory=NULL;
    266                         char* model=NULL;
    267                         int   nsamples;
    268                         _printf0_("Dealing with qmu statistical computation #" << i << "\n");
    269                
    270                         sprintf(string,"md.qmu.statistics.method(%i).name",i);
    271                         iomodel->FindConstant(&name,string);
    272 
    273                         sprintf(string,"md.qmu.statistics.method(%i).fields",i);
    274                         iomodel->FindConstant(&fields,&nfields,string);
    275                         parameters->AddObject(new StringArrayParam(FieldsEnum,fields,nfields));
    276 
    277                         sprintf(string,"md.qmu.statistics.method(%i).steps",i);
    278                         iomodel->FetchData(&steps,&dummy,&nsteps,string);
    279                         parameters->AddObject(new IntVecParam(StepsEnum,steps,nsteps));
    280 
    281                         if (strcmp(name,"Histogram")==0){
    282                                 /*fetch nbins: */
    283                                 sprintf(string,"md.qmu.statistics.method(%i).nbins",i);
    284                                 iomodel->FindConstant(&nbins,string);
    285                                 parameters->AddObject(new IntParam(NbinsEnum,nbins));
    286                                 ComputeHistogram(parameters,results,color,statcomm);
    287                         }
    288                         else if (strcmp(name,"SampleSeries")==0){
    289                                 /*fetch indices: */
    290                                 sprintf(string,"md.qmu.statistics.method(%i).indices",i);
    291                                 iomodel->FetchData(&indices,&dummy,&nindices,string);
    292                                 parameters->AddObject(new IntVecParam(IndicesEnum,indices,nindices));
    293                
    294                                 ComputeSampleSeries(parameters,results,color,statcomm);
    295                         }
    296                         else if (strcmp(name,"MeanVariance")==0){
    297                                 ComputeMeanVariance(parameters,results,color,statcomm);
    298                         }
    299                         else _error_(" error creating qmu statistics methods parameters: unsupported method " << name);
    300                 }
    301         }
    302         //close model file:
    303         fclose(fid);
    304 
    305         /*output results:*/
    306         OutputStatistics(parameters,results,color,statcomm);
    307        
    308         /*all meet here: */
    309         ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); _printf0_("Output file.\n");
    310 
    311         /*Delete ressources:*/
    312         delete parameters;
    313         delete results;
    314 
    315 
    316        
    317 } /*}}}*/

  • issm/branches/trunk-larour-SLPS2020/src/c/main/issm_post.cpp

    r25600 r25615  
    11/*!\file:  issm_post.cpp
    2  * \brief: ISSM POST processing main program
     2 * \brief: ISSM DAKOTA post-processing of statistics
    33 */
    4 /*includes and prototypes:*/
     4
    55#include "./issm.h"
     6#include <sys/stat.h>
    67
    78int main(int argc,char **argv){ /*{{{*/
    8        
    9         char* method=NULL;
    109
    11         int nfields;
    12         char*  string=NULL;
    13         char** fields=NULL;
    14         char*  field=NULL;
    15         int    offset;
    16         char*  stepstring=NULL;
    17         int    step1,step2;
    18         char*  pattern=NULL;
    19         int    nsteps;
    20         int*   steps=NULL;
    21         int    nindices;
    22         int*   indices=NULL;
    23         int    nbins;
    24        
    25         /*datasets*/
    26         Parameters  *parameters   = NULL;
    27         Results  *results   = NULL;
     10        char* dakota_input_file=NULL;
     11        char* dakota_output_file = NULL;
     12        char* dakota_error_file = NULL;
     13        bool statistics=false;
    2814
    29         /*Initialize environment (MPI, PETSC, MUMPS, etc ...)*/
    30         ISSM_MPI_Comm comm=EnvironmentInit(argc,argv);
     15        /*Initialize MPI: */
     16        ISSM_MPI_Init(&argc, &argv); // initialize MPI
    3117
    32         /*First things first, store the communicator, and set it as a global variable: */
    33         IssmComm::SetComm(comm);
    34 
    35         /*Initialize and retrieve parameters:{{{*/
    36         parameters   = new Parameters();
    37         results   = new Results();
    38 
    39         /*Method and Solution:*/
    40         method=argv[1];
    41         parameters->AddObject(new IntParam(SolutionTypeEnum,StatisticsSolutionEnum));
    42         parameters->AddObject(new StringParam(AnalysisTypeEnum,method));
    43 
    44         /*Directory:*/
    45         parameters->AddObject(new StringParam(DirectoryNameEnum,argv[2]));
    46 
    47         /*Model name:*/
    48         parameters->AddObject(new StringParam(InputFileNameEnum,argv[3]));
    49 
    50         /*Number of samples:*/
    51         parameters->AddObject(new IntParam(QmuNsampleEnum,atoi(argv[4])));
    52 
    53         /*Retrieve fields:*/
    54         nfields=atoi(argv[5]);
    55         fields=xNew<char*>(nfields);
    56         for(int i=0;i<nfields;i++){
    57                 string= argv[6+i];
    58                 field=xNew<char>(strlen(string)+1);
    59                 xMemCpy<char>(field,string,(strlen(string)+1));
    60                 fields[i]=field;
    61         }
    62         parameters->AddObject(new StringArrayParam(FieldsEnum,fields,nfields));
    63         offset=6+nfields;
    64 
    65         /*free some memory: */
    66         for(int i=0;i<nfields;i++){
    67                 char* field=fields[i];
    68                 xDelete<char>(field);
    69         }
    70         xDelete<char*>(fields);
    71 
    72         /*Retrieve time steps:*/
    73         stepstring=argv[offset];
    74         pattern=strstr(stepstring,":");
    75         if (pattern==NULL){
    76                 step1=atoi(stepstring);
    77                 step2=step1;
    78         }
    79         else{
    80                 step2=atoi(pattern+1);
    81                 stepstring[pattern-stepstring]='\0';
    82                 step1=atoi(stepstring);
    83         }
    84         nsteps=step2-step1+1;
    85         steps=xNew<int>(nsteps);
    86         for (int i=step1;i<=step2;i++)steps[i-step1]=i;
    87         parameters->AddObject(new IntVecParam(StepsEnum,steps,nsteps));
    88         offset++;
    89 
    90         /*free some memory:*/
    91         xDelete<int>(steps);
    92 
    93         /*}}}*/
    94        
    95         /*Key off method:*/
    96         if(strcmp(method,"MeanVariance")==0){
    97 
    98                 ComputeMeanVariance(parameters,results,0,ISSM_MPI_COMM_WORLD);
    99         }
    100         else if(strcmp(method,"Histogram")==0){
    101                
    102                 /*Retrieve the size of the histogram (number of bins):*/
    103                 nbins=atoi(argv[offset]); offset++;
    104                 parameters->AddObject(new IntParam(NbinsEnum,nbins));
    105                
    106                 ComputeHistogram(parameters,results,0,ISSM_MPI_COMM_WORLD);
    107 
    108         }
    109         else if(strcmp(method,"SampleSeries")==0){
    110 
    111                 /*Retrieve the vertex indices where we'll retrieve our sample time series:*/
    112                 nindices=atoi(argv[offset]); offset++;
    113                 indices=xNew<int>(nindices);
    114                 for (int i=0;i<nindices;i++){
    115                         indices[i]=atoi(argv[offset+i]);
    116                 }
    117                 parameters->AddObject(new IntVecParam(IndicesEnum,indices,nindices));
    118 
    119                 /*free some memory:*/
    120                 xDelete<int>(indices);
    121                
    122                 ComputeSampleSeries(parameters,results,0,ISSM_MPI_COMM_WORLD);
    123         }
    124 
    125         else _error_("unknown method: " << method << "\n");
    126 
    127         /*output results:*/
    128         ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); _printf0_("Output file.\n");
    129         OutputStatistics(parameters,results,0,ISSM_MPI_COMM_WORLD);
    130 
    131         /*Delete ressources:*/
    132         delete parameters;
    133         delete results;
    134 
    135         /*Finalize ISSM:*/
    136         ISSM_MPI_Finalize();
     18        /*Run statistics:*/
     19        DakotaStatistics(argc,argv);
    13720
    13821        /*Return unix success: */

  • issm/branches/trunk-larour-SLPS2020/src/c/modules/QmuStatisticsx/QmuStatisticsx.cpp

    r25605 r25615  
    22 */
    33/*includes and prototypes:*/
     4#include <sys/stat.h>
    45#include "./QmuStatisticsx.h"
    56#include "../OutputResultsx/OutputResultsx.h"
    … …  
    433434                                                IssmDouble mi=*minxs[counter];
    434435                                                int index=(scalar-mi)/(ma-mi)*nbins; if (index==nbins)index--;
     436                                                if(ma==mi)index=0;
     437                                                //_printf_( index << "|" << scalar << "|" << mi << "|" << ma << "|" << nbins << "\n");
    435438                                                localhistogram[index]++;
    436439                                                histogram[counter]=localhistogram;
    … …  
    443446                                                        IssmDouble scalar=doublemat[k];
    444447                                                        int index=(scalar-mi[k])/(ma[k]-mi[k])*nbins; if (index==nbins)index--;
     448                                                        if (mi[k]==ma[k])index=0;
    445449                                                        _assert_(scalar<=ma[k]); _assert_(scalar>=mi[k]); _assert_(index<nbins);
    446450                                                        localhistogram[k*nbins+index]++;
    … …  
    457461                                                IssmDouble mi=*minxs[counter];
    458462                                                int index=(scalar-mi)/(ma-mi)*nbins; if (index==nbins)index=nbins-1;
     463                                                if(ma==mi)index=0;
    459464                                                localhistogram[index]++;
    460465                                        }
    … …  
    466471                                                        IssmDouble scalar=doublemat[k];
    467472                                                        int index=(scalar-mi[k])/(ma[k]-mi[k])*nbins; if (index==nbins)index=nbins-1;
     473                                                        if (mi[k]==ma[k])index=0;
     474
    468475                                                        localhistogram[k*nbins+index]++;
    469476                                                }
    … …  
    496503                                        IssmDouble mi=*minmeans[f];
    497504                                        int index=(timemean-mi)/(ma-mi)*nbins; if (index==nbins)index=nbins-1;
     505                                        if(ma==mi)index=0;
    498506                                        localhistogram[index]++;
    499507                                        timehistogram[f]=localhistogram;
    … …  
    504512                                        IssmDouble mi=*minmeans[f];
    505513                                        int index=(timemean-mi)/(ma-mi)*nbins; if (index==nbins)index=nbins-1;
     514                                        if(ma==mi)index=0;
    506515                                        localhistogram[index]++;
    507516                                }
    … …  
    526535                                                IssmDouble scalar=timemean[k];
    527536                                                int index=(scalar-mi[k])/(ma[k]-mi[k])*nbins; if (index==nbins)index=nbins-1;
     537                                                if (mi[k]==ma[k])index=0;
    528538                                                localhistogram[k*nbins+index]++;
    529539                                        }
    … …  
    539549                                                IssmDouble scalar=timemean[k];
    540550                                                int index=(scalar-mi[k])/(ma[k]-mi[k])*nbins; if (index==nbins)index=nbins-1;
     551                                                if (mi[k]==ma[k])index=0;
     552
    541553                                                localhistogram[k*nbins+index]++;
    542554                                        }
    … …  
    11861198        OutputResultsx(femmodel);
    11871199} /*}}}*/
     1200bool DakotaDirStructure(int argc,char** argv){ /*{{{*/
     1201
     1202        char* input_file;
     1203        FILE* fid;
     1204        IoModel* iomodel=NULL;
     1205        int check;
     1206
     1207        //qmu statistics
     1208        bool statistics    = false;
     1209        int  numdirectories = 0;
     1210
     1211        /*First things first, set the communicator as a global variable: */
     1212        IssmComm::SetComm(MPI_COMM_WORLD);
     1213
     1214        /*Barrier:*/
     1215        ISSM_MPI_Barrier(IssmComm::GetComm());
     1216        _printf0_("Preparing directory structure for model outputs:" << "\n");
     1217
     1218        //open model input file for reading
     1219        input_file=xNew<char>((strlen(argv[2])+strlen(argv[3])+strlen(".bin")+2));
     1220        sprintf(input_file,"%s/%s%s",argv[2],argv[3],".bin");
     1221        fid=fopen(input_file,"rb");
     1222        if (fid==NULL) Cerr << "dirstructure error message: could not open model " << input_file << " to retrieve qmu statistics parameters" << std::endl;
     1223
     1224        //initialize IoModel, but light version, we just need it to fetch one constant:
     1225        iomodel=new IoModel();
     1226        iomodel->fid=fid;
     1227        iomodel->FetchConstants();
     1228
     1229        //early return if statistics not requested:
     1230        iomodel->FindConstant(&statistics,"md.qmu.statistics");
     1231        if(!statistics){
     1232                delete iomodel;
     1233                fclose(fid);
     1234                return false; //important return value!
     1235        }
     1236
     1237        iomodel->FindConstant(&numdirectories,"md.qmu.statistics.ndirectories");
     1238
     1239        /*Ok, we have everything we need to create the directory structure:*/
     1240        if(IssmComm::GetRank()==0){
     1241                for (int i=0;i<numdirectories;i++){
     1242                        char directory[1000];
     1243                        sprintf(directory,"./%i",i+1);
     1244
     1245                        check = mkdir(directory,ACCESSPERMS);
     1246                        if (check) _error_("dirstructure error message: could not create directory " << directory << "\n");
     1247                }
     1248        }
     1249
     1250        //close model file:
     1251        fclose(fid);
     1252
     1253        //return value:
     1254        return true; //statistics computation on!
     1255} /*}}}*/
     1256int DakotaStatistics(int argc,char** argv){ /*{{{*/
     1257
     1258        char* input_file;
     1259        FILE* fid;
     1260        IoModel* iomodel=NULL;
     1261        ISSM_MPI_Comm statcomm;
     1262        int my_rank;
     1263
     1264        //qmu statistics
     1265        bool statistics    = false;
     1266        int  numstatistics = 0;
     1267        int  numdirectories = 0;
     1268        int  nfilesperdirectory = 0;
     1269        char string[1000];
     1270        char* name = NULL;
     1271        char** fields = NULL;
     1272        int    nfields;
     1273        int*   steps=NULL;
     1274        int    nsteps;
     1275        int    nbins;
     1276        int*   indices=NULL;
     1277        int    nindices;
     1278        int    nsamples;
     1279        int    dummy;
     1280        char*  directory=NULL;
     1281        char*  model=NULL;
     1282        Results* results=NULL;
     1283        Parameters* parameters=NULL;
     1284        int color;
     1285
     1286        /*First things first, set the communicator as a global variable: */
     1287        IssmComm::SetComm(MPI_COMM_WORLD);
     1288        my_rank=IssmComm::GetRank();
     1289
     1290        /*Barrier:*/
     1291        ISSM_MPI_Barrier(IssmComm::GetComm());
     1292        _printf0_("Dakota Statistic Computation" << "\n");
     1293
     1294        //open model input file for reading
     1295        input_file=xNew<char>((strlen(argv[2])+strlen(argv[3])+strlen(".bin")+2));
     1296        sprintf(input_file,"%s/%s%s",argv[2],argv[3],".bin");
     1297        fid=fopen(input_file,"rb");
     1298        if (fid==NULL) Cerr << "issm_dakota_statistics error message: could not open model " << input_file << " to retrieve qmu statistics parameters" << std::endl;
     1299
     1300        //initialize IoModel, but light version, we'll need it to fetch constants:
     1301        iomodel=new IoModel();
     1302        iomodel->fid=fid;
     1303        iomodel->FetchConstants();
     1304
     1305        //early return if statistics not requested:
     1306        iomodel->FindConstant(&statistics,"md.qmu.statistics");
     1307        if(!statistics){
     1308                delete iomodel;
     1309                fclose(fid);
     1310                return 0;
     1311        }
     1312
     1313        //create parameters datasets with al the qmu statistics settings we need:
     1314        if(statistics){
     1315
     1316                /*Initialize parameters and results:*/
     1317                results   = new Results();
     1318                parameters=new Parameters();
     1319               
     1320                //solution type:
     1321                parameters->AddObject(new IntParam(SolutionTypeEnum,StatisticsSolutionEnum));
     1322       
     1323                //root  directory
     1324                directory=xNew<char>(strlen(argv[2])+1);
     1325                xMemCpy<char>(directory,argv[2],strlen(argv[2])+1);
     1326                parameters->AddObject(new StringParam(DirectoryNameEnum,directory));
     1327
     1328                //model  name
     1329                model=xNew<char>(strlen(argv[3])+1);
     1330                xMemCpy<char>(model,argv[3],strlen(argv[3])+1);
     1331                parameters->AddObject(new StringParam(InputFileNameEnum,model));
     1332
     1333                //nsamples
     1334                iomodel->FindConstant(&nsamples,"md.qmu.method.params.samples");
     1335                parameters->AddObject(new IntParam(QmuNsampleEnum,nsamples));
     1336
     1337                //ndirectories
     1338                iomodel->FindConstant(&numdirectories,"md.qmu.statistics.ndirectories");
     1339                parameters->AddObject(new IntParam(QmuNdirectoriesEnum,numdirectories));
     1340
     1341                //nfiles per directory
     1342                iomodel->FindConstant(&nfilesperdirectory,"md.qmu.statistics.nfiles_per_directory");
     1343                parameters->AddObject(new IntParam(QmuNfilesPerDirectoryEnum,nfilesperdirectory));
     1344
     1345                //At this point, we don't want to go forward any longer, we want to create an MPI
     1346                //communicator on which to carry out the computations:
     1347                if ((my_rank+1)*nfilesperdirectory>nsamples)color=MPI_UNDEFINED;
     1348                else color=0;
     1349                ISSM_MPI_Comm_split(ISSM_MPI_COMM_WORLD,color, my_rank, &statcomm);
     1350
     1351
     1352                iomodel->FindConstant(&numstatistics,"md.qmu.statistics.numstatistics");
     1353                for (int i=1;i<=numstatistics;i++){
     1354
     1355                        char* directory=NULL;
     1356                        char* model=NULL;
     1357                        int   nsamples;
     1358                        _printf0_("Dealing with qmu statistical computation #" << i << "\n");
     1359               
     1360                        sprintf(string,"md.qmu.statistics.method(%i).name",i);
     1361                        iomodel->FindConstant(&name,string);
     1362
     1363                        sprintf(string,"md.qmu.statistics.method(%i).fields",i);
     1364                        iomodel->FindConstant(&fields,&nfields,string);
     1365                        parameters->AddObject(new StringArrayParam(FieldsEnum,fields,nfields));
     1366
     1367                        sprintf(string,"md.qmu.statistics.method(%i).steps",i);
     1368                        iomodel->FetchData(&steps,&dummy,&nsteps,string);
     1369                        parameters->AddObject(new IntVecParam(StepsEnum,steps,nsteps));
     1370
     1371                        if (strcmp(name,"Histogram")==0){
     1372                                /*fetch nbins: */
     1373                                sprintf(string,"md.qmu.statistics.method(%i).nbins",i);
     1374                                iomodel->FindConstant(&nbins,string);
     1375                                parameters->AddObject(new IntParam(NbinsEnum,nbins));
     1376                                ComputeHistogram(parameters,results,color,statcomm);
     1377                        }
     1378                        else if (strcmp(name,"SampleSeries")==0){
     1379                                /*fetch indices: */
     1380                                sprintf(string,"md.qmu.statistics.method(%i).indices",i);
     1381                                iomodel->FetchData(&indices,&dummy,&nindices,string);
     1382                                parameters->AddObject(new IntVecParam(IndicesEnum,indices,nindices));
     1383               
     1384                                ComputeSampleSeries(parameters,results,color,statcomm);
     1385                        }
     1386                        else if (strcmp(name,"MeanVariance")==0){
     1387                                ComputeMeanVariance(parameters,results,color,statcomm);
     1388                        }
     1389                        else _error_(" error creating qmu statistics methods parameters: unsupported method " << name);
     1390                }
     1391        }
     1392        //close model file:
     1393        fclose(fid);
     1394
     1395        /*output results:*/
     1396        OutputStatistics(parameters,results,color,statcomm);
     1397       
     1398        /*all meet here: */
     1399        ISSM_MPI_Barrier(ISSM_MPI_COMM_WORLD); _printf0_("Output file.\n");
     1400
     1401        /*Delete ressources:*/
     1402        delete parameters;
     1403        delete results;
     1404
     1405
     1406       
     1407} /*}}}*/

  • issm/branches/trunk-larour-SLPS2020/src/c/modules/QmuStatisticsx/QmuStatisticsx.h

    r25599 r25615  
    1212int ComputeHistogram(Parameters* parameters,Results* results,int color, ISSM_MPI_Comm statcomm);
    1313int readdata(IssmDouble** pdoublemat, int* pdoublematsize, IssmDouble* pdouble, FILE* fid,char* field,int step);
     14bool DakotaDirStructure(int argc,char** argv);
     15int DakotaStatistics(int argc,char** argv);
    1416
    1517/* local prototypes: */
Note: See TracChangeset for help on using the changeset viewer.