source: XIOS/dev/XIOS_DEV_CMIP6/src/transformation/domain_algorithm_reorder.cpp @ 1457

Last change on this file since 1457 was 1457, checked in by ymipsl, 6 years ago

Add new domain filter : reorder_domain
Reoder the data along the global domain but works only for rectilinear domain

  • invert_lat : invert the latitute axis
  • shift_lon_fraction : shift the longitude axis of a fration of global size
  • lon_min/lon_max : fixe the range of longitude value (ex : -180:180 or 0:360)

YM

  • Property svn:eol-style set to native
File size: 5.2 KB
Line 
1/*!
2   \file domain_algorithm_reorder.cpp
3   \brief Algorithm for reorder a domain.
4 */
5#include "domain_algorithm_reorder.hpp"
6#include "reorder_domain.hpp"
7#include "domain.hpp"
8#include "grid.hpp"
9#include "grid_transformation_factory_impl.hpp"
10
11namespace xios {
12CGenericAlgorithmTransformation* CDomainAlgorithmReorder::create(CGrid* gridDst, CGrid* gridSrc,
13                                                             CTransformation<CDomain>* transformation,
14                                                             int elementPositionInGrid,
15                                                             std::map<int, int>& elementPositionInGridSrc2ScalarPosition,
16                                                             std::map<int, int>& elementPositionInGridSrc2AxisPosition,
17                                                             std::map<int, int>& elementPositionInGridSrc2DomainPosition,
18                                                             std::map<int, int>& elementPositionInGridDst2ScalarPosition,
19                                                             std::map<int, int>& elementPositionInGridDst2AxisPosition,
20                                                             std::map<int, int>& elementPositionInGridDst2DomainPosition)
21{
22  std::vector<CDomain*> domainListDestP = gridDst->getDomains();
23  std::vector<CDomain*> domainListSrcP  = gridSrc->getDomains();
24
25  CReorderDomain* reorderDomain = dynamic_cast<CReorderDomain*> (transformation);
26  int domainDstIndex = elementPositionInGridDst2DomainPosition[elementPositionInGrid];
27  int domainSrcIndex = elementPositionInGridSrc2DomainPosition[elementPositionInGrid];
28
29  return (new CDomainAlgorithmReorder(domainListDestP[domainDstIndex], domainListSrcP[domainSrcIndex], reorderDomain));
30}
31
32bool CDomainAlgorithmReorder::registerTrans()
33{
34  CGridTransformationFactory<CDomain>::registerTransformation(TRANS_REORDER_DOMAIN, create);
35}
36
37CDomainAlgorithmReorder::CDomainAlgorithmReorder(CDomain* domainDestination, CDomain* domainSource, CReorderDomain* reorderDomain)
38: CDomainAlgorithmTransformation(domainDestination, domainSource)
39{
40  reorderDomain->checkValid(domainSource);
41  if (domainDestination->type !=  CDomain::type_attr::rectilinear)
42  {
43      ERROR("CDomainAlgorithmReorder::CDomainAlgorithmReorder(CDomain* domainDestination, CDomain* domainSource, CReorderDomain* reorderDomain)",
44           << "Domain destination is not rectilinear. This filter work only for rectilinear domain and destination domain with < id = "
45           <<domainDestination->getId() <<" > is of type "<<domainDestination->type<<std::endl);
46  }
47 
48  if (domainDestination == domainSource)
49  {
50    ERROR("CDomainAlgorithmReorder::CDomainAlgorithmReorder(CDomain* domainDestination, CDomain* domainSource, CReorderDomain* reorderDomain)",
51           << "Domain source and domain destination are the same. Please make sure domain destination refers to domain source" << std::endl
52           << "Domain source " <<domainSource->getId() << std::endl
53           << "Domain destination " <<domainDestination->getId() << std::endl);
54  }
55  this->type_ = (ELEMENT_MODIFICATION_WITHOUT_DATA);
56
57  if (!reorderDomain->invert_lat.isEmpty())
58  {
59    CArray<int,1>& j_index=domainDestination->j_index ;
60    int nglo = j_index.numElements() ;
61    int nj_glo =domainDestination->nj_glo ;
62   
63    for (size_t i = 0; i < nglo ; ++i)
64    {
65      j_index(i)=(nj_glo-1)-j_index(i) ;
66    }
67  }
68
69  if (!reorderDomain->shift_lon_fraction.isEmpty())
70  {
71    int ni_glo =domainDestination->ni_glo ;
72    int  offset = ni_glo*reorderDomain->shift_lon_fraction ;
73    CArray<int,1>& i_index=domainDestination->i_index ;
74    int nglo = i_index.numElements() ;
75   
76    for (size_t i = 0; i < nglo ; ++i)
77    {
78      i_index(i)=  (i_index(i)+offset+ni_glo)%ni_glo ;
79    }
80  }
81
82  if (!reorderDomain->min_lon.isEmpty() && !reorderDomain->max_lon.isEmpty())
83  {
84    double min_lon=reorderDomain->min_lon ;
85    double max_lon=reorderDomain->max_lon ;
86    double delta=max_lon-min_lon ;
87   
88    if (!domainDestination->lonvalue_1d.isEmpty() )
89    {
90      CArray<double,1>& lon=domainDestination->lonvalue_1d ;
91      for (int i=0;i<lon.numElements();++i)
92      {
93        while  (lon(i) > max_lon) lon(i)=lon(i)-delta ;
94        while  (lon(i) < min_lon) lon(i)=lon(i)+delta ;
95      }
96    }
97
98    if (!domainDestination->bounds_lon_1d.isEmpty() )
99    {
100      CArray<double,2>& bounds_lon=domainDestination->bounds_lon_1d ;
101      for (int i=0;i<bounds_lon.extent(0);++i)
102      {
103        while  (bounds_lon(0,i) > max_lon) bounds_lon(0,i)=bounds_lon(0,i)-delta ;
104        while  (bounds_lon(1,i) > max_lon) bounds_lon(1,i)=bounds_lon(1,i)-delta ;
105
106        while  (bounds_lon(0,i) < min_lon) bounds_lon(0,i)=bounds_lon(0,i)+delta ;
107        while  (bounds_lon(1,i) < min_lon) bounds_lon(1,i)=bounds_lon(1,i)+delta ;
108      }
109    }
110  }
111   
112 
113}
114
115/*!
116  Compute the index mapping between domain on grid source and one on grid destination
117*/
118void CDomainAlgorithmReorder::computeIndexSourceMapping_(const std::vector<CArray<double,1>* >& dataAuxInputs)
119{
120/*
121  this->transformationMapping_.resize(1);
122  this->transformationWeight_.resize(1);
123
124  TransformationIndexMap& transMap = this->transformationMapping_[0];
125  TransformationWeightMap& transWeight = this->transformationWeight_[0];
126*/
127}
128
129
130}
Note: See TracBrowser for help on using the repository browser.