@@ -40,7 +40,6 @@ void AdjustWindsAtPoles::print(std::ostream& out) const {
40
40
void AdjustWindsAtPoles::execute (context::Context& ctx) const {
41
41
auto & field = ctx.field ();
42
42
43
- ASSERT (field.dimensions () % 2 == 0 );
44
43
ASSERT (field.dimensions () > 0 );
45
44
46
45
const auto missingValue = field.hasMissing () ? field.missingValue () : std::numeric_limits<double >::quiet_NaN ();
@@ -65,22 +64,38 @@ void AdjustWindsAtPoles::execute(context::Context& ctx) const {
65
64
}
66
65
}
67
66
68
-
69
- // set v to zero at pole(s)
70
- for (size_t i = 0 ; i < field.dimensions (); i += 2 ) {
71
- auto & u = field.direct (i);
72
- auto & v = field.direct (i + 1 );
73
- ASSERT (u.size () == N);
74
- ASSERT (v.size () == N);
75
-
76
- for (const auto & indices : {north, south}) {
77
- for (const auto & idx : indices) {
78
- if (u[idx] == missingValue || v[idx] == missingValue) {
79
- u[idx] = missingValue;
80
- v[idx] = missingValue;
67
+ if (field.dimensions () % 2 == 0 ) {
68
+ // two-component mode: set v to zero at pole(s)
69
+ for (size_t i = 0 ; i < field.dimensions (); i += 2 ) {
70
+ auto & u = field.direct (i);
71
+ auto & v = field.direct (i + 1 );
72
+ ASSERT (u.size () == N);
73
+ ASSERT (v.size () == N);
74
+
75
+ for (const auto & indices : {north, south}) {
76
+ for (const auto & idx : indices) {
77
+ if (u[idx] == missingValue || v[idx] == missingValue) {
78
+ u[idx] = missingValue;
79
+ v[idx] = missingValue;
80
+ }
81
+ else {
82
+ v[idx] = 0 ;
83
+ }
81
84
}
82
- else {
83
- v[idx] = 0 ;
85
+ }
86
+ }
87
+ }
88
+ else {
89
+ // single-component mode: set u, v to zero at pole(s)
90
+ for (size_t i = 0 ; i < field.dimensions (); ++i) {
91
+ auto & v = field.direct (i);
92
+ ASSERT (v.size () == N);
93
+
94
+ for (const auto & indices : {north, south}) {
95
+ for (const auto & idx : indices) {
96
+ if (v[idx] != missingValue) {
97
+ v[idx] = 0 .;
98
+ }
84
99
}
85
100
}
86
101
}
0 commit comments