- Timestamp:
- 01/23/19 10:31:44 (5 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
XIOS/dev/branch_openmp/src/transformation/domain_algorithm_zoom.cpp
r1460 r1642 1 /*!2 \file domain_algorithm_zoom.cpp3 \author Ha NGUYEN4 \since 02 Jul 20155 \date 02 Jul 20156 7 \brief Algorithm for zooming on an domain.8 */9 1 #include "domain_algorithm_zoom.hpp" 10 2 #include "zoom_domain.hpp" … … 12 4 #include "grid.hpp" 13 5 #include "grid_transformation_factory_impl.hpp" 6 #include "attribute_template.hpp" 14 7 15 8 namespace xios { … … 23 16 std::map<int, int>& elementPositionInGridDst2AxisPosition, 24 17 std::map<int, int>& elementPositionInGridDst2DomainPosition) 18 TRY 25 19 { 26 20 std::vector<CDomain*> domainListDestP = gridDst->getDomains(); … … 33 27 return (new CDomainAlgorithmZoom(domainListDestP[domainDstIndex], domainListSrcP[domainSrcIndex], zoomDomain)); 34 28 } 29 CATCH 35 30 36 31 bool CDomainAlgorithmZoom::registerTrans() 32 TRY 37 33 { 38 34 CGridTransformationFactory<CDomain>::registerTransformation(TRANS_ZOOM_DOMAIN, create); 39 35 } 36 CATCH 40 37 41 38 CDomainAlgorithmZoom::CDomainAlgorithmZoom(CDomain* domainDestination, CDomain* domainSource, CZoomDomain* zoomDomain) 42 39 : CDomainAlgorithmTransformation(domainDestination, domainSource) 40 TRY 43 41 { 44 42 zoomDomain->checkValid(domainSource); … … 67 65 << "Zoom size is " << zoomNj_ ); 68 66 } 69 } 67 68 // Calculate the size of local domain 69 int ind, indLocSrc, indLocDest, iIdxSrc, jIdxSrc, destIBegin = -1, destJBegin = -1, niDest = 0, njDest = 0, ibeginDest, jbeginDest ; 70 int indGloDest, indGloSrc, niGloSrc = domainSrc_->ni_glo, iSrc, jSrc, nvertex = 0; 71 for (int j = 0; j < domainSrc_->nj.getValue(); j++) 72 { 73 for (int i = 0; i < domainSrc_->ni.getValue(); i++) 74 { 75 ind = j*domainSrc_->ni + i; 76 iIdxSrc = domainSrc_->i_index(ind); 77 if ((iIdxSrc >= zoomIBegin_) && (iIdxSrc <= zoomIEnd_)) 78 { 79 jIdxSrc = domainSrc_->j_index(ind); 80 if ((jIdxSrc >= zoomJBegin_) && (jIdxSrc <= zoomJEnd_)) 81 { 82 if ((niDest == 0) && (njDest == 0)) 83 { 84 destIBegin = i; 85 destJBegin = j; 86 } 87 if (i == destIBegin) ++njDest; 88 } 89 if (j == destJBegin) ++niDest; 90 91 } 92 } 93 } 94 ibeginDest = destIBegin + domainSrc_->ibegin - zoomIBegin_; 95 jbeginDest = destJBegin + domainSrc_->jbegin - zoomJBegin_; 96 domainDest_->ni_glo.setValue(zoomNi_); 97 domainDest_->nj_glo.setValue(zoomNj_); 98 domainDest_->ni.setValue(niDest); 99 domainDest_->nj.setValue(njDest); 100 if ( (niDest==0) || (njDest==0)) 101 { 102 domainDest_->ibegin.setValue(0); 103 domainDest_->jbegin.setValue(0); 104 } 105 else 106 { 107 domainDest_->ibegin.setValue(ibeginDest); 108 domainDest_->jbegin.setValue(jbeginDest); 109 } 110 domainDest_->i_index.resize(niDest*njDest); 111 domainDest_->j_index.resize(niDest*njDest); 112 113 domainDest_->data_ni.setValue(niDest); 114 domainDest_->data_nj.setValue(njDest); 115 domainDest_->data_ibegin.setValue(0); 116 domainDest_->data_jbegin.setValue(0); 117 domainDest_->data_i_index.resize(niDest*njDest); 118 domainDest_->data_j_index.resize(niDest*njDest); 119 120 domainDest_->domainMask.resize(niDest*njDest); 121 122 if (!domainSrc_->lonvalue_1d.isEmpty()) 123 { 124 if (domainDest_->type == CDomain::type_attr::rectilinear) 125 { 126 domainDest_->lonvalue_1d.resize(niDest); 127 domainDest_->latvalue_1d.resize(njDest); 128 } 129 else if (domainDest_->type == CDomain::type_attr::unstructured) 130 { 131 domainDest_->lonvalue_1d.resize(niDest); 132 domainDest_->latvalue_1d.resize(niDest); 133 } 134 else if (domainDest_->type == CDomain::type_attr::curvilinear) 135 { 136 domainDest_->lonvalue_1d.resize(niDest*njDest); 137 domainDest_->latvalue_1d.resize(niDest*njDest); 138 } 139 } 140 else if (!domainSrc_->lonvalue_2d.isEmpty()) 141 { 142 domainDest_->lonvalue_2d.resize(niDest,njDest); 143 domainDest_->latvalue_2d.resize(niDest,njDest); 144 } 145 146 if (domainSrc_->hasBounds) 147 { 148 nvertex = domainSrc_->nvertex; 149 domainDest_->nvertex.setValue(nvertex); 150 if (!domainSrc_->bounds_lon_1d.isEmpty()) 151 { 152 if (domainDest_->type == CDomain::type_attr::rectilinear) 153 { 154 domainDest_->bounds_lon_1d.resize(nvertex, niDest); 155 domainDest_->bounds_lat_1d.resize(nvertex, njDest); 156 } 157 else if (domainDest_->type == CDomain::type_attr::unstructured) 158 { 159 domainDest_->bounds_lon_1d.resize(nvertex, niDest); 160 domainDest_->bounds_lat_1d.resize(nvertex, niDest); 161 } 162 else if (domainDest_->type == CDomain::type_attr::curvilinear) 163 { 164 domainDest_->bounds_lon_1d.resize(nvertex, niDest*njDest); 165 domainDest_->bounds_lat_1d.resize(nvertex, niDest*njDest); 166 } 167 } 168 else if (!domainSrc_->bounds_lon_2d.isEmpty()) 169 { 170 domainDest_->bounds_lon_2d.resize(nvertex, niDest, njDest); 171 domainDest_->bounds_lat_2d.resize(nvertex, niDest, njDest); 172 } 173 } 174 if (domainSrc_->hasArea) domainDest_->area.resize(niDest,njDest); 175 176 this->transformationMapping_.resize(1); 177 this->transformationWeight_.resize(1); 178 TransformationIndexMap& transMap = this->transformationMapping_[0]; 179 TransformationWeightMap& transWeight = this->transformationWeight_[0]; 180 181 for (int iDest = 0; iDest < niDest; iDest++) 182 { 183 iSrc = iDest + destIBegin; 184 for (int jDest = 0; jDest < njDest; jDest++) 185 { 186 jSrc = jDest + destJBegin; 187 ind = jSrc * domainSrc_->ni + iSrc; 188 iIdxSrc = domainSrc_->i_index(ind); 189 jIdxSrc = domainSrc_->j_index(ind); 190 indLocDest = jDest*niDest + iDest; 191 indGloDest = (jDest + jbeginDest)*zoomNi_ + (iDest + ibeginDest); 192 indLocSrc = (jDest+destJBegin)*domainSrc_->ni + (iDest+destIBegin); 193 indGloSrc = (jIdxSrc )* niGloSrc + iIdxSrc; 194 domainDest_->i_index(indLocDest) = iDest + ibeginDest; // i_index contains global positions 195 domainDest_->j_index(indLocDest) = jDest + jbeginDest; // i_index contains global positions 196 domainDest_->data_i_index(indLocDest) = (domainSrc_->data_dim == 1) ? indLocDest : iDest; // data_i_index contains local positions 197 domainDest_->data_j_index(indLocDest) = (domainSrc_->data_dim == 1) ? 0 :jDest; // data_i_index contains local positions 198 domainDest_->domainMask(indLocDest) = domainSrc_->domainMask(indLocSrc); 199 200 if (domainSrc_->hasArea) 201 domainDest_->area(iDest,jDest) = domainSrc_->area(iSrc,jSrc); 202 203 if (domainSrc_->hasLonLat) 204 { 205 if (!domainSrc_->latvalue_1d.isEmpty()) 206 { 207 if (domainDest_->type == CDomain::type_attr::rectilinear) 208 { 209 domainDest_->latvalue_1d(jDest) = domainSrc_->latvalue_1d(jSrc); 210 } 211 else 212 { 213 domainDest_->lonvalue_1d(indLocDest) = domainSrc_->lonvalue_1d(ind); 214 domainDest_->latvalue_1d(indLocDest) = domainSrc_->latvalue_1d(ind); 215 } 216 } 217 else if (!domainSrc_->latvalue_2d.isEmpty()) 218 { 219 domainDest_->lonvalue_2d(iDest,jDest) = domainSrc_->lonvalue_2d(iSrc,jSrc); 220 domainDest_->latvalue_2d(iDest,jDest) = domainSrc_->latvalue_2d(iSrc,jSrc); 221 } 222 } 223 224 if (domainSrc_->hasBounds) 225 { 226 if (!domainSrc_->bounds_lon_1d.isEmpty()) 227 { 228 if (domainDest_->type == CDomain::type_attr::rectilinear) 229 { 230 for (int n = 0; n < nvertex; ++n) 231 domainDest_->bounds_lat_1d(n,jDest) = domainSrc_->bounds_lat_1d(n,jSrc); 232 } 233 else 234 { 235 for (int n = 0; n < nvertex; ++n) 236 { 237 domainDest_->bounds_lon_1d(n,indLocDest) = domainSrc_->bounds_lon_1d(n,ind); 238 domainDest_->bounds_lat_1d(n,indLocDest) = domainSrc_->bounds_lat_1d(n,ind); 239 } 240 } 241 } 242 else if (!domainSrc_->bounds_lon_2d.isEmpty()) 243 { 244 for (int n = 0; n < nvertex; ++n) 245 { 246 domainDest_->bounds_lon_2d(n,iDest,jDest) = domainSrc_->bounds_lon_2d(n,iSrc,jSrc); 247 domainDest_->bounds_lat_2d(n,iDest,jDest) = domainSrc_->bounds_lat_2d(n,iSrc,jSrc); 248 } 249 } 250 251 } 252 253 transMap[indGloDest].push_back(indGloSrc); 254 transWeight[indGloDest].push_back(1.0); 255 } 256 257 if (domainSrc_->hasLonLat && !domainSrc_->latvalue_1d.isEmpty()) 258 { 259 if (domainDest_->type == CDomain::type_attr::rectilinear) 260 { 261 domainDest_->lonvalue_1d(iDest) = domainSrc_->lonvalue_1d(iSrc); 262 } 263 } 264 265 if (domainSrc_->hasBounds && !domainSrc_->bounds_lon_1d.isEmpty()) 266 { 267 if (domainDest_->type == CDomain::type_attr::rectilinear) 268 { 269 for (int n = 0; n < nvertex; ++n) 270 domainDest_->bounds_lon_1d(n,iDest) = domainSrc_->bounds_lon_1d(n,iSrc); 271 } 272 } 273 } 274 domainDest_->computeLocalMask(); 275 } 276 CATCH 70 277 71 278 /*! … … 74 281 void CDomainAlgorithmZoom::computeIndexSourceMapping_(const std::vector<CArray<double,1>* >& dataAuxInputs) 75 282 { 76 77 int niGlob = domainSrc_->ni_glo.getValue(); 78 int njGlob = domainSrc_->nj_glo.getValue(); 79 80 this->transformationMapping_.resize(1); 81 this->transformationWeight_.resize(1); 82 83 TransformationIndexMap& transMap = this->transformationMapping_[0]; 84 TransformationWeightMap& transWeight = this->transformationWeight_[0]; 85 86 int domainGlobalIndex; 87 int iglob ; 88 int jglob ; 89 const CArray<int,1>& i_index = domainSrc_->i_index.getValue() ; 90 const CArray<int,1>& j_index = domainSrc_->j_index.getValue() ; 91 92 int nglo = i_index.numElements() ; 93 for (size_t i = 0; i < nglo ; ++i) 94 { 95 iglob=i_index(i) ; jglob=j_index(i) ; 96 if (iglob>=zoomIBegin_ && iglob<=zoomIEnd_ && jglob>=zoomJBegin_ && jglob<=zoomJEnd_) 97 { 98 domainGlobalIndex = jglob*niGlob + iglob; 99 transMap[domainGlobalIndex].push_back(domainGlobalIndex); 100 transWeight[domainGlobalIndex].push_back(1.0); 101 } 102 } 103 updateZoom(); 104 } 105 106 /*! 107 After a zoom on domain, it should be certain that (global) zoom begin and (global) zoom size are updated 108 */ 109 void CDomainAlgorithmZoom::updateZoom() 110 { 111 domainDest_->global_zoom_ibegin = zoomIBegin_; 112 domainDest_->global_zoom_jbegin = zoomJBegin_; 113 domainDest_->global_zoom_ni = zoomNi_; 114 domainDest_->global_zoom_nj = zoomNj_; 115 } 116 117 /*! 118 Update mask on domain 119 Because only zoomed region on domain is not masked, the remaining must be masked to make sure 120 correct index be extracted 121 */ 122 // void CDomainAlgorithmZoom::updateDomainDestinationMask() 123 // { 124 // int niMask = domainDest_->ni.getValue(); 125 // int iBeginMask = domainDest_->ibegin.getValue(); 126 // int njMask = domainDest_->nj.getValue(); 127 // int jBeginMask = domainDest_->jbegin.getValue(); 128 // int niGlob = domainDest_->ni_glo.getValue(); 129 // int globalIndexMask = 0; 130 131 // TransformationIndexMap& transMap = this->transformationMapping_[0]; 132 // TransformationIndexMap::const_iterator ite = (transMap).end(); 133 // for (int j = 0; j < njMask; ++j) 134 // { 135 // for (int i = 0; i < niMask; ++i) 136 // { 137 // globalIndexMask = (j+jBeginMask) * niGlob + (i + iBeginMask); 138 // if (transMap.find(globalIndexMask) == ite) 139 // (domainDest_->mask_1d)(i+j*niMask) = false; 140 // } 141 // } 142 // } 143 144 } 283 } 284 285 286 }
Note: See TracChangeset
for help on using the changeset viewer.