40#include "matpack_data.h"
94 as.insert(as.begin(), name);
101 const Tensor4 dummy = af.
data;
104 af.
resize(nf + 1, dummy.npages(), dummy.nrows(), dummy.ncols());
108 af.
data(Range((prepend && nf > 1) ? 1 : 0, nf - 1), joker, joker, joker) =
132 ArrayOfLagrangeLogInterpolation& lag_p,
134 ConstVectorView p_grid_out,
135 ConstVectorView p_grid_in,
136 const Index& interp_order,
140 out2 <<
" Interpolation order: " << interp_order <<
"\n";
143 ing_max = p_grid_out.nelem() - 1;
145 "Atmospheric field to p_grid_out", p_grid_in, p_grid_out, interp_order);
147 Index nelem_in_range = ing_max - ing_min + 1;
150 if (nelem_in_range > 0) {
151 lag_p = my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(
152 p_grid_out, p_grid_in, interp_order, 0.5);
159 Tensor3& atmtensor_out,
161 const Tensor3& atmtensor_in_orig,
162 const Vector& p_grid_new,
163 const Vector& p_grid_old,
164 const Index& interp_order,
169 "p_grid_old is supposed to be the p_grid associated with the "
170 "atmospheric field.\n"
171 "However, it is not as their sizes are inconsistent.\n")
173 const Tensor3* atmtensor_in_pnt;
174 Tensor3 atmtensor_in_copy;
176 if (&atmtensor_in_orig == &atmtensor_out) {
177 atmtensor_in_copy = atmtensor_in_orig;
178 atmtensor_in_pnt = &atmtensor_in_copy;
180 atmtensor_in_pnt = &atmtensor_in_orig;
182 const Tensor3& atmtensor_in = *atmtensor_in_pnt;
185 atmtensor_out.resize(
186 p_grid_new.nelem(), atmtensor_in.nrows(), atmtensor_in.ncols());
188 ArrayOfLagrangeLogInterpolation lag_p;
191 Index ing_min, ing_max;
204 (ing_max - ing_min + 1 != p_grid_new.nelem()),
205 "New grid seems not to be sufficiently covered by old grid.\n")
207 for (Index i = 0; i < atmtensor_in.nrows(); i++)
208 for (Index j = 0; j < atmtensor_in.ncols(); j++)
209 reinterp(atmtensor_out(joker, i, j), atmtensor_in(joker, i, j), itw, lag_p);
214 Tensor4& atmtensor_out,
216 const Tensor4& atmtensor_in_orig,
217 const Vector& p_grid_new,
218 const Vector& p_grid_old,
219 const Index& interp_order,
221 const Tensor4* atmtensor_in_pnt;
222 Tensor4 atmtensor_in_copy;
224 if (&atmtensor_in_orig == &atmtensor_out) {
225 atmtensor_in_copy = atmtensor_in_orig;
226 atmtensor_in_pnt = &atmtensor_in_copy;
228 atmtensor_in_pnt = &atmtensor_in_orig;
230 const Tensor4& atmtensor_in = *atmtensor_in_pnt;
233 atmtensor_out.resize(atmtensor_in.nbooks(),
235 atmtensor_in.nrows(),
236 atmtensor_in.ncols());
238 ArrayOfLagrangeLogInterpolation lag_p;
241 Index ing_min, ing_max;
254 (ing_max - ing_min + 1 != p_grid_new.nelem()),
255 "New grid seems not to be sufficiently covered by old grid.\n")
257 for (Index
b = 0;
b < atmtensor_in.nbooks();
b++)
258 for (Index i = 0; i < atmtensor_in.nrows(); i++)
259 for (Index j = 0; j < atmtensor_in.ncols(); j++)
260 reinterp(atmtensor_out(
b, joker, i, j),
261 atmtensor_in(
b, joker, i, j),
278 const Vector& lon_grid,
285 if (lon_grid.empty()) {
288 if (lat_grid.empty())
307 const Vector& p_grid
_U_,
308 const Vector& lat_grid,
309 const Vector& lon_grid,
315 field_out = gfraw_in.
data;
322 const Vector& p_grid,
323 const Vector& lat_grid,
324 const Vector& lon_grid,
333 field_out = gfraw_in.
data;
340 const Vector& p_grid,
341 const Vector& lat_grid,
342 const Vector& lon_grid,
351 field_out = gfraw_in.
data;
358 const Vector& p_grid,
359 const Vector& lat_grid,
360 const Vector& lon_grid,
364 if (!gfraw_in.
nelem()) {
366 out1 <<
" Warning: gfraw_in is empty, proceeding anyway\n";
367 field_out.resize(0, 0, 0, 0);
369 field_out.resize(gfraw_in.
nelem(),
371 gfraw_in[0].data.nrows(),
372 gfraw_in[0].data.ncols());
375 for (Index i = 0; i < gfraw_in.
nelem(); i++) {
378 "p_grid",
"gfield.p_grid", p_grid, gfraw_in[i].get_numeric_grid(0));
381 lat_grid, lon_grid, 1, 2, gfraw_in[i]);
383 field_out(i, joker, joker, joker) = gfraw_in[i].data;
396 if (&gfraw_in_orig == &gfraw_out) {
397 gfraw_in_copy = gfraw_in_orig;
398 gfraw_in_pnt = &gfraw_in_copy;
400 gfraw_in_pnt = &gfraw_in_orig;
408 "Can't expand data because number of Latitudes and Longitudes is greater than 1");
414 if (gfraw_in.
data.nrows() == 1 && gfraw_in.
data.ncols() != 1) {
420 for (Index j = 0; j < gfraw_in.
data.ncols(); j++)
421 gfraw_out.
data(joker, j) = gfraw_in.
data(0, j);
422 }
else if (gfraw_in.
data.nrows() != 1 && gfraw_in.
data.ncols() == 1) {
428 for (Index j = 0; j < gfraw_in.
data.nrows(); j++)
429 gfraw_out.
data(j, joker) = gfraw_in.
data(j, 0);
439 gfraw_out.
data = gfraw_in.
data(0, 0);
452 if (&gfraw_in_orig == &gfraw_out) {
453 gfraw_in_copy = gfraw_in_orig;
454 gfraw_in_pnt = &gfraw_in_copy;
456 gfraw_in_pnt = &gfraw_in_orig;
464 "Can't expand data because number of Latitudes and Longitudes is greater than 1");
472 if (gfraw_in.
data.nrows() == 1 && gfraw_in.
data.ncols() != 1) {
476 gfraw_out.
resize(gfraw_in.
data.npages(), 2, gfraw_in.
data.ncols());
478 for (Index i = 0; i < gfraw_in.
data.npages(); i++)
479 for (Index j = 0; j < gfraw_in.
data.ncols(); j++)
480 gfraw_out.
data(i, joker, j) = gfraw_in.
data(i, 0, j);
481 }
else if (gfraw_in.
data.nrows() != 1 && gfraw_in.
data.ncols() == 1) {
485 gfraw_out.
resize(gfraw_in.
data.npages(), gfraw_in.
data.nrows(), 2);
487 for (Index i = 0; i < gfraw_in.
data.npages(); i++)
488 for (Index j = 0; j < gfraw_in.
data.nrows(); j++)
489 gfraw_out.
data(i, j, joker) = gfraw_in.
data(i, j, 0);
497 gfraw_out.
resize(gfraw_in.
data.npages(), 2, 2);
499 for (Index i = 0; i < gfraw_in.
data.npages(); i++)
500 gfraw_out.
data(i, joker, joker) = gfraw_in.
data(i, 0, 0);
513 if (&gfraw_in_orig == &gfraw_out) {
514 gfraw_in_copy = gfraw_in_orig;
515 gfraw_in_pnt = &gfraw_in_copy;
517 gfraw_in_pnt = &gfraw_in_orig;
525 "Can't expand data because number of Latitudes and Longitudes is greater than 1");
537 if (gfraw_in.
data.nrows() == 1 && gfraw_in.
data.ncols() != 1) {
542 gfraw_in.
data.npages(),
544 gfraw_in.
data.ncols());
546 for (Index k = 0; k < gfraw_in.
data.nbooks(); k++)
547 for (Index i = 0; i < gfraw_in.
data.npages(); i++)
548 for (Index j = 0; j < gfraw_in.
data.ncols(); j++)
549 gfraw_out.
data(k, i, joker, j) = gfraw_in.
data(k, i, 0, j);
550 }
else if (gfraw_in.
data.nrows() != 1 && gfraw_in.
data.ncols() == 1) {
555 gfraw_in.
data.npages(),
556 gfraw_in.
data.nrows(),
559 for (Index k = 0; k < gfraw_in.
data.nbooks(); k++)
560 for (Index i = 0; i < gfraw_in.
data.npages(); i++)
561 for (Index j = 0; j < gfraw_in.
data.nrows(); j++)
562 gfraw_out.
data(k, i, j, joker) = gfraw_in.
data(k, i, j, 0);
570 gfraw_out.
resize(gfraw_in.
data.nbooks(), gfraw_in.
data.npages(), 2, 2);
572 for (Index k = 0; k < gfraw_in.
data.nbooks(); k++)
573 for (Index i = 0; i < gfraw_in.
data.npages(); i++)
574 gfraw_out.
data(k, i, joker, joker) = gfraw_in.
data(k, i, 0, 0);
584 gfraw_out.resize(gfraw_in.
nelem());
586 for (Index i = 0; i < gfraw_in.
nelem(); i++)
611 ArrayOfLagrangeLogInterpolation& lag_p,
615 const Index p_grid_index,
616 ConstVectorView p_grid,
617 const Index& interp_order,
618 const Index& zeropadding,
624 out2 <<
" Interpolation order: " << interp_order <<
"\n";
629 gfraw_out.
set_grid(p_grid_index, Vector{p_grid});
633 if (in_p_grid[0] < p_grid[p_grid.nelem() - 1] ||
634 in_p_grid[in_p_grid.nelem() - 1] > p_grid[0]) {
636 ing_max = ing_min - 1;
640 "Raw field to p_grid",
646 ing_max = p_grid.nelem() - 1;
648 "Raw field to p_grid", in_p_grid, p_grid, interp_order);
650 Index nelem_in_range = ing_max - ing_min + 1;
653 if (nelem_in_range > 0) {
654 lag_p = my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid[Range(ing_min, nelem_in_range)], in_p_grid, interp_order, 0.5);
663 const Vector& p_grid,
666 const Index& interp_order,
667 const Index& zeropadding,
672 if (&gfraw_in_orig == &gfraw_out) {
673 gfraw_in_copy = gfraw_in_orig;
674 gfraw_in_pnt = &gfraw_in_copy;
676 gfraw_in_pnt = &gfraw_in_orig;
680 const Index p_grid_index = 0;
684 p_grid.nelem(), gfraw_in.
data.nrows(), gfraw_in.
data.ncols());
690 ArrayOfLagrangeLogInterpolation lag_p;
693 Index ing_min, ing_max;
708 if (ing_max - ing_min < 0)
710 else if (ing_max - ing_min + 1 != p_grid.nelem()) {
712 for (Index i = 0; i < gfraw_in.
data.nrows(); i++)
713 for (Index j = 0; j < gfraw_in.
data.ncols(); j++) {
718 reinterp(gfraw_out.
data(Range(ing_min, ing_max - ing_min + 1), i, j),
719 gfraw_in.
data(joker, i, j),
724 for (Index i = 0; i < gfraw_in.
data.nrows(); i++)
725 for (Index j = 0; j < gfraw_in.
data.ncols(); j++)
727 gfraw_out.
data(joker, i, j), gfraw_in.
data(joker, i, j), itw, lag_p);
734 const Vector& p_grid,
737 const Index& interp_order,
738 const Index& zeropadding,
743 if (&gfraw_in_orig == &gfraw_out) {
744 gfraw_in_copy = gfraw_in_orig;
745 gfraw_in_pnt = &gfraw_in_copy;
747 gfraw_in_pnt = &gfraw_in_orig;
751 const Index p_grid_index = 1;
756 gfraw_in.
data.nrows(),
757 gfraw_in.
data.ncols());
765 ArrayOfLagrangeLogInterpolation lag_p;
768 Index ing_min, ing_max;
783 if (ing_max - ing_min < 0)
785 else if (ing_max - ing_min + 1 != p_grid.nelem()) {
787 for (Index
b = 0;
b < gfraw_in.
data.nbooks();
b++)
788 for (Index i = 0; i < gfraw_in.
data.nrows(); i++)
789 for (Index j = 0; j < gfraw_in.
data.ncols(); j++) {
794 reinterp(gfraw_out.
data(
b, Range(ing_min, ing_max - ing_min), i, j),
795 gfraw_in.
data(
b, joker, i, j),
800 for (Index
b = 0;
b < gfraw_in.
data.nbooks();
b++)
801 for (Index i = 0; i < gfraw_in.
data.nrows(); i++)
802 for (Index j = 0; j < gfraw_in.
data.ncols(); j++)
803 reinterp(gfraw_out.
data(
b, joker, i, j),
804 gfraw_in.
data(
b, joker, i, j),
813 const Vector& p_grid,
816 const Index& interp_order,
817 const Index& zeropadding,
819 agfraw_out.resize(agfraw_in.
nelem());
821 for (Index i = 0; i < agfraw_in.
nelem(); i++) {
849 ArrayOfLagrangeCyclic0to360Interpolation& lag_lon,
853 const Index lat_grid_index,
854 const Index lon_grid_index,
855 ConstVectorView lat_true,
856 ConstVectorView lon_true,
857 const Index& interp_order,
862 "The new latitude grid is not allowed to be empty.");
864 "The new longitude grid is not allowed to be empty.");
870 "Raw data has to be true 3D data (nlat>1 and nlon>1).");
872 out2 <<
" Interpolation order: " << interp_order <<
"\n";
878 gfraw_out.
set_grid(lat_grid_index, Vector{lat_true});
881 gfraw_out.
set_grid(lon_grid_index, Vector{lon_true});
886 "Raw field to lat_grid, 3D case", in_lat_grid, lat_true, interp_order);
889 lag_lat = my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_true, in_lat_grid, interp_order);
890 lag_lon = my_interp::lagrange_interpolation_list<LagrangeCyclic0to360Interpolation>(lon_true, in_lon_grid, interp_order, 0.5);
898 const Vector& lat_true,
899 const Vector& lon_true,
902 const Index& interp_order,
905 "The new latitude grid is not allowed to be empty.");
907 "The new longitude grid is not allowed to be empty.");
912 if (&gfraw_in_orig == &gfraw_out) {
913 gfraw_in_copy = gfraw_in_orig;
914 gfraw_in_pnt = &gfraw_in_copy;
916 gfraw_in_pnt = &gfraw_in_orig;
920 const Index lat_grid_index = 0;
921 const Index lon_grid_index = 1;
925 "Raw data has to be true 3D data (nlat>1 and nlon>1).\n"
926 "Use GriddedFieldLatLonExpand to convert 1D or 2D data to 3D!\n")
929 gfraw_out.
resize(lat_true.nelem(), lon_true.nelem());
931 ArrayOfLagrangeInterpolation lag_lat;
932 ArrayOfLagrangeCyclic0to360Interpolation lag_lon;
936 const Vector& in_lat_grid =
938 const Vector& in_lon_grid =
941 if (is_lon_cyclic(in_lon_grid)) {
942 for (Index lat = 0; lat < in_lat_grid.nelem(); lat++) {
944 gfraw_in.
data(lat, in_lon_grid.nelem() - 1),
946 "Data values at 0 and 360 degrees for a cyclic longitude grid must match: \n"
947 ,
"Mismatch at latitude index : ", lat,
" ("
948 , in_lat_grid[lat],
" degrees)\n"
949 ,
"Value at 0 degrees longitude : ", gfraw_in.
data(lat, 0)
951 ,
"Value at 360 degrees longitude: "
952 , gfraw_in.
data(lat, in_lon_grid.nelem() - 1),
"\n"
954 , gfraw_in.
data(lat, in_lon_grid.nelem() - 1) -
955 gfraw_in.
data(lat, 0)
974 reinterp(gfraw_out.
data, gfraw_in.
data, itw, lag_lat, lag_lon);
981 const Vector& lat_true,
982 const Vector& lon_true,
985 const Index& interp_order,
988 "The new latitude grid is not allowed to be empty.");
990 "The new longitude grid is not allowed to be empty.");
995 if (&gfraw_in_orig == &gfraw_out) {
996 gfraw_in_copy = gfraw_in_orig;
997 gfraw_in_pnt = &gfraw_in_copy;
999 gfraw_in_pnt = &gfraw_in_orig;
1003 const Index lat_grid_index = 1;
1004 const Index lon_grid_index = 2;
1008 "Raw data has to be true 3D data (nlat>1 and nlon>1).\n"
1009 "Use GriddedFieldLatLonExpand to convert 1D or 2D data to 3D!\n")
1012 gfraw_out.
resize(gfraw_in.
data.npages(), lat_true.nelem(), lon_true.nelem());
1016 ArrayOfLagrangeInterpolation lag_lat;
1017 ArrayOfLagrangeCyclic0to360Interpolation lag_lon;
1022 const Vector& in_lat_grid =
1024 const Vector& in_lon_grid =
1027 if (is_lon_cyclic(in_lon_grid)) {
1028 for (Index g0 = 0; g0 < in_grid0.nelem(); g0++)
1029 for (Index lat = 0; lat < in_lat_grid.nelem(); lat++) {
1031 gfraw_in.
data(g0, lat, 0),
1032 gfraw_in.
data(g0, lat, in_lon_grid.nelem() - 1),
1034 "Data values at 0 and 360 degrees for a cyclic longitude grid must match: \n"
1035 ,
"Mismatch at 1st grid index : " , g0 ,
" (" , in_grid0[g0]
1037 ,
" at latitude index : " , lat ,
" ("
1038 , in_lat_grid[lat] ,
" degrees)\n"
1039 ,
"Value at 0 degrees longitude : " , gfraw_in.
data(g0, lat, 0)
1041 ,
"Value at 360 degrees longitude: "
1042 , gfraw_in.
data(g0, lat, in_lon_grid.nelem() - 1) ,
"\n"
1044 , gfraw_in.
data(g0, lat, in_lon_grid.nelem() - 1) -
1045 gfraw_in.
data(g0, lat, 0)
1064 for (Index i = 0; i < gfraw_in.
data.npages(); i++)
1065 reinterp(gfraw_out.
data(i, joker, joker),
1066 gfraw_in.
data(i, joker, joker),
1076 const Vector& lat_true,
1077 const Vector& lon_true,
1080 const Index& interp_order,
1083 "The new latitude grid is not allowed to be empty.");
1085 "The new longitude grid is not allowed to be empty.");
1090 if (&gfraw_in_orig == &gfraw_out) {
1091 gfraw_in_copy = gfraw_in_orig;
1092 gfraw_in_pnt = &gfraw_in_copy;
1094 gfraw_in_pnt = &gfraw_in_orig;
1098 const Index lat_grid_index = 2;
1099 const Index lon_grid_index = 3;
1103 "Raw data has to be true 3D data (nlat>1 and nlon>1).\n"
1104 "Use GriddedFieldLatLonExpand to convert 1D or 2D data to 3D!\n")
1108 gfraw_in.
data.npages(),
1116 ArrayOfLagrangeInterpolation lag_lat;
1117 ArrayOfLagrangeCyclic0to360Interpolation lag_lon;
1135 const Vector& in_lat_grid =
1137 const Vector& in_lon_grid =
1140 if (is_lon_cyclic(in_lon_grid)) {
1141 for (Index g0 = 0; g0 < in_grid0.nelem(); g0++)
1142 for (Index g1 = 0; g1 < in_grid1.nelem(); g1++)
1143 for (Index lat = 0; lat < in_lat_grid.nelem(); lat++) {
1145 gfraw_in.
data(g0, g1, lat, 0),
1146 gfraw_in.
data(g0, g1, lat, in_lon_grid.nelem() - 1),
1148 "Data values at 0 and 360 degrees for a cyclic longitude grid must match: \n"
1149 "Mismatch at 1st grid index : ", g0,
" ("
1150 , in_grid0[g0],
")\n"
1151 ,
" at 2nd grid index : ", g1,
" ("
1152 , in_grid1[g1],
")\n"
1153 ,
" at latitude index : ", lat,
" ("
1154 , in_lat_grid[lat],
" degrees)\n"
1155 ,
"Value at 0 degrees longitude : "
1156 , gfraw_in.
data(g0, g1, lat, 0),
"\n"
1157 ,
"Value at 360 degrees longitude: "
1158 , gfraw_in.
data(g0, g1, lat, in_lon_grid.nelem() - 1),
"\n"
1160 , gfraw_in.
data(g0, g1, lat, in_lon_grid.nelem() - 1) -
1161 gfraw_in.
data(g0, g1, lat, 0)
1168 for (Index i = 0; i < gfraw_in.
data.nbooks(); i++)
1169 for (Index j = 0; j < gfraw_in.
data.npages(); j++)
1170 reinterp(gfraw_out.
data(i, j, joker, joker),
1171 gfraw_in.
data(i, j, joker, joker),
1181 const Vector& lat_true,
1182 const Vector& lon_true,
1185 const Index& interp_order,
1187 agfraw_out.resize(agfraw_in.
nelem());
1189 for (Index i = 0; i < agfraw_in.
nelem(); i++) {
1219 ArrayOfLagrangeInterpolation& lag_p,
1222 const Index z_grid_index,
1223 ConstVectorView z_grid,
1224 const Index& interp_order,
1225 const Index& zeropadding,
1231 out2 <<
" Interpolation order: " << interp_order <<
"\n";
1236 if (in_z_grid[0] > z_grid[z_grid.nelem() - 1] ||
1237 in_z_grid[in_z_grid.nelem() - 1] < z_grid[0]) {
1239 ing_max = ing_min - 1;
1243 "Raw field to z_grid",
1249 ing_max = z_grid.nelem() - 1;
1251 "Raw field to p_grid", in_z_grid, z_grid, interp_order);
1254 Index nelem_in_range = ing_max - ing_min + 1;
1257 if (nelem_in_range > 0) {
1258 lag_p = my_interp::lagrange_interpolation_list<LagrangeInterpolation>(z_grid[Range(ing_min, nelem_in_range)], in_z_grid, interp_order);
1267 const Vector& p_grid,
1268 const Vector& lat_grid,
1269 const Vector& lon_grid,
1270 const Tensor3& z_field,
1273 const Index& interp_order,
1274 const Index& zeropadding,
1278 z_field.nrows() == lat_grid.nelem()) &&
1279 z_field.ncols() == lon_grid.nelem()),
1280 "*z_field* must be of the same size as *p_grid*, *lat_grid*, and *lon_grid* in *GriddedFieldZToPRegrid*.");
1289 ARTS_USER_ERROR_IF (lat_grid.nelem() != lat_in.nelem() || lon_grid.nelem() != lon_in.nelem(),
1290 "Gridding of field to regrid is bad.\n*GriddedFieldZToPRegrid* requires latitude and longitude to be on the same grid as *z_field*.");
1291 for (Index ii = 0; ii < lat_grid.nelem(); ii++)
1293 "Gridding of field to regrid is bad.\n*GriddedFieldZToPRegrid* requires latitude and longitude of the gridded field to be the same as for *z_field*.");
1294 for (Index ii = 0; ii < lon_grid.nelem(); ii++)
1296 "Gridding of field to regrid is bad.\n*GriddedFieldZToPRegrid* requires latitude and longitude of the gridded field to be the same as for *z_field*.");
1302 if (&gfraw_in_orig == &gfraw_out) {
1303 gfraw_in_copy = gfraw_in_orig;
1304 gfraw_in_pnt = &gfraw_in_copy;
1306 gfraw_in_pnt = &gfraw_in_orig;
1312 gfraw_out.
resize(p_grid.nelem(), lat_grid.nelem(), lon_grid.nelem());
1319 gfraw_out.
data = 0.;
1321 ArrayOfLagrangeInterpolation lag_p;
1324 Index ing_min, ing_max;
1326 for (Index lat_index = 0; lat_index < lat_grid.nelem(); lat_index++) {
1327 for (Index lon_index = 0; lon_index < lon_grid.nelem(); lon_index++) {
1328 const Vector z_out{z_field(joker, lat_index, lon_index)};
1341 if (ing_max - ing_min >= 0) {
1344 if (ing_max - ing_min + 1 != z_out.nelem()) {
1345 r = Range(ing_min, ing_max - ing_min + 1);
1348 reinterp(gfraw_out.
data(r, lat_index, lon_index),
1349 gfraw_in.
data(joker, lat_index, lon_index),
1362 const Index& atmosphere_dim,
1369 "Atmospheric dimension must be 1.")
1371 const Index np = im.nrows();
1372 const Index nf = im.ncols() - 1;
1378 "Cannot extract fields from Matrix.\n"
1379 "*field_names* must have one element less than there are\n"
1385 for (Index f = 0; f < field_names.
nelem(); f++) {
1386 fn_upper = field_names[f];
1389 if (fn_upper !=
"IGNORE") {
1396 Index nf_1 = f_1.size();
1398 for (Index f = 0; f < nf_1; f++) {
1399 field_names_1[f] = field_names[f_1[f]];
1406 af.
set_grid(GFIELD4_FIELD_NAMES, field_names_1);
1408 af.
set_grid(GFIELD4_P_GRID, Vector{im(Range(joker), 0)});
1410 af.
set_grid(GFIELD4_LAT_GRID, Vector());
1411 af.
set_grid(GFIELD4_LON_GRID, Vector());
1413 af.
resize(nf_1, np, 1, 1);
1414 for (Index f = 0; f < nf_1; f++)
1415 af.
data(f, Range(joker), 0, 0) = im(Range(joker), f_1[f] + 1);
1425 const Numeric& value,
1426 const Index& prepend,
1434 if (condensibles.
nelem()) {
1435 const Tensor4& vmrs = af.
data;
1437 Tensor3 condensible_sum(vmrs.npages(), vmrs.nrows(), vmrs.ncols(), 1.);
1438 for (Index
c = 0;
c < condensibles.
nelem();
c++) {
1439 bool species_found =
false;
1440 for (Index i = 0; !species_found && i < species.
nelem(); i++) {
1441 if (species[i] == condensibles[
c]) {
1442 condensible_sum -= vmrs(i, joker, joker, joker);
1443 species_found =
true;
1447 "Condensible species \"", condensibles[
c],
"\" not found "
1450 condensible_sum *= value;
1453 af.
data(0, joker, joker, joker) = condensible_sum;
1455 af.
data(nf - 1, joker, joker, joker) = condensible_sum;
1459 af.
data(0, joker, joker, joker) = value;
1461 af.
data(nf - 1, joker, joker, joker) = value;
1472 const Index& prepend,
1477 ConstVectorView af_p_grid =
1479 ConstVectorView af_lat_grid =
1481 ConstVectorView af_lon_grid =
1489 atm_fields_compact, new_n_fields, name, prepend, verbosity);
1491 const Index insert_pos = (prepend) ? 0 : new_n_fields - 1;
1497 "species p_grid to atm_fields_compact p_grid", sp_p_grid, af_p_grid);
1500 p2gridpos(p_gridpos, sp_p_grid, af_p_grid);
1502 if (sp_lat_grid.nelem() > 1) {
1508 gridpos(lat_gridpos, sp_lat_grid, af_lat_grid);
1510 if (sp_lon_grid.nelem() > 1) {
1515 gridpos(lon_gridpos, sp_lon_grid, af_lon_grid);
1522 af_p_grid.nelem(), af_lat_grid.nelem(), af_lon_grid.nelem());
1523 interp(newfield, itw, species.
data, p_gridpos, lat_gridpos, lon_gridpos);
1525 atm_fields_compact.
data(insert_pos, joker, joker, joker) = newfield;
1528 Tensor3 itw(p_gridpos.
nelem(), lat_gridpos.
nelem(), 4);
1531 Matrix newfield(af_p_grid.nelem(), af_lat_grid.nelem());
1533 newfield, itw, species.
data(joker, joker, 0), p_gridpos, lat_gridpos);
1535 atm_fields_compact.
data(insert_pos, joker, joker, 0) = newfield;
1538 Matrix itw(p_gridpos.
nelem(), 2);
1541 Vector newfield(af_p_grid.nelem());
1542 interp(newfield, itw, species.
data(joker, 0, 0), p_gridpos);
1544 atm_fields_compact.
data(insert_pos, joker, 0, 0) = newfield;
1553 const Numeric& threshold,
1556 Tensor4View afd = atm_fields_compact.
data;
1561 for (Index i = 0; i < afd.nbooks(); i++)
1562 if (atm_fields_compact.
get_string_grid(GFIELD4_FIELD_NAMES)[i] !=
"T" &&
1564 for (Index j = 0; j < afd.npages(); j++)
1565 for (Index k = 0; k < afd.nrows(); k++)
1566 for (Index l = 0; l < afd.ncols(); l++)
1567 if (afd(i, j, k, l) < threshold) afd(i, j, k, l) = 0.0;
1583 sp_name_grid[0] = name;
1585 atm_fields_compact.
set_grid(0, sp_name_grid);
1586 atm_fields_compact.
set_grid(1, Vector{sp_p_grid});
1587 atm_fields_compact.
set_grid(2, Vector{sp_lat_grid});
1588 atm_fields_compact.
set_grid(3, Vector{sp_lon_grid});
1590 atm_fields_compact.
data.resize(
1591 1, sp_p_grid.nelem(), sp_lat_grid.nelem(), sp_lon_grid.nelem());
1593 atm_fields_compact.
data(0, joker, joker, joker) = field.
data;
1602 const Numeric& value,
1603 const Index& prepend,
1606 for (Index i = 0; i < batch_atm_fields_compact.
nelem(); i++) {
1623 const Index& prepend,
1625 const Index nelem = batch_atm_fields_compact.
nelem();
1628 bool failed =
false;
1632#pragma omp parallel for if (!arts_omp_in_parallel() && \
1633 nelem >= arts_omp_get_max_threads())
1634 for (Index i = 0; i < nelem; i++) {
1637 batch_atm_fields_compact[i], name, species, prepend, verbosity);
1638 }
catch (
const std::exception& e) {
1639#pragma omp critical(batch_atm_fields_compactAddSpecies_fail)
1641 fail_msg = e.what();
1655 const Numeric& threshold,
1657 for (Index i = 0; i < batch_atm_fields_compact.
nelem(); i++) {
1659 batch_atm_fields_compact[i], threshold, verbosity);
1667 const Index& atmosphere_dim,
1669 const ArrayOfMatrix& am,
1673 const Index amnelem = am.nelem();
1676 "No elements in atmospheric scenario batch.\n"
1677 "Check, whether any batch atmosphere file has been read!")
1683 batch_atm_fields_compact.resize(amnelem);
1686 bool failed =
false;
1689#pragma omp parallel for if (!arts_omp_in_parallel() && \
1690 amnelem >= arts_omp_get_max_threads())
1691 for (Index i = 0; i < amnelem; ++i) {
1693 if (failed)
continue;
1708 }
catch (
const std::exception& e) {
1709#pragma omp critical(batch_atm_fields_compactFromArrayOfMatrix_fail)
1711 fail_msg = e.what();
1728 Tensor4& particle_bulkprop_field,
1734 const Index& atmosphere_dim,
1736 const Numeric& p_min,
1738 const Index& check_gridnames,
1746 c.get_numeric_grid(GFIELD4_P_GRID),
1747 c.get_numeric_grid(GFIELD4_LAT_GRID),
1748 c.get_numeric_grid(GFIELD4_LON_GRID));
1751 if (check_gridnames == 1) {
1757 const Index nf =
c.get_grid_size(GFIELD4_FIELD_NAMES);
1758 const Index np =
c.get_grid_size(GFIELD4_P_GRID);
1761 Index nlat =
c.get_grid_size(GFIELD4_LAT_GRID);
1762 Index nlon =
c.get_grid_size(GFIELD4_LON_GRID);
1763 if (nlat == 0) nlat++;
1764 if (nlon == 0) nlon++;
1767 p_grid =
c.get_numeric_grid(GFIELD4_P_GRID);
1768 lat_grid =
c.get_numeric_grid(GFIELD4_LAT_GRID);
1769 lon_grid =
c.get_numeric_grid(GFIELD4_LON_GRID);
1773 bool search_toa =
true;
1774 while (search_toa && l > 0) {
1775 if (p_grid[l - 1] < p_min)
1781 "At least one atmospheric level with pressure larger p_min (=",
1783 "is needed, but none is found.")
1784 const Index npn = l + 1;
1785 p_grid = p_grid[Range(0, npn)];
1787 const Index nsa = abs_species.
nelem();
1791 "There must be at least one absorption species.")
1806 const String as_type =
"abs_species";
1807 const String ss_type =
"scat_species";
1811 t_field.resize(npn, nlat, nlon);
1812 for (Index i = 0; i < nf; ++i) {
1813 if (
c.get_string_grid(GFIELD4_FIELD_NAMES)[i] ==
"T") {
1815 "Only one temperature ('T') field allowed,\n"
1816 "but found at least 2.")
1818 t_field =
c.data(i, Range(0, npn), Range(joker), Range(joker));
1822 "One temperature ('T') field required, but none found")
1826 z_field.resize(npn, nlat, nlon);
1827 for (Index i = 0; i < nf; ++i) {
1828 if (
c.get_string_grid(GFIELD4_FIELD_NAMES)[i] ==
"z") {
1830 "Only one altitude ('z') field allowed,\n"
1831 "but found at least 2.")
1833 z_field =
c.data(i, Range(0, npn), Range(joker), Range(joker));
1837 "One altitude ('z') field required, but none found")
1841 vmr_field.resize(nsa, npn, nlat, nlon);
1842 for (Index j = 0; j < nsa; ++j) {
1843 const String as_name = Species::toShortName(abs_species[j][0].Spec());
1848 while (!found && i < nf) {
1850 species_type,
c.get_string_grid(GFIELD4_FIELD_NAMES)[i], delim);
1852 if (species_type == as_type) {
1854 species_name,
c.get_string_grid(GFIELD4_FIELD_NAMES)[i], delim);
1855 if (species_name == as_name) {
1857 vmr_field(j, Range(joker), Range(joker), Range(joker)) =
1858 c.data(i, Range(0, npn), Range(joker), Range(joker));
1864 "No field for absorption species '", as_name,
"' found.")
1868 std::vector<Index> Idx;
1871 for (Index i = 0; i < nf; ++i) {
1873 species_type,
c.get_string_grid(GFIELD4_FIELD_NAMES)[i], delim);
1875 if (species_type == ss_type) {
1880 const Index nsp = Idx.size();
1883 particle_bulkprop_field.resize(nsp, npn, nlat, nlon);
1884 particle_bulkprop_field = NAN;
1885 particle_bulkprop_names.resize(nsp);
1889 for (Index j = 0; j < nsp; ++j) {
1894 scat_type,
c.get_string_grid(GFIELD4_FIELD_NAMES)[Idx[j]], delim);
1897 species_name,
c.get_string_grid(GFIELD4_FIELD_NAMES)[Idx[j]], delim);
1899 particle_bulkprop_field(j, Range(joker), Range(joker), Range(joker)) =
1900 c.data(Idx[j], Range(0, npn), Range(joker), Range(joker));
1902 particle_bulkprop_names[j] = species_name + delim + scat_type;
1908 Index& atmosphere_dim,
1915 out2 <<
" Sets the atmospheric dimensionality to 1.\n";
1916 out3 <<
" atmosphere_dim = 1\n";
1917 out3 <<
" lat_grid is set to be an empty vector\n";
1918 out3 <<
" lon_grid is set to be an empty vector\n";
1927 Index& atmosphere_dim,
1933 out2 <<
" Sets the atmospheric dimensionality to 2.\n";
1934 out3 <<
" atmosphere_dim = 2\n";
1935 out3 <<
" lon_grid is set to be an empty vector\n";
1943 Index& atmosphere_dim,
1950 out2 <<
" Sets the atmospheric dimensionality to 3.\n";
1951 out3 <<
" atmosphere_dim = 3\n";
1965 const Vector& p_grid,
1966 const Vector& lat_grid,
1967 const Vector& lon_grid,
1973 const Vector& nlte_energies,
1974 const Index& atmosphere_dim,
1976 const Index& interp_order,
1977 const Index& vmr_zeropadding,
1978 const Index& vmr_nonegative,
1979 const Index& nlte_when_negative,
1983 const Vector& tfr_p_grid =
1985 const Vector& tfr_lat_grid =
1987 const Vector& tfr_lon_grid =
1989 const Vector& zfr_p_grid =
1991 const Vector& zfr_lat_grid =
1993 const Vector& zfr_lon_grid =
1996 out2 <<
" Interpolation order: " << interp_order <<
"\n";
2007 nlte_field.
vib_energy = nlte_ids.
nelem() == nlte_field_raw.
nelem() ? nlte_energies : Vector(0);
2010 if (atmosphere_dim == 1) {
2012 "Temperature data (T_field) has wrong dimension "
2016 "Altitude data (z_field) has wrong dimension "
2022 temp_gfield3, p_grid, t_field_raw, interp_order, 0, verbosity);
2023 t_field = temp_gfield3.
data;
2026 temp_gfield3, p_grid, z_field_raw, interp_order, 0, verbosity);
2027 z_field = temp_gfield3.
data;
2037 }
catch (
const std::runtime_error& e) {
2040 "Note that you can explicitly set vmr_zeropadding "
2041 "to 1 in the method call.")
2044 vmr_field, p_grid, lat_grid, lon_grid, temp_agfield3, verbosity);
2047 if (nlte_ids.
nelem() == nlte_field_raw.
nelem()) {
2049 temp_agfield3, p_grid, nlte_field_raw, interp_order, 0, verbosity);
2051 nlte_field.
value, p_grid, lat_grid, lon_grid, temp_agfield3, verbosity);
2054 nlte_field.
value.resize(0, 0, 0, 0);
2059 else if (atmosphere_dim == 2) {
2061 "Raw data has wrong dimension (1D). "
2062 "You have to use \n"
2063 "AtmFieldsCalcExpand1D instead of AtmFieldsCalc.");
2066 t_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2067 z_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2069 vmr_field_raw.
nelem(), p_grid.nelem(), lat_grid.nelem(), 1);
2070 if (nlte_ids.
nelem() == nlte_field_raw.
nelem())
2071 nlte_field.
value.resize(
2072 nlte_field_raw.
nelem(), p_grid.nelem(), lat_grid.nelem(), 1);
2074 nlte_field.
value.resize(0, 0, 0, 0);
2081 "Raw temperature to p_grid, 2D case", tfr_p_grid, p_grid, interp_order);
2088 auto lag_p=my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid, tfr_p_grid, interp_order, 0.5);
2089 auto lag_lat=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_grid, tfr_lat_grid, interp_order);
2093 reinterp(t_field(joker, joker, 0),
2094 t_field_raw.
data(joker, joker, 0),
2104 "Raw z to p_grid, 2D case", zfr_p_grid, p_grid, interp_order);
2106 "Raw z to lat_grid, 2D case", zfr_lat_grid, lat_grid, interp_order);
2109 lag_p=my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid, zfr_p_grid, interp_order, 0.5);
2110 lag_lat=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_grid, zfr_lat_grid, interp_order);
2114 reinterp(z_field(joker, joker, 0),
2115 z_field_raw.
data(joker, joker, 0),
2122 for (Index gas_i = 0; gas_i < vmr_field_raw.
nelem(); gas_i++) {
2123 ARTS_USER_ERROR_IF(!(vmr_field_raw[gas_i].get_numeric_grid(GFIELD3_LAT_GRID).nelem() !=
2125 vmr_field_raw[gas_i].get_numeric_grid(GFIELD3_LON_GRID).nelem() ==
2127 "VMR data of the ", gas_i,
" the species has "
2128 "wrong dimension (1D or 3D). \n")
2133 var_string(
"Raw VMR[", gas_i,
"] to p_grid, 2D case"),
2134 vmr_field_raw[gas_i].get_numeric_grid(GFIELD3_P_GRID),
2138 var_string(
"Raw VMR[", gas_i,
"] to lat_grid, 2D case"),
2139 vmr_field_raw[gas_i].get_numeric_grid(GFIELD3_LAT_GRID),
2143 lag_p=my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid, vmr_field_raw[gas_i].get_numeric_grid(GFIELD3_P_GRID), interp_order, 0.5);
2144 lag_lat=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_grid, vmr_field_raw[gas_i].get_numeric_grid(GFIELD3_LAT_GRID), interp_order);
2148 reinterp(vmr_field(gas_i, joker, joker, 0),
2149 vmr_field_raw[gas_i].data(joker, joker, 0),
2156 for (Index qi_i = 0; qi_i < nlte_field_raw.
nelem(); qi_i++) {
2157 ARTS_USER_ERROR_IF (!(nlte_field_raw[qi_i].get_numeric_grid(GFIELD3_LAT_GRID).nelem() !=
2159 nlte_field_raw[qi_i].get_numeric_grid(GFIELD3_LON_GRID).nelem() ==
2161 "NLTE data of the ", qi_i,
" temperature field has "
2162 "wrong dimension (1D or 3D). \n")
2167 var_string(
"Raw NLTE[", qi_i,
"] to p_grid, 2D case"),
2168 nlte_field_raw[qi_i].get_numeric_grid(GFIELD3_P_GRID),
2172 var_string(
"Raw NLTE[", qi_i,
"] to lat_grid, 2D case"),
2173 nlte_field_raw[qi_i].get_numeric_grid(GFIELD3_LAT_GRID),
2177 lag_p=my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid, nlte_field_raw[qi_i].get_numeric_grid(GFIELD3_P_GRID), interp_order, 0.5);
2178 lag_lat=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_grid, nlte_field_raw[qi_i].get_numeric_grid(GFIELD3_LAT_GRID), interp_order);
2182 if (nlte_ids.
nelem() == nlte_field_raw.
nelem())
2183 reinterp(nlte_field.
value(qi_i, joker, joker, 0),
2184 nlte_field_raw[qi_i].data(joker, joker, 0),
2189 nlte_field.
value.resize(0, 0, 0, 0);
2195 else if (atmosphere_dim == 3) {
2197 "Raw data has wrong dimension. You have to use \n"
2198 "AtmFieldsCalcExpand1D instead of AtmFieldsCalc.");
2203 temp_gfield3, lat_grid, lon_grid, t_field_raw, interp_order, verbosity);
2205 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2206 t_field = temp_gfield3.
data;
2209 temp_gfield3, lat_grid, lon_grid, z_field_raw, interp_order, verbosity);
2211 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2212 z_field = temp_gfield3.
data;
2228 }
catch (
const std::runtime_error& e) {
2231 "Note that you can explicitly set vmr_zeropadding "
2232 "to 1 in the method call.")
2235 vmr_field, p_grid, lat_grid, lon_grid, temp_agfield3, verbosity);
2237 if (nlte_field_raw.
nelem()) {
2247 temp_agfield3, p_grid, temp_agfield3, interp_order, 0, verbosity);
2248 }
catch (
const std::runtime_error& e) {
2250 "Note that you can explicitly set vmr_zeropadding "
2251 "to 1 in the method call.")
2254 if (nlte_ids.
nelem() == nlte_field_raw.
nelem())
2256 nlte_field.
value, p_grid, lat_grid, lon_grid, temp_agfield3, verbosity);
2258 nlte_field.
value.resize(0, 0, 0, 0);
2267 if (vmr_nonegative) {
2268 for (Index ib = 0; ib < vmr_field.nbooks(); ib++) {
2269 for (Index ip = 0; ip < vmr_field.npages(); ip++) {
2270 for (Index ir = 0; ir < vmr_field.nrows(); ir++) {
2271 for (Index ic = 0; ic < vmr_field.ncols(); ic++) {
2272 if (vmr_field(ib, ip, ir, ic) < 0) {
2273 vmr_field(ib, ip, ir, ic) = 0;
2282 if (nlte_when_negative != -1) {
2283 if (nlte_field_raw.
nelem()) {
2284 for (Index ib = 0; ib < nlte_field.
value.nbooks(); ib++) {
2285 for (Index ip = 0; ip < nlte_field.
value.npages(); ip++) {
2286 for (Index ir = 0; ir < nlte_field.
value.nrows(); ir++) {
2287 for (Index ic = 0; ic < nlte_field.
value.ncols(); ic++) {
2288 if (nlte_field.
value(ib, ip, ir, ic) < 0) {
2290 nlte_field.
value(ib, ip, ir, ic) =
2291 nlte_when_negative == 1 ? t_field(ip, ir, ic) : 0;
2308 Tensor3& mag_u_field,
2309 Tensor3& mag_v_field,
2310 Tensor3& mag_w_field,
2312 const Tensor3& z_field,
2313 const Vector& lat_grid,
2314 const Vector& lon_grid,
2315 const Vector& refellipsoid,
2319 ARTS_USER_ERROR_IF(z_field.nrows() not_eq lat_grid.nelem(),
"Must have same lat_grid.nelem() [", lat_grid.nelem(),
"] as z_field.nrows() [", z_field.nrows(),
']')
2320 ARTS_USER_ERROR_IF(z_field.ncols() not_eq lon_grid.nelem(),
"Must have same lon_grid.nelem() [", lon_grid.nelem(),
"] as z_field.ncols() [", z_field.ncols(),
']')
2321 ARTS_USER_ERROR_IF(refellipsoid.nelem() not_eq 2 or refellipsoid[0] < 1 or refellipsoid[1] >= 1 or refellipsoid[1] < 0,
"Bad refellipsoid: ", refellipsoid)
2325 mag_u_field = mag.
u;
2326 mag_v_field = mag.v;
2327 mag_w_field = mag.w;
2332 Tensor3& mag_u_field,
2333 Tensor3& mag_v_field,
2334 Tensor3& mag_w_field,
2336 const Vector& p_grid,
2337 const Vector& lat_grid,
2338 const Vector& lon_grid,
2342 const Index& atmosphere_dim,
2344 const Index& interp_order,
2348 const Vector& ufr_p_grid =
2350 const Vector& ufr_lat_grid =
2352 const Vector& ufr_lon_grid =
2354 const Vector& vfr_p_grid =
2356 const Vector& vfr_lat_grid =
2358 const Vector& vfr_lon_grid =
2360 const Vector& wfr_p_grid =
2362 const Vector& wfr_lat_grid =
2364 const Vector& wfr_lon_grid =
2367 out2 <<
" Interpolation order: " << interp_order <<
"\n";
2376 if (atmosphere_dim == 1) {
2378 "Magnetic u field data has wrong dimension (2D or 3D).\n");
2380 "Magnetic v field data has wrong dimension (2D or 3D).\n");
2382 "Magnetic w field data has wrong dimension (2D or 3D).\n");
2387 temp_gfield3, p_grid, mag_u_field_raw, interp_order, 0, verbosity);
2388 mag_u_field = temp_gfield3.
data;
2391 temp_gfield3, p_grid, mag_v_field_raw, interp_order, 0, verbosity);
2392 mag_v_field = temp_gfield3.
data;
2395 temp_gfield3, p_grid, mag_w_field_raw, interp_order, 0, verbosity);
2396 mag_w_field = temp_gfield3.
data;
2400 else if (atmosphere_dim == 2) {
2402 "Raw data has wrong dimension (1D). You have to use \n"
2403 "MagFieldsCalcExpand1D instead of MagFieldsCalc.");
2405 "Raw data has wrong dimension (1D). You have to use \n"
2406 "MagFieldsCalcExpand1D instead of MagFieldsCalc.");
2408 "Raw data has wrong dimension (1D). You have to use \n"
2409 "MagFieldsCalcExpand1D instead of MagFieldsCalc.");
2412 mag_u_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2413 mag_v_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2414 mag_w_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2421 "Raw u field to p_grid, 2D case", ufr_p_grid, p_grid, interp_order);
2428 auto lag_p = my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(ufr_p_grid, p_grid, interp_order, 0.5);
2429 auto lag_lat = my_interp::lagrange_interpolation_list<LagrangeInterpolation>(ufr_lat_grid, lat_grid, interp_order);
2433 reinterp(mag_u_field(joker, joker, 0),
2434 mag_u_field_raw.
data(joker, joker, 0),
2444 "Raw v field to p_grid, 2D case", vfr_p_grid, p_grid, interp_order);
2451 lag_p = my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(vfr_p_grid, p_grid, interp_order, 0.5);
2452 lag_lat = my_interp::lagrange_interpolation_list<LagrangeInterpolation>(vfr_lat_grid, lat_grid, interp_order);
2456 reinterp(mag_v_field(joker, joker, 0),
2457 mag_v_field_raw.
data(joker, joker, 0),
2467 "Raw w field to p_grid, 2D case", wfr_p_grid, p_grid, interp_order);
2474 lag_p = my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(wfr_p_grid, p_grid, interp_order, 0.5);
2475 lag_lat = my_interp::lagrange_interpolation_list<LagrangeInterpolation>(wfr_lat_grid, lat_grid, interp_order);
2479 reinterp(mag_w_field(joker, joker, 0),
2480 mag_w_field_raw.
data(joker, joker, 0),
2488 else if (atmosphere_dim == 3) {
2490 "Raw data has wrong dimension. You have to use \n"
2491 "MagFieldsCalcExpand1D instead of MagFieldsCalc.");
2493 "Raw data has wrong dimension. You have to use \n"
2494 "MagFieldsCalcExpand1D instead of MagFieldsCalc.");
2496 "Raw data has wrong dimension. You have to use \n"
2497 "MagFieldsCalcExpand1D instead of MagFieldsCalc.");
2508 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2509 mag_u_field = temp_gfield3.
data;
2518 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2519 mag_v_field = temp_gfield3.
data;
2528 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2529 mag_w_field = temp_gfield3.
data;
2540 Tensor3& mag_u_field,
2541 Tensor3& mag_v_field,
2542 Tensor3& mag_w_field,
2544 const Vector& lat_grid,
2545 const Vector& lon_grid,
2546 const Tensor3& z_field,
2551 const Index& interp_order,
2552 const Numeric& extrapolation_factor,
2554 const auto nalt = z_field.npages();
2555 const auto nlat = z_field.nrows();
2556 const auto nlon = z_field.ncols();
2559 for (
auto& gf3 : {mag_u_field_raw, mag_v_field_raw, mag_w_field_raw}) {
2561 gf3.get_grid_name(1) not_eq
"Latitude" or
2562 gf3.get_grid_name(2) not_eq
"Longitude",
2564 "Grids must be Altitude, Latitude, Longitude, but are: ",
2565 gf3.get_grid_name(0),
", ", gf3.get_grid_name(1),
", ",
2566 gf3.get_grid_name(2),
'\n')
2572 u, lat_grid, lon_grid, mag_u_field_raw, interp_order, verbosity);
2574 v, lat_grid, lon_grid, mag_v_field_raw, interp_order, verbosity);
2576 w, lat_grid, lon_grid, mag_w_field_raw, interp_order, verbosity);
2579 mag_u_field.resize(nalt, nlat, nlon);
2580 mag_v_field.resize(nalt, nlat, nlon);
2581 mag_w_field.resize(nalt, nlat, nlon);
2582 for (Index ilat = 0; ilat < nlat; ilat++) {
2583 for (Index ilon = 0; ilon < nlon; ilon++) {
2585 u.get_numeric_grid(0),
2586 z_field(joker, ilat, ilon),
2588 extrapolation_factor,
2590 auto lag=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(z_field(joker, ilat, ilon),
u.get_numeric_grid(0), interp_order, extrapolation_factor);
2592 reinterp(mag_u_field(joker, ilat, ilon),
u.data(joker, ilat, ilon), itw, lag);
2595 v.get_numeric_grid(0),
2596 z_field(joker, ilat, ilon),
2598 extrapolation_factor,
2601 lag=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(z_field(joker, ilat, ilon),
v.get_numeric_grid(0), interp_order, extrapolation_factor);
2603 reinterp(mag_v_field(joker, ilat, ilon),
v.data(joker, ilat, ilon), itw, lag);
2606 w.get_numeric_grid(0),
2607 z_field(joker, ilat, ilon),
2609 extrapolation_factor,
2612 lag=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(z_field(joker, ilat, ilon),
w.get_numeric_grid(0), interp_order, extrapolation_factor);
2614 reinterp(mag_w_field(joker, ilat, ilon),
w.data(joker, ilat, ilon), itw, lag);
2621 Tensor3& wind_u_field,
2622 Tensor3& wind_v_field,
2623 Tensor3& wind_w_field,
2625 const Vector& p_grid,
2626 const Vector& lat_grid,
2627 const Vector& lon_grid,
2631 const Index& atmosphere_dim,
2633 const Index& interp_order,
2637 const Vector& ufr_p_grid =
2639 const Vector& ufr_lat_grid =
2641 const Vector& ufr_lon_grid =
2643 const Vector& vfr_p_grid =
2645 const Vector& vfr_lat_grid =
2647 const Vector& vfr_lon_grid =
2649 const Vector& wfr_p_grid =
2651 const Vector& wfr_lat_grid =
2653 const Vector& wfr_lon_grid =
2656 out2 <<
" Interpolation order: " << interp_order <<
"\n";
2665 if (atmosphere_dim == 1) {
2667 "Wind u field data has wrong dimension (2D or 3D).\n");
2669 "Wind v field data has wrong dimension (2D or 3D).\n");
2671 "Wind w field data has wrong dimension (2D or 3D).\n");
2676 temp_gfield3, p_grid, wind_u_field_raw, interp_order, 0, verbosity);
2677 wind_u_field = temp_gfield3.
data;
2680 temp_gfield3, p_grid, wind_v_field_raw, interp_order, 0, verbosity);
2681 wind_v_field = temp_gfield3.
data;
2684 temp_gfield3, p_grid, wind_w_field_raw, interp_order, 0, verbosity);
2685 wind_w_field = temp_gfield3.
data;
2689 else if (atmosphere_dim == 2) {
2691 "Raw data has wrong dimension (1D). You have to use \n"
2692 "WindFieldsCalcExpand1D instead of WindFieldsCalc.");
2694 "Raw data has wrong dimension (1D). You have to use \n"
2695 "WindFieldsCalcExpand1D instead of WindFieldsCalc.");
2697 "Raw data has wrong dimension (1D). You have to use \n"
2698 "WindFieldsCalcExpand1D instead of WindFieldsCalc.");
2701 wind_u_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2702 wind_v_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2703 wind_w_field.resize(p_grid.nelem(), lat_grid.nelem(), 1);
2710 "Raw u field to p_grid, 2D case", ufr_p_grid, p_grid, interp_order);
2717 auto lag_p=my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid, ufr_p_grid, interp_order, 0.5);
2718 auto lag_lat=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_grid, ufr_lat_grid, interp_order);
2722 reinterp(wind_u_field(joker, joker, 0),
2723 wind_u_field_raw.
data(joker, joker, 0),
2733 "Raw v field to p_grid, 2D case", vfr_p_grid, p_grid, interp_order);
2740 lag_p=my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid, vfr_p_grid, interp_order, 0.5);
2741 lag_lat=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_grid, vfr_lat_grid, interp_order);
2745 reinterp(wind_v_field(joker, joker, 0),
2746 wind_v_field_raw.
data(joker, joker, 0),
2756 "Raw w field to p_grid, 2D case", wfr_p_grid, p_grid, interp_order);
2763 lag_p=my_interp::lagrange_interpolation_list<LagrangeLogInterpolation>(p_grid, wfr_p_grid, interp_order, 0.5);
2764 lag_lat=my_interp::lagrange_interpolation_list<LagrangeInterpolation>(lat_grid, wfr_lat_grid, interp_order);
2768 reinterp(wind_w_field(joker, joker, 0),
2769 wind_w_field_raw.
data(joker, joker, 0),
2777 else if (atmosphere_dim == 3) {
2779 "Raw data has wrong dimension. You have to use \n"
2780 "WindFieldsCalcExpand1D instead of WindFieldsCalc.");
2782 "Raw data has wrong dimension. You have to use \n"
2783 "WindFieldsCalcExpand1D instead of WindFieldsCalc.");
2785 "Raw data has wrong dimension. You have to use \n"
2786 "WindFieldsCalcExpand1D instead of WindFieldsCalc.");
2797 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2798 wind_u_field = temp_gfield3.
data;
2807 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2808 wind_v_field = temp_gfield3.
data;
2817 temp_gfield3, p_grid, temp_gfield3, interp_order, 0, verbosity);
2818 wind_w_field = temp_gfield3.
data;
2832 const Vector& p_grid,
2833 const Vector& lat_grid,
2834 const Vector& lon_grid,
2840 const Vector& nlte_energies,
2841 const Index& atmosphere_dim,
2842 const Index& interp_order,
2843 const Index& vmr_zeropadding,
2844 const Index& vmr_nonegative,
2845 const Index& nlte_when_negative,
2851 "This function is intended for 2D and 3D. For 1D, use *AtmFieldsCalc*.");
2855 Tensor3 t_temp, z_temp;
2879 const Index np = p_grid.nelem();
2880 const Index nlat = lat_grid.nelem();
2881 Index nlon = lon_grid.nelem();
2882 if (atmosphere_dim == 2) {
2885 const Index nspecies = vmr_temp.nbooks();
2889 t_field.resize(np, nlat, nlon);
2890 z_field.resize(np, nlat, nlon);
2891 vmr_field.resize(nspecies, np, nlat, nlon);
2892 if (nlte_field_raw.
nelem()) {
2894 nlte_field.
value.resize(nlte_field_raw.
nelem(), np, nlat, nlon);
2895 nlte_field.
levels = nlte_ids;
2901 for (Index ilon = 0; ilon < nlon; ilon++) {
2902 for (Index ilat = 0; ilat < nlat; ilat++) {
2903 for (Index ip = 0; ip < np; ip++) {
2904 t_field(ip, ilat, ilon) = t_temp(ip, 0, 0);
2905 z_field(ip, ilat, ilon) = z_temp(ip, 0, 0);
2906 for (Index is = 0; is < nspecies; is++) {
2907 vmr_field(is, ip, ilat, ilon) = vmr_temp(is, ip, 0, 0);
2909 for (Index is = 0; is < nlte_field_raw.
nelem(); is++) {
2910 nlte_field.
value(is, ip, ilat, ilon) = nlte_temp.
value(is, ip, 0, 0);
2920 Tensor3& mag_v_field,
2921 Tensor3& mag_w_field,
2922 const Vector& p_grid,
2923 const Vector& lat_grid,
2924 const Vector& lon_grid,
2928 const Index& atmosphere_dim,
2929 const Index& interp_order,
2935 "This function is intended for 2D and 3D. For 1D, use *MagFieldsCalc*.");
2939 Tensor3 mag_u_field_temp, mag_v_field_temp, mag_w_field_temp;
2954 const Index np = p_grid.nelem();
2955 const Index nlat = lat_grid.nelem();
2956 Index nlon = lon_grid.nelem();
2957 if (atmosphere_dim == 2) {
2965 mag_u_field.resize(np, nlat, nlon);
2966 mag_v_field.resize(np, nlat, nlon);
2967 mag_w_field.resize(np, nlat, nlon);
2969 for (Index ilon = 0; ilon < nlon; ilon++) {
2970 for (Index ilat = 0; ilat < nlat; ilat++) {
2971 for (Index ip = 0; ip < np; ip++) {
2972 mag_u_field(ip, ilat, ilon) = mag_u_field_temp(ip, 0, 0);
2973 mag_v_field(ip, ilat, ilon) = mag_v_field_temp(ip, 0, 0);
2974 mag_w_field(ip, ilat, ilon) = mag_w_field_temp(ip, 0, 0);
2982 Tensor3& wind_v_field,
2983 Tensor3& wind_w_field,
2984 const Vector& p_grid,
2985 const Vector& lat_grid,
2986 const Vector& lon_grid,
2990 const Index& atmosphere_dim,
2991 const Index& interp_order,
2997 "This function is intended for 2D and 3D. For 1D, use *WindFieldsCalc*.");
3001 Tensor3 wind_u_field_temp, wind_v_field_temp, wind_w_field_temp;
3016 const Index np = p_grid.nelem();
3017 const Index nlat = lat_grid.nelem();
3018 Index nlon = lon_grid.nelem();
3019 if (atmosphere_dim == 2) {
3027 wind_u_field.resize(np, nlat, nlon);
3028 wind_v_field.resize(np, nlat, nlon);
3029 wind_w_field.resize(np, nlat, nlon);
3031 for (Index ilon = 0; ilon < nlon; ilon++) {
3032 for (Index ilat = 0; ilat < nlat; ilat++) {
3033 for (Index ip = 0; ip < np; ip++) {
3034 wind_u_field(ip, ilat, ilon) = wind_u_field_temp(ip, 0, 0);
3035 wind_v_field(ip, ilat, ilon) = wind_v_field_temp(ip, 0, 0);
3036 wind_w_field(ip, ilat, ilon) = wind_w_field_temp(ip, 0, 0);
3046 const Vector& p_grid,
3047 const Vector& lat_grid,
3048 const Vector& lon_grid,
3049 const Index& atmosphere_dim,
3050 const Index& chk_vmr_nan,
3056 const Index np = p_grid.nelem();
3057 const Index nlat = lat_grid.nelem();
3058 const Index nlon =
max(Index(1), lon_grid.nelem());
3059 const Index nspecies = vmr_field.nbooks();
3061 const bool chknan = chk_vmr_nan;
3064 "No use in calling this method for 1D.");
3065 chk_atm_field(
"t_field", t_field, 1, p_grid, Vector(0), Vector(0));
3066 chk_atm_field(
"z_field", z_field, 1, p_grid, Vector(0), Vector(0));
3078 Tensor3 t_temp = t_field, z_temp = z_field;
3079 Tensor4 vmr_temp = vmr_field;
3082 t_field.resize(np, nlat, nlon);
3083 z_field.resize(np, nlat, nlon);
3084 vmr_field.resize(nspecies, np, nlat, nlon);
3086 for (Index ilon = 0; ilon < nlon; ilon++) {
3087 for (Index ilat = 0; ilat < nlat; ilat++) {
3088 for (Index ip = 0; ip < np; ip++) {
3089 t_field(ip, ilat, ilon) = t_temp(ip, 0, 0);
3090 z_field(ip, ilat, ilon) = z_temp(ip, 0, 0);
3091 for (Index is = 0; is < nspecies; is++) {
3092 vmr_field(is, ip, ilat, ilon) = vmr_temp(is, ip, 0, 0);
3109 if (atmosphere_dim == 1) {
3114 "Invalid of *ilat*. It must be >= 0 and less than "
3115 "length of *lat_grid*.");
3117 if (atmosphere_dim == 2) {
3119 vtmp = t_field(joker, ilat, 0);
3120 t_field.resize(t_field.npages(), 1, 1);
3121 t_field(joker, 0, 0) = vtmp;
3122 vtmp = z_field(joker, ilat, 0);
3123 z_field.resize(z_field.npages(), 1, 1);
3124 z_field(joker, 0, 0) = vtmp;
3126 mtmp = vmr_field(joker, joker, ilat, 0);
3127 vmr_field.resize(vmr_field.nbooks(), vmr_field.npages(), 1, 1);
3128 vmr_field(joker, joker, 0, 0) = mtmp;
3129 }
else if (atmosphere_dim == 3) {
3131 "Invalid of *ilon*. It must be >= 0 and less than "
3132 "length of *lon_grid*.");
3134 vtmp = t_field(joker, ilat, ilon);
3135 t_field.resize(t_field.npages(), 1, 1);
3136 t_field(joker, 0, 0) = vtmp;
3137 vtmp = z_field(joker, ilat, ilon);
3138 z_field.resize(z_field.npages(), 1, 1);
3139 z_field(joker, 0, 0) = vtmp;
3141 mtmp = vmr_field(joker, joker, ilat, ilon);
3142 vmr_field.resize(vmr_field.nbooks(), vmr_field.npages(), 1, 1);
3143 vmr_field(joker, joker, 0, 0) = mtmp;
3159 Index& atmfields_checked,
3160 Index& atmgeom_checked,
3161 Index& cloudbox_checked,
3163 const Vector& lat_grid,
3164 const Vector& lon_grid,
3165 const Index& atmosphere_dim,
3167 const Numeric& p_step,
3168 const Index& interp_order,
3181 chk_atm_field(
"t_field", t_field, atmosphere_dim, p_grid, lat_grid, lon_grid);
3184 chk_atm_field(
"z_field", z_field, atmosphere_dim, p_grid, lat_grid, lon_grid);
3207 AtmFieldPRegrid(z_field, z_field, p_grid, p_old, interp_order, verbosity);
3208 AtmFieldPRegrid(t_field, t_field, p_grid, p_old, interp_order, verbosity);
3209 AtmFieldPRegrid(vmr_field, vmr_field, p_grid, p_old, interp_order, verbosity);
3219 Vector& nlte_vibrational_energies,
3227 String tmp_basename = basename;
3228 if (basename.length() && basename[basename.length() - 1] !=
'/')
3229 tmp_basename +=
".";
3232 String file_name = tmp_basename +
"t.xml";
3235 out3 <<
"Temperature field read from file: " << file_name <<
"\n";
3238 file_name = tmp_basename +
"z.xml";
3241 out3 <<
"Altitude field read from file: " << file_name <<
"\n";
3244 vmr_field_raw.resize(0);
3247 for (Index i = 0; i < abs_species.
nelem(); i++) {
3249 file_name = tmp_basename +
3250 String(Species::toShortName(abs_species[i].
Species())) +
".xml";
3254 vmr_field_raw.push_back(vmr_field_data);
3258 file_name, vmr_field_raw[vmr_field_raw.
nelem() - 1], verbosity);
3261 out3 <<
" " << Species::toShortName(abs_species[i].
Species())
3262 <<
" profile read from file: " << file_name <<
"\n";
3266 nlte_field_raw.resize(0);
3267 nlte_quantum_identifiers.resize(0);
3268 nlte_vibrational_energies.resize(0);
3281 String tmp_basename = basename;
3282 if (basename.length() && basename[basename.length() - 1] !=
'/')
3283 tmp_basename +=
".";
3286 String file_name = tmp_basename +
"mag_u.xml";
3289 out3 <<
"Bu field read from file: " << file_name <<
"\n";
3292 file_name = tmp_basename +
"mag_v.xml";
3295 out3 <<
"Bv field read from file: " << file_name <<
"\n";
3298 file_name = tmp_basename +
"mag_w.xml";
3301 out3 <<
"Bw field read from file: " << file_name <<
"\n";
3314 String tmp_basename = basename;
3315 if (basename.length() && basename[basename.length() - 1] !=
'/')
3316 tmp_basename +=
".";
3319 String file_name = tmp_basename +
"wind_u.xml";
3322 out3 <<
"Wind u field read from file: " << file_name <<
"\n";
3325 file_name = tmp_basename +
"wind_v.xml";
3328 out3 <<
"Wind v field read from file: " << file_name <<
"\n";
3331 file_name = tmp_basename +
"wind_w.xml";
3334 out3 <<
"Wind w field read from file: " << file_name <<
"\n";
3344 Vector& nlte_vibrational_energies,
3349 const Index& expect_vibrational_energies,
3353 String tmp_basename = basename;
3354 if (basename.length() && basename[basename.length() - 1] !=
'/')
3355 tmp_basename +=
".";
3358 String file_name = tmp_basename +
"t.xml";
3361 out3 <<
"Temperature field read from file: " << file_name <<
"\n";
3364 file_name = tmp_basename +
"z.xml";
3367 out3 <<
"Altitude field read from file: " << file_name <<
"\n";
3370 vmr_field_raw.resize(0);
3373 for (Index i = 0; i < abs_species.
nelem(); i++) {
3375 file_name = tmp_basename +
3376 String(Species::toShortName(abs_species[i].
Species())) +
".xml";
3380 vmr_field_raw.push_back(vmr_field_data);
3384 file_name, vmr_field_raw[vmr_field_raw.
nelem() - 1], verbosity);
3387 out3 <<
" " << Species::toShortName(abs_species[i].
Species())
3388 <<
" profile read from file: " << file_name <<
"\n";
3392 file_name = tmp_basename +
"nlte.xml";
3395 out3 <<
"NLTE field array read from file: " << file_name <<
"\n";
3398 file_name = tmp_basename +
"qi.xml";
3401 out3 <<
"NLTE identifier array read from file: " << file_name <<
"\n";
3403 if (expect_vibrational_energies) {
3405 file_name = tmp_basename +
"ev.xml";
3408 out3 <<
"NLTE energy levels array read from file: " << file_name <<
"\n";
3411 nlte_vibrational_energies.resize(0);
3415 (nlte_field_raw.
nelem() != nlte_vibrational_energies.nelem() and
3416 0 != nlte_vibrational_energies.nelem()),
3417 "The quantum identifers and the NLTE temperature fields\n"
3418 "are of different lengths. This should not be the case.\n"
3419 "please check the qi.xml and t_nlte.xml files under\n",
3421 "to correct this error.\n")
3426 const Vector& lat_grid,
3427 const Vector& lon_grid,
3429 const Index& interp_order,
3430 const Index& set_lowest_altitude_to_zero,
3434 out3 <<
"Reading GriddedField2 surface altitude from " << filename <<
"\n";
3438 out3 <<
"Surface altitude field interpolated back to lat_grid and lon_grid\n";
3445 z_surface = z_surface_field.
data;
3446 if (set_lowest_altitude_to_zero) {
3447 z_surface -=
min(z_surface);
3453 const Vector& lat_grid,
3454 const Vector& lon_grid,
3455 const Numeric& altitude,
3458 out3 <<
"Setting surface to constant altitude of " << altitude <<
" m\n";
3459 z_surface = Matrix(lat_grid.nelem() ? lat_grid.nelem() : 1, lon_grid.nelem() ? lon_grid.nelem() : 1, altitude);
3464 const Index& atmosphere_dim,
3465 const Vector& p_grid,
3466 const Vector& lat_grid,
3467 const Vector& lon_grid,
3468 const Tensor3& z_field,
3469 const Vector& rtp_pos,
3470 const Tensor3& field,
3498 out3 <<
" Result = " << outvalue <<
"\n";
3504 Index& atmfields_checked,
3505 Index& atmgeom_checked,
3506 Index& cloudbox_checked,
3508 const Vector& p_grid_old,
3516 "The old and new grids (p_grid and p_grid_old) are not allowed\n"
3517 "to be identical (pointing to same memory space).\n"
3518 "But they are doing in your case.")
3523 atmfields_checked = 0;
3524 atmgeom_checked = 0;
3525 cloudbox_checked = 0;
3530 // Nothing to do if nfill=0
3532 // Allocate new size for p_grid
3533 const Index n0 = p_grid_old.nelem();
3534 p_grid.resize((n0 - 1) * (1 + nfill) + 1);
3537 p_grid[0] = p_grid_old[0];
3539 for (Index i = 1; i < n0; i++) {
3542 pnew, 2 + nfill, p_grid_old[i - 1], p_grid_old[i], verbosity);
3543 for (Index j = 1; j < nfill + 2; j++) {
3545 p_grid[iout] = pnew[j];
3551/* Workspace method: Doxygen documentation will be auto-generated */
3552void p_gridRefine( // WS Output:
3554 Index& atmfields_checked,
3555 Index& atmgeom_checked,
3556 Index& cloudbox_checked,
3558 const Vector& p_grid_old,
3559 // Control Parameters:
3560 const Numeric& p_step10,
3562 // Check that p_grid and p_grid_old are not the same variable (pointing to the
3563 // same memory space). this as p_grid will be overwritten, but we will need
3564 // both data later on for data regridding.
3565 ARTS_USER_ERROR_IF (&p_grid == &p_grid_old,
3566 "The old and
new grids (p_grid and p_grid_old) are not allowed\n
"
3567 "to be identical (pointing to same memory space).\n
"
3568 "But they are doing in your
case.
")
3570 // as we manipoulate the overall vertical grid (but not simultaneously the
3571 // atmospheric fields), we reset all atmfields related checked WSV to
3572 // unchecked, forcing the user to do the checks again.
3573 atmfields_checked = 0;
3574 atmgeom_checked = 0;
3575 cloudbox_checked = 0;
3577 // Check the keyword argument:
3578 ARTS_USER_ERROR_IF (p_step10 <= 0,
3579 "The keyword argument p_step must be >0.
")
3581 // For consistency with other code around arts (e.g., correlation
3582 // lengths in atmlab), p_step is given as log10(p[Pa]). However, we
3583 // convert it here to the natural log:
3584 const Numeric p_step = log(pow(10.0, p_step10));
3586 // Now starting modification of p_grid
3588 // We will need the log of the pressure grid:
3589 Vector log_p_old(p_grid_old.nelem());
3590 transform(log_p_old, log, p_grid_old);
3592 // const Numeric epsilon = 0.01 * p_step; // This is the epsilon that
3593 // // we use for comparing p grid spacings.
3595 // Construct p_grid (new)
3596 // ----------------------
3598 ArrayOfNumeric log_p_new; // We take log_p_new as an array of Numeric, so
3599 // that we can easily build it up by appending new
3600 // elements to the end.
3602 // Check whether there are pressure levels that are further apart
3603 // (in log(p)) than p_step, and insert additional levels if
3606 log_p_new.push_back(log_p_old[0]);
3608 for (Index i = 1; i < log_p_old.nelem(); ++i) {
3610 log_p_old[i - 1] - log_p_old[i]; // The grid is descending.
3612 const Numeric dp_by_p_step = dp / p_step;
3613 // cout << "dp_by_p_step:
" << dp_by_p_step << "\n
";
3615 // How many times does p_step fit into dp?
3616 const Index n = (Index)ceil(dp_by_p_step);
3617 // n is the number of intervals that we want to have in the
3618 // new grid. The number of additional points to insert is
3619 // n-1. But we have to insert the original point as well.
3620 // cout << n << "\n
";
3622 const Numeric ddp = dp / (Numeric)n;
3623 // cout << "ddp:
" << ddp << "\n
";
3625 for (Index j = 1; j <= n; ++j)
3626 log_p_new.push_back(log_p_old[i - 1] - (Numeric)j * ddp);
3629 // Copy ArrayOfNumeric to proper vector.
3630 Vector log_p(log_p_new.nelem());
3631 for (Index i = 0; i < log_p_new.nelem(); ++i) log_p[i] = log_p_new[i];
3633 // Copy the new grid to abs_p, removing the log:
3634 p_grid.resize(log_p.nelem());
3635 transform(p_grid, exp, log_p);
3638/* Workspace method: Doxygen documentation will be auto-generated */
3639void p_gridFromZRaw( //WS Output
3642 const GriddedField3& z_field_raw,
3643 const Index& no_negZ,
3645 // original version excludes negative z. not clear, why this is. maybe is
3646 // currently a convention somehwere in ARTS (DOIT?). negative z seem, however,
3647 // to work fine for clear-sky cases. so we make the negative z exclude an
3648 // option (for consistency until unclarities solved, default: do exclude)
3649 Vector p_grid_raw = z_field_raw.get_numeric_grid(GFIELD3_P_GRID);
3652 if (is_increasing(z_field_raw.data(joker, 0, 0))) {
3655 while (z_field_raw.data(i, 0, 0) < 0.0) i++;
3657 p_grid = p_grid_raw[Range(i, joker)];
3658 } else if (is_decreasing(z_field_raw.data(joker, 0, 0))) {
3659 i = z_field_raw.data.npages() - 1;
3661 while (z_field_raw.data(i, 0, 0) < 0.0) i--;
3663 p_grid = reverse(p_grid_raw);
3665 ARTS_USER_ERROR ("z_field_raw needs to be monotonous, but
this is not the
case.\n
")
3669/* Workspace method: Doxygen documentation will be auto-generated */
3670void lat_gridFromZRaw( //WS Output
3673 const GriddedField3& z_field_raw,
3675 lat_grid = z_field_raw.get_numeric_grid(GFIELD3_LAT_GRID);
3678/* Workspace method: Doxygen documentation will be auto-generated */
3679void lon_gridFromZRaw( //WS Output
3682 const GriddedField3& z_field_raw,
3684 lon_grid = z_field_raw.get_numeric_grid(GFIELD3_LON_GRID);
3687/* Workspace method: Doxygen documentation will be auto-generated */
3688void atm_gridsFromZRaw( //WS Output
3693 const GriddedField3& z_field_raw,
3694 const Index& no_negZ,
3695 const Verbosity& v) {
3696 p_gridFromZRaw(p_grid, z_field_raw, no_negZ, v);
3697 lat_gridFromZRaw(lat_grid, z_field_raw, v);
3698 lon_gridFromZRaw(lon_grid, z_field_raw, v);
3701/* Workspace method: Doxygen documentation will be auto-generated */
3702void lat_gridFromRawField( //WS Output
3705 const GriddedField3& field_raw,
3707 lat_grid = field_raw.get_numeric_grid(GFIELD3_LAT_GRID);
3710/* Workspace method: Doxygen documentation will be auto-generated */
3711void lon_gridFromRawField( //WS Output
3714 const GriddedField3& field_raw,
3716 lon_grid = field_raw.get_numeric_grid(GFIELD3_LON_GRID);
3719/* Workspace method: Doxygen documentation will be auto-generated */
3720void wind_u_fieldIncludePlanetRotation(Tensor3& wind_u_field,
3721 const Index& atmosphere_dim,
3722 const Vector& p_grid,
3723 const Vector& lat_grid,
3724 const Vector& lon_grid,
3725 const Vector& refellipsoid,
3726 const Tensor3& z_field,
3727 const Numeric& planet_rotation_period,
3729 ARTS_USER_ERROR_IF (atmosphere_dim < 3,
3730 "No need to use
this method
for 1D and 2D.
");
3732 const Index np = p_grid.nelem();
3733 const Index na = lat_grid.nelem();
3734 const Index no = lon_grid.nelem();
3736 chk_atm_field("z_field
", z_field, atmosphere_dim, p_grid, lat_grid, lon_grid);
3737 if (wind_u_field.npages() > 0) {
3738 chk_atm_field("wind_u_field
",
3745 wind_u_field.resize(np, na, no);
3749 const Numeric k1 = 2 * Constant::pi / planet_rotation_period;
3751 for (Index a = 0; a < na; a++) {
3752 const Numeric k2 = k1 * Conversion::cosd(lat_grid[a]);
3753 const Numeric re = refell2r(refellipsoid, lat_grid[a]);
3755 for (Index o = 0; o < no; o++) {
3756 for (Index p = 0; p < np; p++) {
3757 wind_u_field(p, a, o) += k2 * (re + z_field(p, a, o));
3763// A small help function
3764void z2g(Numeric& g, const Numeric& r, const Numeric& g0, const Numeric& z) {
3765 const Numeric x = r / (r + z);
3769/* Workspace method: Doxygen documentation will be auto-generated */
3770void z_fieldFromHSE(Workspace& ws,
3772 const Index& atmosphere_dim,
3773 const Vector& p_grid,
3774 const Vector& lat_grid,
3775 const Vector& lon_grid,
3776 const Vector& lat_true,
3777 const Vector& lon_true,
3778 const ArrayOfArrayOfSpeciesTag& abs_species,
3779 const Tensor3& t_field,
3780 const Tensor4& vmr_field,
3781 const Vector& refellipsoid,
3782 const Matrix& z_surface,
3783 const Index& atmfields_checked,
3784 const Agenda& g0_agenda,
3785 const Numeric& molarmass_dry_air,
3786 const Numeric& p_hse,
3787 const Numeric& z_hse_accuracy,
3788 const Verbosity& verbosity) {
3789 ARTS_USER_ERROR_IF (atmfields_checked != 1,
3790 "The atmospheric fields must be flagged to have
"
3791 "passed
a consistency check (atmfields_checked=1).
");
3793 // Some general variables
3794 const Index np = p_grid.nelem();
3795 const Index nlat = t_field.nrows();
3796 const Index nlon = t_field.ncols();
3798 const Index firstH2O = find_first_species(
3799 abs_species, Species::fromShortName("H2O
"));
3803 out1 << "No water vapour tag group in *abs_species*.\n
"
3804 << "Be aware that
this leads to significant variations in atmospheres\n
"
3805 << "that contain considerable amounts of water vapour (e.g. Earth)!\n
";
3808 ARTS_USER_ERROR_IF (p_hse > p_grid[0] || p_hse < p_grid[np - 1],
3809 "The value of *p_hse* must be inside the range of *p_grid*:
"
3810 " p_hse =
", p_hse, " Pa\n
"
3811 " p_grid =
", p_grid[np-1],
3812 " -
", p_grid[0], " Pa\n
")
3814 ARTS_USER_ERROR_IF (z_hse_accuracy <= 0,
3815 "The value of *z_hse_accuracy* must be > 0.
");
3817 chk_latlon_true(atmosphere_dim, lat_grid, lat_true, lon_true);
3819 // Determine interpolation weights for p_hse
3821 ArrayOfGridPos gp(1);
3823 p2gridpos(gp, p_grid, Vector(1, p_hse));
3824 interpweights(itw, gp);
3826 // // Molecular weight of water vapour
3827 const Numeric mw = 18.016;
3829 // mw/molarmass_dry_air matches eps in Eq. 3.14 in Wallace&Hobbs:
3830 const Numeric k = 1 - mw / molarmass_dry_air;
3832 // Gas constant for 1kg dry air:
3833 const Numeric rd = 1e3 * GAS_CONSTANT / molarmass_dry_air;
3837 for (Index ilat = 0; ilat < nlat; ilat++) {
3838 // The reference ellipsoid is already adjusted to internal 1D or 2D
3839 // views, and lat_grid is the relevant grid for *refellipsoid*, also
3840 // for 2D. On the other hand, extraction of g0 requires that the true
3841 // position is determined.
3843 // Radius of reference ellipsoid
3845 if (atmosphere_dim == 1) {
3846 re = refellipsoid[0];
3848 re = refell2r(refellipsoid, lat_grid[ilat]);
3851 for (Index ilon = 0; ilon < nlon; ilon++) {
3852 // Determine true latitude and longitude
3854 Vector pos(atmosphere_dim); // pos[0] can be a dummy value
3855 if (atmosphere_dim >= 2) {
3856 pos[1] = lat_grid[ilat];
3857 if (atmosphere_dim == 3) {
3858 pos[2] = lon_grid[ilon];
3862 lat, lon, atmosphere_dim, lat_grid, lat_true, lon_true, pos);
3866 g0_agendaExecute(ws, g0, lat, lon, g0_agenda);
3868 // Determine altitude for p_hse
3870 interp(z_hse, itw, z_field(joker, ilat, ilon), gp);
3872 Numeric z_acc = 2 * z_hse_accuracy;
3874 while (z_acc > z_hse_accuracy) {
3878 for (Index ip = 0; ip < np - 1; ip++) {
3879 // Calculate average g
3881 z2g(g2, re, g0, z_field(ip, ilat, ilon));
3883 const Numeric g1 = g2;
3884 z2g(g2, re, g0, z_field(ip + 1, ilat, ilon));
3886 const Numeric g = (g1 + g2) / 2.0;
3888 //Average water vapour VMR
3893 hm = 0.5 * (vmr_field(firstH2O, ip, ilat, ilon) +
3894 vmr_field(firstH2O, ip + 1, ilat, ilon));
3897 // Average virtual temperature (no liquid water)
3898 // (cf. e.g. Eq. 3.16 in Wallace&Hobbs)
3900 (1 / (2 * (1 - hm * k))) *
3901 (t_field(ip, ilat, ilon) + t_field(ip + 1, ilat, ilon));
3903 // Change in vertical altitude from ip to ip+1
3904 // (cf. e.g. Eq. 3.24 in Wallace&Hobbs)
3905 const Numeric dz = rd * (tv / g) * log(p_grid[ip] / p_grid[ip + 1]);
3908 Numeric znew = z_field(ip, ilat, ilon) + dz;
3909 z_acc = max(z_acc, fabs(znew - z_field(ip + 1, ilat, ilon)));
3910 z_field(ip + 1, ilat, ilon) = znew;
3915 interp(z_tmp, itw, z_field(joker, ilat, ilon), gp);
3916 z_field(joker, ilat, ilon) -= z_tmp[0] - z_hse[0];
3921 // Check that there is no gap between the surface and lowest pressure
3923 // (This code is copied from *basics_checkedCalc*. Make this to an internal
3924 // function if used in more places.)
3925 for (Index row = 0; row < z_surface.nrows(); row++) {
3926 for (Index col = 0; col < z_surface.ncols(); col++) {
3927 ARTS_USER_ERROR_IF (z_surface(row, col) < z_field(0, row, col) ||
3928 z_surface(row, col) >= z_field(z_field.npages() - 1, row, col),
3929 "The surface altitude (*z_surface*) cannot be outside
"
3930 "of the altitudes in *z_field*.
")
3935/* Workspace method: Doxygen documentation will be auto-generated */
3936void vmr_fieldSetConstant(Tensor4& vmr_field,
3937 const ArrayOfArrayOfSpeciesTag& abs_species,
3938 const String& species,
3939 const Numeric& vmr_value,
3942 chk_if_in_range("vmr_value
", vmr_value, 0, 1);
3944 ARTS_USER_ERROR_IF (abs_species.nelem() != vmr_field.nbooks(),
3945 "Size of *vmr_field* and length of *abs_species*
do not agree.
");
3947 // Find position for this species.
3948 const ArrayOfSpeciesTag tag(species);
3949 Index si = chk_contains("species
", abs_species, tag);
3952 vmr_field(si, joker, joker, joker) = vmr_value;
3955/* Workspace method: Doxygen documentation will be auto-generated */
3956void vmr_fieldSetAllConstant(Tensor4& vmr_field,
3957 const ArrayOfArrayOfSpeciesTag& abs_species,
3958 const Vector& vmr_values,
3959 const Verbosity& verbosity) {
3962 const Index nspecies = abs_species.nelem();
3964 ARTS_USER_ERROR_IF (vmr_values.nelem() not_eq nspecies,
3965 "Not same number of vmr_values as abs_species.
");
3967 out3 << "Setting all
" << nspecies << " species to constant VMR\n
";
3969 for (Index i = 0; i < nspecies; i++) {
3970 const ArrayOfSpeciesTag& a_abs_species = abs_species[i];
3971 const String species_tag_name = a_abs_species.Name();
3972 vmr_fieldSetConstant(
3973 vmr_field, abs_species, species_tag_name, vmr_values[i], verbosity);
3978/* Workspace method: Doxygen documentation will be auto-generated */
3979void vmr_fieldSetRh(Workspace& ws,
3981 const ArrayOfArrayOfSpeciesTag& abs_species,
3982 const Tensor3& t_field,
3983 const Vector& p_grid,
3984 const Agenda& water_p_eq_agenda,
3986 const Numeric& vmr_threshold,
3990 const Index ih2o = find_first_species(abs_species, Species::fromShortName("H2O
"));
3991 ARTS_USER_ERROR_IF (ih2o < 0, "There is no H2O species in *abs_species*.
")
3993 // Calculate partial pressure matching selected RH
3995 water_p_eq_agendaExecute(ws, p_rh, t_field, water_p_eq_agenda);
3999 for (Index p=0; p<vmr_field.npages(); ++p) {
4000 for (Index a=0; a<vmr_field.nrows(); ++a) {
4001 for (Index o=0; o<vmr_field.ncols(); ++o) {
4002 if (vmr_field(ih2o, p, a, o) > vmr_threshold) {
4003 vmr_field(ih2o, p, a, o) = p_rh(p, a, o) / p_grid[p];
4010/* Workspace method: Doxygen documentation will be auto-generated */
4011void nlte_fieldLteExternalPartitionFunction(
4013 EnergyLevelMap& nlte_field,
4014 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
4015 const ArrayOfQuantumIdentifier& nlte_quantum_identifiers,
4016 const Tensor3& t_field,
4017 const Verbosity& verbosity) {
4021 const Index nn = nlte_quantum_identifiers.nelem(), np = t_field.npages(),
4022 nlat = t_field.nrows(), nlon = t_field.ncols();
4023 if (nn == 0) return;
4025 Tensor4 nlte_tensor4(nn, np, nlat, nlon);
4027 ArrayOfIndex checked(nn, 0);
4029 for (Index in = 0; in < nn; in++) {
4030 const QuantumIdentifier& qi = nlte_quantum_identifiers[in];
4031 Tensor3View lte = nlte_tensor4(in, joker, joker, joker);
4033 for (auto& abs_lines : abs_lines_per_species) {
4034 for (auto& band : abs_lines) {
4035 for (Index k=0; k<band.NumLines(); k++) {
4036 const Quantum::Number::StateMatch lt(qi, band.lines[k].localquanta, band.quantumidentity);
4037 if (lt == Quantum::Number::StateMatchType::Level and lt.low) {
4038 band.population = Absorption::PopulationType::NLTE;
4040 if (not checked[in]) {
4043 for (Index ip = 0; ip < np; ip++) {
4044 for (Index ilat = 0; ilat < nlat; ilat++) {
4045 for (Index ilon = 0; ilon < nlon; ilon++) {
4046 lte(ip, ilat, ilon) =
4047 boltzman_factor(t_field(ip, ilat, ilon), band.lines[k].E0) *
4048 band.lines[k].glow / single_partition_function(
4049 t_field(ip, ilat, ilon), band.Isotopologue());
4056 if (lt == Quantum::Number::StateMatchType::Level and lt.upp) {
4057 band.population = Absorption::PopulationType::NLTE;
4059 if (not checked[in]) {
4061 for (Index ip = 0; ip < np; ip++) {
4062 for (Index ilat = 0; ilat < nlat; ilat++) {
4063 for (Index ilon = 0; ilon < nlon; ilon++) {
4064 lte(ip, ilat, ilon) =
4065 boltzman_factor(t_field(ip, ilat, ilon), band.lines[k].E0 + h*band.lines[k].F0) *
4066 band.lines[k].gupp / single_partition_function(
4067 t_field(ip, ilat, ilon), band.Isotopologue());
4078 for (Index in = 0; in < nn; in++) {
4079 if (not checked[in]) {
4080 out2 << "Did not find match among lines
for:
"
4081 << nlte_quantum_identifiers[in] << "\n
";
4085 nlte_field = EnergyLevelMap(nlte_tensor4, nlte_quantum_identifiers);
4088/* Workspace method: Doxygen documentation will be auto-generated */
4089void nlte_fieldLteInternalPartitionFunction(
4091 EnergyLevelMap& nlte_field,
4092 ArrayOfArrayOfAbsorptionLines& abs_lines_per_species,
4093 const ArrayOfQuantumIdentifier& nlte_quantum_identifiers,
4094 const Tensor3& t_field,
4095 const Verbosity& verbosity) {
4099 const Index nn = nlte_quantum_identifiers.nelem(), np = t_field.npages(),
4100 nlat = t_field.nrows(), nlon = t_field.ncols();
4101 if (nn == 0) return;
4103 // Find where they are positioned and how many different molecules there are for the NLTE fields
4104 ArrayOfIndex part_fun_pos(nn, 0);
4106 for (Index in = 1; in < nn; in++) {
4108 for (Index ix = 0; ix < in; ix++) {
4109 if (nlte_quantum_identifiers[in].Species() ==
4110 nlte_quantum_identifiers[ix].Species() and
4111 nlte_quantum_identifiers[in].Isotopologue() ==
4112 nlte_quantum_identifiers[ix].Isotopologue()) {
4113 part_fun_pos[in] = part_fun_pos[ix];
4119 part_fun_pos[in] = x;
4124 Tensor4 part_fun(x, np, nlat, nlon, 0.0);
4125 Tensor4 nlte_tensor4(nn, np, nlat, nlon, 0);
4127 ArrayOfIndex checked(nn, 0);
4129 for (Index in = 0; in < nn; in++) {
4130 const QuantumIdentifier& qi = nlte_quantum_identifiers[in];
4131 Tensor3View lte = nlte_tensor4(in, joker, joker, joker);
4133 for (auto& abs_lines : abs_lines_per_species) {
4134 for (auto& band : abs_lines) {
4135 for (Index k=0; k<band.NumLines(); k++) {
4136 const Quantum::Number::StateMatch lt(qi, band.lines[k].localquanta, band.quantumidentity);
4137 if (lt == Quantum::Number::StateMatchType::Level and lt.low) {
4138 band.population = Absorption::PopulationType::NLTE;
4140 if (not checked[in]) {
4143 for (Index ip = 0; ip < np; ip++) {
4144 for (Index ilat = 0; ilat < nlat; ilat++) {
4145 for (Index ilon = 0; ilon < nlon; ilon++) {
4146 lte(ip, ilat, ilon) =
4147 boltzman_factor(t_field(ip, ilat, ilon), band.lines[k].E0) *
4149 part_fun(part_fun_pos[in], ip, ilat, ilon) +=
4150 lte(ip, ilat, ilon);
4157 if (lt == Quantum::Number::StateMatchType::Level and lt.upp) {
4158 band.population = Absorption::PopulationType::NLTE;
4160 if (not checked[in]) {
4163 for (Index ip = 0; ip < np; ip++) {
4164 for (Index ilat = 0; ilat < nlat; ilat++) {
4165 for (Index ilon = 0; ilon < nlon; ilon++) {
4166 lte(ip, ilat, ilon) =
4167 boltzman_factor(t_field(ip, ilat, ilon),
4168 band.lines[k].E0 + h * band.lines[k].F0) *
4170 part_fun(part_fun_pos[in], ip, ilat, ilon) +=
4171 lte(ip, ilat, ilon);
4182 for (Index in = 0; in < nn; in++) {
4183 if (not checked[in]) {
4184 out2 << "Did not find match among lines
for:
"
4185 << nlte_quantum_identifiers[in] << "\n
";
4189 for (Index in = 0; in < nn; in++) {
4191 nlte_tensor4(in, joker, joker, joker) /=
4192 part_fun(part_fun_pos[in], joker, joker, joker);
4196 nlte_field = EnergyLevelMap(nlte_tensor4, nlte_quantum_identifiers);
Declarations required for the calculation of absorption coefficients.
Declarations for agendas.
base max(const Array< base > &x)
Max function.
base min(const Array< base > &x)
Min function.
The global header file for ARTS.
Constants of physical expressions as constexpr.
Index nelem() const ARTS_NOEXCEPT
void resize(const GriddedField2 &gf)
Make this GriddedField2 the same size as the given one.
bool checksize() const final
Consistency check.
void resize(const GriddedField3 &gf)
Make this GriddedField3 the same size as the given one.
bool checksize() const final
Consistency check.
void resize(const GriddedField4 &gf)
Make this GriddedField4 the same size as the given one.
Index get_grid_size(Index i) const
Get the size of a grid.
const ArrayOfString & get_string_grid(Index i) const
Get a string grid.
void set_grid_name(Index i, const String &s)
Set grid name.
void set_grid(Index i, const Vector &g)
Set a numeric grid.
const Vector & get_numeric_grid(Index i) const
Get a numeric grid.
const String & get_grid_name(Index i) const
Get grid name.
void toupper()
Convert to upper case.
void parse_atmcompact_speciesname(String &species_name, const String &field_name, const String &delim)
void parse_atmcompact_speciestype(String &species_type, const String &field_name, const String &delim)
void parse_atmcompact_scattype(String &scat_type, const String &field_name, const String &delim)
Internal cloudbox functions.
#define ARTS_ASSERT(condition,...)
#define ARTS_USER_ERROR(...)
std::string var_string(Args &&... args)
#define ARTS_USER_ERROR_IF(condition,...)
Implementation of gridded fields.
void gridpos(ArrayOfGridPos &gp, ConstVectorView old_grid, ConstVectorView new_grid, const Numeric &extpolfac)
Set up a grid position Array.
void interpweights(VectorView itw, const GridPos &tc)
Red 1D interpolation weights.
Numeric interp(ConstVectorView itw, ConstVectorView a, const GridPos &tc)
Red 1D Interpolate.
Header file for interpolation.cc.
Constains various line scaling functions.
void AtmosphereSet2D(Index &atmosphere_dim, Vector &lon_grid, const Verbosity &verbosity)
WORKSPACE METHOD: AtmosphereSet2D.
void MagRawRead(GriddedField3 &mag_u_field_raw, GriddedField3 &mag_v_field_raw, GriddedField3 &mag_w_field_raw, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: MagRawRead.
void GriddedFieldZToPRegridHelper(Index &ing_min, Index &ing_max, ArrayOfLagrangeInterpolation &lag_p, Matrix &itw, const GriddedField &gfraw_in, const Index z_grid_index, ConstVectorView z_grid, const Index &interp_order, const Index &zeropadding, const Verbosity &verbosity)
Calculate grid positions and interpolations weights for GriddedFieldZToPRegrid.
void atm_fields_compactCreateFromField(GriddedField4 &atm_fields_compact, const String &name, const GriddedField3 &field, const Verbosity &)
WORKSPACE METHOD: atm_fields_compactCreateFromField.
void AtmRawRead(GriddedField3 &t_field_raw, GriddedField3 &z_field_raw, ArrayOfGriddedField3 &vmr_field_raw, ArrayOfGriddedField3 &nlte_field_raw, ArrayOfQuantumIdentifier &nlte_quantum_identifiers, Vector &nlte_vibrational_energies, const ArrayOfArrayOfSpeciesTag &abs_species, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: AtmRawRead.
void GriddedFieldPRegrid(GriddedField3 &gfraw_out, const Vector &p_grid, const GriddedField3 &gfraw_in_orig, const Index &interp_order, const Index &zeropadding, const Verbosity &verbosity)
WORKSPACE METHOD: GriddedFieldPRegrid.
void WindFieldsCalcExpand1D(Tensor3 &wind_u_field, Tensor3 &wind_v_field, Tensor3 &wind_w_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField3 &wind_u_field_raw, const GriddedField3 &wind_v_field_raw, const GriddedField3 &wind_w_field_raw, const Index &atmosphere_dim, const Index &interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: WindFieldsCalcExpand1D.
void atm_fields_compactFromMatrix(GriddedField4 &af, const Index &atmosphere_dim, const Matrix &im, const ArrayOfString &field_names, const Verbosity &)
WORKSPACE METHOD: atm_fields_compactFromMatrix.
void AtmosphereSet1D(Index &atmosphere_dim, Vector &lat_grid, Vector &lon_grid, const Verbosity &verbosity)
WORKSPACE METHOD: AtmosphereSet1D.
void AtmFieldsExpand1D(Tensor3 &t_field, Tensor3 &z_field, Tensor4 &vmr_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Index &atmosphere_dim, const Index &chk_vmr_nan, const Verbosity &)
WORKSPACE METHOD: AtmFieldsExpand1D.
void AtmFieldsExtract1D(Index &atmosphere_dim, Vector &lat_grid, Vector &lon_grid, Tensor3 &t_field, Tensor3 &z_field, Tensor4 &vmr_field, const Index &ilat, const Index &ilon, const Verbosity &verbosity)
WORKSPACE METHOD: AtmFieldsExtract1D.
void GriddedFieldPRegridHelper(Index &ing_min, Index &ing_max, ArrayOfLagrangeLogInterpolation &lag_p, Matrix &itw, GriddedField &gfraw_out, const GriddedField &gfraw_in, const Index p_grid_index, ConstVectorView p_grid, const Index &interp_order, const Index &zeropadding, const Verbosity &verbosity)
Calculate grid positions and interpolations weights for GriddedFieldPRegrid.
void AtmFieldsAndParticleBulkPropFieldFromCompact(Vector &p_grid, Vector &lat_grid, Vector &lon_grid, Tensor3 &t_field, Tensor3 &z_field, Tensor4 &vmr_field, Tensor4 &particle_bulkprop_field, ArrayOfString &particle_bulkprop_names, const ArrayOfArrayOfSpeciesTag &abs_species, const GriddedField4 &atm_fields_compact, const Index &atmosphere_dim, const String &delim, const Numeric &p_min, const Index &check_gridnames, const Verbosity &)
WORKSPACE METHOD: AtmFieldsAndParticleBulkPropFieldFromCompact.
void MagFieldsCalcIGRF(Tensor3 &mag_u_field, Tensor3 &mag_v_field, Tensor3 &mag_w_field, const Tensor3 &z_field, const Vector &lat_grid, const Vector &lon_grid, const Vector &refellipsoid, const Time &time, const Verbosity &)
WORKSPACE METHOD: MagFieldsCalcIGRF.
void MagFieldsCalc(Tensor3 &mag_u_field, Tensor3 &mag_v_field, Tensor3 &mag_w_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField3 &mag_u_field_raw, const GriddedField3 &mag_v_field_raw, const GriddedField3 &mag_w_field_raw, const Index &atmosphere_dim, const Index &interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: MagFieldsCalc.
void GriddedFieldZToPRegrid(GriddedField3 &gfraw_out, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const GriddedField3 &gfraw_in_orig, const Index &interp_order, const Index &zeropadding, const Verbosity &verbosity)
WORKSPACE METHOD: GriddedFieldZToPRegrid.
void GriddedFieldLatLonRegrid(GriddedField2 &gfraw_out, const Vector &lat_true, const Vector &lon_true, const GriddedField2 &gfraw_in_orig, const Index &interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: GriddedFieldLatLonRegrid.
void GriddedFieldLatLonExpand(GriddedField2 &gfraw_out, const GriddedField2 &gfraw_in_orig, const Verbosity &)
WORKSPACE METHOD: GriddedFieldLatLonExpand.
void AtmFieldPRegridHelper(Index &ing_min, Index &ing_max, ArrayOfLagrangeLogInterpolation &lag_p, Matrix &itw, ConstVectorView p_grid_out, ConstVectorView p_grid_in, const Index &interp_order, const Verbosity &verbosity)
Calculate grid positions and interpolations weights for AtmFieldPRegrid.
void p_gridDensify(Vector &p_grid, Index &atmfields_checked, Index &atmgeom_checked, Index &cloudbox_checked, const Vector &p_grid_old, const Index &nfill, const Verbosity &verbosity)
WORKSPACE METHOD: p_gridDensify.
void InterpAtmFieldToPosition(Numeric &outvalue, const Index &atmosphere_dim, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const Vector &rtp_pos, const Tensor3 &field, const Verbosity &verbosity)
WORKSPACE METHOD: InterpAtmFieldToPosition.
void batch_atm_fields_compactAddConstant(ArrayOfGriddedField4 &batch_atm_fields_compact, const String &name, const Numeric &value, const Index &prepend, const ArrayOfString &condensibles, const Verbosity &verbosity)
WORKSPACE METHOD: batch_atm_fields_compactAddConstant.
void FieldFromGriddedFieldCheckLatLonHelper(const Vector &lat_grid, const Vector &lon_grid, const Index ilat, const Index ilon, const GriddedField &gfield)
Check for correct grid dimensions.
void batch_atm_fields_compactAddSpecies(ArrayOfGriddedField4 &batch_atm_fields_compact, const String &name, const GriddedField3 &species, const Index &prepend, const Verbosity &verbosity)
WORKSPACE METHOD: batch_atm_fields_compactAddSpecies.
void AtmFieldsCalcExpand1D(Tensor3 &t_field, Tensor3 &z_field, Tensor4 &vmr_field, EnergyLevelMap &nlte_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField3 &t_field_raw, const GriddedField3 &z_field_raw, const ArrayOfGriddedField3 &vmr_field_raw, const ArrayOfGriddedField3 &nlte_field_raw, const ArrayOfQuantumIdentifier &nlte_ids, const Vector &nlte_energies, const Index &atmosphere_dim, const Index &interp_order, const Index &vmr_zeropadding, const Index &vmr_nonegative, const Index &nlte_when_negative, const Verbosity &verbosity)
WORKSPACE METHOD: AtmFieldsCalcExpand1D.
void atm_fields_compactAddSpecies(GriddedField4 &atm_fields_compact, const String &name, const GriddedField3 &species, const Index &prepend, const Verbosity &verbosity)
WORKSPACE METHOD: atm_fields_compactAddSpecies.
void MagFieldsCalcExpand1D(Tensor3 &mag_u_field, Tensor3 &mag_v_field, Tensor3 &mag_w_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField3 &mag_u_field_raw, const GriddedField3 &mag_v_field_raw, const GriddedField3 &mag_w_field_raw, const Index &atmosphere_dim, const Index &interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: MagFieldsCalcExpand1D.
constexpr Numeric GAS_CONSTANT
void AtmFieldsRefinePgrid(Vector &p_grid, Tensor3 &t_field, Tensor3 &z_field, Tensor4 &vmr_field, Index &atmfields_checked, Index &atmgeom_checked, Index &cloudbox_checked, const Vector &lat_grid, const Vector &lon_grid, const Index &atmosphere_dim, const Numeric &p_step, const Index &interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: AtmFieldsRefinePgrid.
void MagFieldsFromAltitudeRawCalc(Tensor3 &mag_u_field, Tensor3 &mag_v_field, Tensor3 &mag_w_field, const Vector &lat_grid, const Vector &lon_grid, const Tensor3 &z_field, const GriddedField3 &mag_u_field_raw, const GriddedField3 &mag_v_field_raw, const GriddedField3 &mag_w_field_raw, const Index &interp_order, const Numeric &extrapolation_factor, const Verbosity &verbosity)
WORKSPACE METHOD: MagFieldsFromAltitudeRawCalc.
void WindFieldsCalc(Tensor3 &wind_u_field, Tensor3 &wind_v_field, Tensor3 &wind_w_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField3 &wind_u_field_raw, const GriddedField3 &wind_v_field_raw, const GriddedField3 &wind_w_field_raw, const Index &atmosphere_dim, const Index &interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: WindFieldsCalc.
void AtmWithNLTERawRead(GriddedField3 &t_field_raw, GriddedField3 &z_field_raw, ArrayOfGriddedField3 &vmr_field_raw, ArrayOfGriddedField3 &nlte_field_raw, ArrayOfQuantumIdentifier &nlte_quantum_identifiers, Vector &nlte_vibrational_energies, const ArrayOfArrayOfSpeciesTag &abs_species, const String &basename, const Index &expect_vibrational_energies, const Verbosity &verbosity)
WORKSPACE METHOD: AtmWithNLTERawRead.
void GriddedFieldLatLonRegridHelper(ArrayOfLagrangeInterpolation &lag_lat, ArrayOfLagrangeCyclic0to360Interpolation &lag_lon, Tensor4 &itw, GriddedField &gfraw_out, const GriddedField &gfraw_in, const Index lat_grid_index, const Index lon_grid_index, ConstVectorView lat_true, ConstVectorView lon_true, const Index &interp_order, const Verbosity &verbosity)
Calculate grid positions and interpolations weights for GriddedFieldLatLonRegrid.
void atm_fields_compactAddConstant(GriddedField4 &af, const String &name, const Numeric &value, const Index &prepend, const ArrayOfString &condensibles, const Verbosity &verbosity)
WORKSPACE METHOD: atm_fields_compactAddConstant.
void atm_fields_compactCleanup(GriddedField4 &atm_fields_compact, const Numeric &threshold, const Verbosity &)
WORKSPACE METHOD: atm_fields_compactCleanup.
void atm_fields_compactExpand(GriddedField4 &af, Index &nf, const String &name, const Index &prepend, const Verbosity &)
atm_fields_compactExpand
void AtmosphereSet3D(Index &atmosphere_dim, Vector &lat_true, Vector &lon_true, const Verbosity &verbosity)
WORKSPACE METHOD: AtmosphereSet3D.
void FieldFromGriddedField(Matrix &field_out, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField2 &gfraw_in, const Verbosity &)
WORKSPACE METHOD: FieldFromGriddedField.
void z_surfaceConstantAltitude(Matrix &z_surface, const Vector &lat_grid, const Vector &lon_grid, const Numeric &altitude, const Verbosity &verbosity)
WORKSPACE METHOD: z_surfaceConstantAltitude.
constexpr Numeric EPSILON_LON_CYCLIC
Data value accuracy requirement for values at 0 and 360 deg if longitudes are cyclic.
void batch_atm_fields_compactCleanup(ArrayOfGriddedField4 &batch_atm_fields_compact, const Numeric &threshold, const Verbosity &verbosity)
WORKSPACE METHOD: batch_atm_fields_compactCleanup.
void AtmFieldsCalc(Tensor3 &t_field, Tensor3 &z_field, Tensor4 &vmr_field, EnergyLevelMap &nlte_field, const Vector &p_grid, const Vector &lat_grid, const Vector &lon_grid, const GriddedField3 &t_field_raw, const GriddedField3 &z_field_raw, const ArrayOfGriddedField3 &vmr_field_raw, const ArrayOfGriddedField3 &nlte_field_raw, const ArrayOfQuantumIdentifier &nlte_ids, const Vector &nlte_energies, const Index &atmosphere_dim, const Index &interp_order, const Index &vmr_zeropadding, const Index &vmr_nonegative, const Index &nlte_when_negative, const Verbosity &verbosity)
WORKSPACE METHOD: AtmFieldsCalc.
void batch_atm_fields_compactFromArrayOfMatrix(ArrayOfGriddedField4 &batch_atm_fields_compact, const Index &atmosphere_dim, const ArrayOfMatrix &am, const ArrayOfString &field_names, const Verbosity &verbosity)
WORKSPACE METHOD: batch_atm_fields_compactFromArrayOfMatrix.
void WindRawRead(GriddedField3 &wind_u_field_raw, GriddedField3 &wind_v_field_raw, GriddedField3 &wind_w_field_raw, const String &basename, const Verbosity &verbosity)
WORKSPACE METHOD: WindRawRead.
void AtmFieldPRegrid(Tensor3 &atmtensor_out, const Tensor3 &atmtensor_in_orig, const Vector &p_grid_new, const Vector &p_grid_old, const Index &interp_order, const Verbosity &verbosity)
WORKSPACE METHOD: AtmFieldPRegrid.
void z_surfaceFromFileAndGrid(Matrix &z_surface, const Vector &lat_grid, const Vector &lon_grid, const String &filename, const Index &interp_order, const Index &set_lowest_altitude_to_zero, const Verbosity &verbosity)
WORKSPACE METHOD: z_surfaceFromFileAndGrid.
void p_gridRefine(Vector &p_grid, Index &atmfields_checked, Index &atmgeom_checked, Index &cloudbox_checked, const Vector &p_grid_old, const Numeric &p_step10, const Verbosity &)
WORKSPACE METHOD: p_gridRefine.
Declarations having to do with the four output streams.
my_basic_string< char > String
The String type for ARTS.
constexpr Numeric ideal_gas_constant
Ideal gas constant [J/mol K].
constexpr Index GFIELD3_LON_GRID
Global constant, Index of the longitude grid in GriddedField3.
constexpr Index GFIELD4_LAT_GRID
Global constant, Index of the latitude grid in GriddedField4.
constexpr Index GFIELD4_P_GRID
Global constant, Index of the pressure grid in GriddedField4.
constexpr Index GFIELD4_FIELD_NAMES
Global constant, Index of the field names in GriddedField4.
constexpr Index GFIELD3_P_GRID
Global constant, Index of the pressure grid in GriddedField3.
constexpr Index GFIELD4_LON_GRID
Global constant, Index of the longitude grid in GriddedField4.
constexpr Index GFIELD3_LAT_GRID
Global constant, Index of the latitude grid in GriddedField3.
MagneticField compute(const Tensor3 &z_field, const Vector &lat_grid, const Vector &lon_grid, const Time &time, const Vector &ell)
Computes the magnetic field based on IGRF13 coefficients.
Array< QuantumIdentifier > ArrayOfQuantumIdentifier
Declaration of functions in rte.cc.
void p2gridpos(ArrayOfGridPos &gp, ConstVectorView old_pgrid, ConstVectorView new_pgrid, const Numeric &extpolfac)
Calculates grid positions for pressure values.
void interp_atmfield_by_gp(VectorView x, const Index &atmosphere_dim, ConstTensor3View x_field, const ArrayOfGridPos &gp_p, const ArrayOfGridPos &gp_lat, const ArrayOfGridPos &gp_lon)
Interpolates an atmospheric field given the grid positions.
void rte_pos2gridpos(GridPos &gp_p, GridPos &gp_lat, GridPos &gp_lon, const Index &atmosphere_dim, ConstVectorView p_grid, ConstVectorView lat_grid, ConstVectorView lon_grid, ConstTensor3View z_field, ConstVectorView rte_pos)
Converts a geographical position (rte_pos) to grid positions for p, lat and lon.
Header file for special_interp.cc.
ArrayOfQuantumIdentifier levels
void ThrowIfNotOK() const ARTS_NOEXCEPT
Structure to store a grid position.
Magnetic field for the east (u), north (v), and up (w) components as the ENU-coordinate system.
Class to handle time in ARTS.
This file contains basic functions to handle XML data files.
void xml_read_from_file(const String &filename, T &type, const Verbosity &verbosity)
Reads data from XML file.