Changeset 2478 for XIOS3/branches/xios-3.0-beta/src/io/onetcdf4.cpp
- Timestamp:
- 03/21/23 10:51:27 (16 months ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
XIOS3/branches/xios-3.0-beta/src/io/onetcdf4.cpp
r2264 r2478 314 314 std::vector<StdSize> dimsizes; 315 315 int dimSize = dim.size(); 316 317 int grpid = this->getCurrentGroup(); 318 319 std::vector<StdString>::const_iterator it = dim.begin(), end = dim.end(); 320 321 for (int idx = 0; it != end; it++, ++idx) 322 { 323 const StdString& dimid = *it; 324 dimids.push_back(this->getDimension(dimid)); 325 } 326 327 CNetCdfInterface::defVar(grpid, name, type, dimids.size(), &dimids[0], varid); 328 329 return varid; 330 } 331 332 //--------------------------------------------------------------- 333 334 int CONetCDF4::addChunk(CField* field, nc_type type, 335 const std::vector<StdString>& dim, int compressionLevel) 336 { 337 const StdString& name = field->getFieldOutputName(); 338 int varid = 0; 339 std::vector<StdSize> dimsizes; 340 int dimSize = dim.size(); 316 341 317 342 StdSize size; 318 343 StdSize totalSize; 319 StdSize maxSize = 1024 * 1024 * 256; // == 2GB/8 if output double 344 345 // default chunk size 346 StdSize maxSize = 1024 * 1024 * 20; // = 20 Mo (exp using large NEMO like grids) 347 StdSize targetSize = maxSize; 348 if (!field->chunking_blocksize_target.isEmpty()) 349 { 350 targetSize = field->chunking_blocksize_target.getValue()*1024*1024; 351 } 352 StdSize recordSize = 1; //sizeof(type); 353 if (field->prec.isEmpty()) recordSize *= 4; 354 else recordSize *= field->prec; 320 355 321 356 int grpid = this->getCurrentGroup(); … … 326 361 { 327 362 const StdString& dimid = *it; 328 dimids.push_back(this->getDimension(dimid));329 363 CNetCdfInterface::inqDimLen(grpid, this->getDimension(dimid), size); 330 364 if (size == NC_UNLIMITED) size = 1; 365 recordSize *= size; 331 366 dimsizes.push_back(size); 332 367 } 333 334 CNetCdfInterface::defVar(grpid, name, type, dimids.size(), &dimids[0], varid); 368 double chunkingRatio = (double)recordSize / (double)targetSize; 369 370 varid = this->getVariable(name); 335 371 336 372 // The classic NetCDF format does not support chunking nor fill parameters 337 373 if (!useClassicFormat) 338 374 { 339 // set chunksize : size of one record 340 // but must not be > 2GB (netcdf or HDF5 problem) 341 totalSize = 1; 342 for (vector<StdSize>::reverse_iterator it = dimsizes.rbegin(); it != dimsizes.rend(); ++it) 375 // Browse field's elements (domain, axis) and NetCDF meta-data to retrieve corresponding chunking coefficients 376 // in the file order ! 377 CGrid* grid = field->getGrid(); 378 std::vector<CDomain*> domains = grid->getDomains(); 379 std::vector<CAxis*> axis = grid->getAxis(); 380 381 std::vector<double> userChunkingWeights; // store chunking coefficients defined by users 382 std::vector<StdString>::const_reverse_iterator itId = dim.rbegin(); // NetCDF is fed using dim, see this::addVariable() 383 int elementIdx(0); 384 int outerAxis(-1), outerDomJ(-1), outerDomI(-1); 385 for (vector<StdSize>::reverse_iterator itDim = dimsizes.rbegin(); itDim != dimsizes.rend(); ++itDim, ++itId, ++elementIdx) 343 386 { 344 totalSize *= *it; 345 if (totalSize >= maxSize) *it = 1; 387 bool isAxis(false); 388 for ( auto ax : axis ) 389 { 390 StdString axisDim; 391 // Rebuild axis name from axis before comparing 392 if (ax->dim_name.isEmpty()) axisDim = ax->getAxisOutputName(); 393 else axisDim=ax->dim_name.getValue(); 394 if (axisDim == *itId) 395 { 396 if (!ax->chunking_weight.isEmpty()) userChunkingWeights.push_back( ax->chunking_weight.getValue() ); 397 else userChunkingWeights.push_back( 0. ); 398 isAxis = true; 399 outerAxis = elementIdx; // Going backward in dimsizes, overwriting will keep outer 400 break; 401 } 402 } // end scanning axis 403 bool isDomain(false); 404 for ( auto dom : domains ) 405 { 406 StdString axisDim; 407 // Rebuild axis I name from domain before comparing 408 if (!dom->dim_i_name.isEmpty()) axisDim = dom->dim_i_name.getValue(); 409 else 410 { 411 if (dom->type==CDomain::type_attr::curvilinear) axisDim="x"; 412 if (dom->type==CDomain::type_attr::unstructured) axisDim="cell"; 413 if (dom->type==CDomain::type_attr::rectilinear) axisDim="lon"; 414 } 415 if (axisDim == *itId) 416 { 417 if (!dom->chunking_weight_i.isEmpty()) userChunkingWeights.push_back( dom->chunking_weight_i.getValue() ); 418 else userChunkingWeights.push_back( 0. ); 419 isDomain = true; 420 outerDomI = elementIdx; // Going backward in dimsizes, overwriting will keep outer 421 break; 422 } 423 // Rebuild axis J name from domain before comparing 424 if (!dom->dim_j_name.isEmpty()) axisDim = dom->dim_j_name.getValue(); 425 else { 426 if (dom->type==CDomain::type_attr::curvilinear) axisDim="y"; 427 if (dom->type==CDomain::type_attr::rectilinear) axisDim="lat"; 428 } 429 if (axisDim == *itId) 430 { 431 if (!dom->chunking_weight_j.isEmpty()) userChunkingWeights.push_back( dom->chunking_weight_j.getValue() ); 432 else userChunkingWeights.push_back( 0. ); 433 outerDomJ = elementIdx; // Going backward in dimsizes, overwriting will keep outer 434 isDomain = true; 435 break; 436 } 437 } // end scanning domain 438 // No chunking applied on scalars or time 439 if ((!isAxis)&&(!isDomain)) 440 { 441 userChunkingWeights.push_back( 0. ); 442 } 443 } 444 445 double sumChunkingWeights(0); 446 for (int i=0;i<userChunkingWeights.size();i++) sumChunkingWeights += userChunkingWeights[i]; 447 if ( (!sumChunkingWeights) && (chunkingRatio > 1) ) 448 { 449 if (outerAxis!=-1) userChunkingWeights[outerAxis] = 1; // chunk along outer axis 450 else if (outerDomJ!=-1) userChunkingWeights[outerDomJ] = 1; // if no axis ? -> along j 451 else if (outerDomI!=-1) userChunkingWeights[outerDomI] = 1; // if no j -> along i 452 else {;} 453 } 454 455 int countChunkingDims(0); // number of dimensions on which chunking is operated : algo uses pow(value, 1/countChunkingDims) 456 double normalizingWeight(0); // use to relativize chunking coefficients for all dimensions 457 for (int i=0;i<userChunkingWeights.size();i++) 458 if (userChunkingWeights[i]>0.) 459 { 460 countChunkingDims++; 461 normalizingWeight = userChunkingWeights[i]; 462 } 463 464 std::vector<double> chunkingRatioPerDims; // will store coefficients used to compute chunk size 465 double productRatios = 1; // last_coeff = pow( shrink_ratio / (product of all ratios), 1/countChunkingDims ) 466 for (int i=0;i<userChunkingWeights.size();i++) 467 { 468 chunkingRatioPerDims.push_back( userChunkingWeights[i] / normalizingWeight ); 469 if (chunkingRatioPerDims[i]) productRatios *= chunkingRatioPerDims[i]; 470 } 471 for (int i=0;i<userChunkingWeights.size();i++) 472 { 473 chunkingRatioPerDims[i] *= pow( chunkingRatio / productRatios, 1./countChunkingDims ); 474 } 475 476 std::vector<double>::iterator itChunkingRatios = chunkingRatioPerDims.begin(); 477 //itId = dim.rbegin(); 478 double correctionFromPreviousDim = 1.; 479 for (vector<StdSize>::reverse_iterator itDim = dimsizes.rbegin(); itDim != dimsizes.rend(); ++itDim, ++itChunkingRatios, ++itId) 480 { 481 *itChunkingRatios *= correctionFromPreviousDim; 482 correctionFromPreviousDim = 1; 483 if (*itChunkingRatios > 1) // else target larger than size ! 484 { 485 StdSize dimensionSize = *itDim; 486 //info(0) << *itId << " " << *itDim << " " << *itChunkingRatios << " " << (*itDim)/(*itChunkingRatios) << endl; 487 *itDim = ceil( *itDim / ceil(*itChunkingRatios) ); 488 correctionFromPreviousDim = *itChunkingRatios/ ((double)dimensionSize/(*itDim)); 489 } 346 490 } 347 491 int storageType = (0 == dimSize) ? NC_CONTIGUOUS : NC_CHUNKED;
Note: See TracChangeset
for help on using the changeset viewer.