1 // Copyright (C) 2013-2018 Alexander Barth <barth.alexander@gmail.com>
2 //
3 // This program is free software; you can redistribute it and/or modify it under
4 // the terms of the GNU General Public License as published by the Free Software
5 // Foundation; either version 2 of the License, or (at your option) any later
6 // version.
7 //
8 // This program is distributed in the hope that it will be useful, but WITHOUT
9 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
10 // FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
11 // details.
12 //
13 // You should have received a copy of the GNU General Public License along with
14 // this program; if not, see <http://www.gnu.org/licenses/>.
15
16
17 #include <octave/oct.h>
18 #include <octave/ov-cell.h>
19
20 #include <netcdf.h>
21
22 #include <string>
23 #include <map>
24 #include <iostream>
25 #include <algorithm>
26 #include <vector>
27 #include <inttypes.h>
28
29 #include "config.h"
30
31 // Octave constructor to hold an array of ints
32 #define NETCDF_INT_ARRAY int32NDArray
33 std::map<std::string, octave_value> netcdf_constants;
34
init()35 void init() {
36 #include "netcdf_constants.h"
37 }
38
check_err(int status)39 void check_err(int status)
40 {
41 if (status != NC_NOERR) error("%s",nc_strerror(status));
42 }
43
44 // convert name to upper-case and add "NC_" prefix if it is missing
normalize_ncname(std::string name)45 std::string normalize_ncname(std::string name) {
46 std::string prefix = "NC_";
47 std::string ncname = name;
48 // to upper case
49 std::transform(ncname.begin(), ncname.end(),ncname.begin(), ::toupper);
50
51 // add prefix if it is missing
52 if (ncname.substr(0, prefix.size()) != prefix) {
53 ncname = prefix + ncname;
54 }
55 return ncname;
56 }
57
netcdf_get_constant(octave_value ov)58 octave_value netcdf_get_constant(octave_value ov)
59 {
60 if (netcdf_constants.empty())
61 {
62 init();
63 }
64
65 if (ov.is_scalar_type())
66 {
67 return ov.scalar_value();
68 }
69
70 std::string name = ov.string_value();
71 name = normalize_ncname(name);
72 std::map<std::string, octave_value>::const_iterator cst = netcdf_constants.find(name);
73
74 if (cst != netcdf_constants.end ())
75 {
76 return cst->second;
77 }
78 else
79 {
80 error("unknown netcdf constant: %s",name.c_str());
81 return octave_value();
82 }
83 }
84
to_size_t(octave_value ov)85 size_t to_size_t(octave_value ov) {
86 size_t sz;
87 sz = static_cast<uint64_t>(ov.uint64_scalar_value());
88 return sz;
89 }
90
to_size_t_vector(octave_value ov,int len,size_t * vec)91 void to_size_t_vector(octave_value ov, int len, size_t *vec) {
92 uint64NDArray tmp = ov.uint64_array_value();
93
94 for (int i=0; i<len; i++)
95 {
96 vec[i] = static_cast<uint64_t>(tmp(len-i-1));
97 }
98
99 }
100
to_ptrdiff_t_vector(octave_value ov,int len,ptrdiff_t * vec)101 void to_ptrdiff_t_vector(octave_value ov, int len, ptrdiff_t *vec) {
102 int64NDArray tmp = ov.int64_array_value();
103
104 for (int i=0; i<len; i++)
105 {
106 vec[i] = static_cast<int64_t>(tmp(len-i-1));
107 }
108
109 }
110
start_count_stride(int ncid,int varid,octave_value_list args,int len,int ndims,size_t * start,size_t * count,ptrdiff_t * stride)111 void start_count_stride(int ncid, int varid, octave_value_list args,int len,
112 int ndims, size_t* start,size_t* count,ptrdiff_t* stride)
113 {
114 OCTAVE_LOCAL_BUFFER (int, dimids, ndims);
115 check_err(nc_inq_vardimid (ncid, varid, dimids));
116
117 // default values for start, count and stride
118 // i.e. all variable is loaded
119
120 for (int i=0; i<ndims; i++) {
121 check_err(nc_inq_dimlen(ncid,dimids[i],&(count[i])));
122 start[i] = 0;
123 //cout << "count def " << count[i] << " " << i << endl;
124 stride[i] = 1;
125 }
126
127 // start argument
128
129 if (len > 2)
130 {
131 if (args(2).numel() != ndims)
132 {
133 error("number of elements of argument %s should match the number "
134 "of dimension of the netCDF variable",
135 "start");
136 }
137
138 to_size_t_vector(args(2),ndims,start);
139
140 // if start is specified, the default for count is 1 (how odd!)
141 for (int i=0; i<ndims; i++)
142 {
143 count[i] = 1;
144 }
145 }
146
147 // count argument
148
149 if (len > 3)
150 {
151
152 if (args(3).numel() != ndims)
153 {
154 error("number of elements of argument %s should match the number "
155 "of dimension of the netCDF variable",
156 "count");
157 }
158
159 to_size_t_vector(args(3),ndims,count);
160 }
161
162 // stride argument
163
164 if (len > 4)
165 {
166 if (args(4).numel() != ndims)
167 {
168 error("number of elements of argument %s should match the number "
169 "of dimension of the netCDF variable",
170 "stride");
171 }
172
173 to_ptrdiff_t_vector(args(4),ndims,stride);
174 }
175 }
176
177
178 DEFUN_DLD(netcdf_getConstant, args,,
179 "-*- texinfo -*-\n\
180 @deftypefn {Loadable Function} {@var{value} =} netcdf_getConstant(@var{name}) \n\
181 Returns the value of a NetCDF constant called @var{name}.\n\
182 @seealso{netcdf_getConstantNames}\n\
183 @end deftypefn")
184 {
185 if (args.length() != 1) {
186 print_usage ();
187 return octave_value();
188 }
189
190 return netcdf_get_constant(args(0));
191 }
192
193
194 DEFUN_DLD(netcdf_getConstantNames, args,,
195 "-*- texinfo -*-\n\
196 @deftypefn {Loadable Function} {@var{value} =} netcdf_getConstantNames() \n\
197 Returns a list of all constant names.\n\
198 @end deftypefn\n\
199 @seealso{netcdf_getConstant}\n")
200 {
201
202 if (args.length() != 0)
203 {
204 print_usage ();
205 return octave_value();
206 }
207
208 if (netcdf_constants.empty())
209 {
210 init();
211 }
212
213 Cell c = Cell (dim_vector(1,netcdf_constants.size()));
214
215 int i = 0;
216 for (std::map<std::string, octave_value>::const_iterator p = netcdf_constants.begin ();
217 p != netcdf_constants.end (); p++) {
218 c(i++) = octave_value(p->first);
219 }
220
221 return octave_value(c);
222
223 }
224
225
226 DEFUN_DLD(netcdf_inqLibVers, args,,
227 "-*- texinfo -*-\n\
228 @deftypefn {Loadable Function} {@var{vers} =} netcdf_inqLibVers() \n\
229 Returns the version of the NetCDF library.\n\
230 @end deftypefn\n\
231 @seealso{netcdf_open}\n")
232 {
233 if (args.length() != 0)
234 {
235 print_usage ();
236 return octave_value ();
237 }
238
239 return octave_value(std::string(nc_inq_libvers()));
240 }
241
242 DEFUN_DLD(netcdf_setDefaultFormat, args,,
243 "-*- texinfo -*-\n\
244 @deftypefn {Loadable Function} {@var{old_format} =} netcdf_setDefaultFormat(@var{format}) \n\
245 Sets the default format of the NetCDF library and returns the previous default format (as a numeric value). @var{format} can be \n\
246 \"format_classic\", \"format_64bit\", \"format_netcdf4\" or \"format_netcdf4_classic\". \n\
247 @end deftypefn\n\
248 @seealso{netcdf_open}\n")
249 {
250 if (args.length() != 1)
251 {
252 print_usage ();
253 return octave_value ();
254 }
255
256 int format = netcdf_get_constant(args(0)).int_value();
257 int old_format;
258
259 check_err(nc_set_default_format(format, &old_format));
260
261 return octave_value(old_format);
262 }
263
264
265 // int nc_set_chunk_cache(size_t size, size_t nelems, float preemption);
266
267 DEFUN_DLD(netcdf_setChunkCache, args,,
268 "-*- texinfo -*-\n\
269 @deftypefn {Loadable Function} {} netcdf_setChunkCache(@var{size}, @var{nelems}, @var{preemption}) \n\
270 Sets the default chunk cache settins in the HDF5 library. The settings applies to all files which are subsequently opened or created.\n\
271 @end deftypefn\n\
272 @seealso{netcdf_getChunkCache}\n")
273 {
274 if (args.length() != 3)
275 {
276 print_usage ();
277 return octave_value();
278 }
279
280 size_t size = to_size_t(args(0));
281 size_t nelems = to_size_t(args(1));
282 float preemption = args(2).scalar_value();
283
284 check_err(nc_set_chunk_cache(size, nelems, preemption));
285
286 return octave_value();
287 }
288
289
290 // int nc_get_chunk_cache(size_t *sizep, size_t *nelemsp, float *preemptionp);
291
292 DEFUN_DLD(netcdf_getChunkCache, args,,
293 "-*- texinfo -*-\n\
294 @deftypefn {Loadable Function} {[@var{size}, @var{nelems}, @var{preemption}] =} netcdf_getChunkCache() \n\
295 Gets the default chunk cache settins in the HDF5 library. \n\
296 @end deftypefn\n\
297 @seealso{netcdf_setChunkCache}\n")
298 {
299 if (args.length() != 0)
300 {
301 print_usage ();
302 return octave_value ();
303 }
304
305 size_t size;
306 size_t nelems;
307 float preemption;
308
309 check_err(nc_get_chunk_cache(&size, &nelems, &preemption));
310 octave_value_list retval;
311 retval(0) = octave_value(size);
312 retval(1) = octave_value(nelems);
313 retval(2) = octave_value(preemption);
314
315 return retval;
316 }
317
318
319
320 DEFUN_DLD(netcdf_create, args,,
321 "-*- texinfo -*-\n\
322 @deftypefn {Loadable Function} {@var{ncid} =} netcdf_create(@var{filename},@var{mode}) \n\
323 Creates the file named @var{filename} in the mode @var{mode} which can have the \n\
324 following values: \n\
325 \"clobber\" (overwrite existing files), \n\
326 \"noclobber\" (prevent to overwrite existing files) \n\
327 \"64bit_offset\" (use the 64bit-offset format), \n\
328 \"netcdf4\" (use the NetCDF4, i.e. HDF5 format) or \n\
329 \"share\" (concurrent reading of the dataset). \n\
330 @var{mode} can also be the numeric value return by netcdf_getConstant. In the later-case it can be combined with a bitwise-or. \n\
331 @end deftypefn\n\
332 Example: \n\
333 @example \n\
334 mode = bitor(netcdf.getConstant(\"classic_model\"), ...\n\
335 netcdf.getConstant(\"netcdf4\")); \n\
336 ncid = netcdf.create(\"test.nc\",mode); \n\
337 @end example \n\
338 @seealso{netcdf_close}\n")
339 {
340
341 if (args.length() != 2)
342 {
343 print_usage ();
344 return octave_value ();
345 }
346
347 std::string filename = args(0).string_value();
348 int mode = netcdf_get_constant(args(1)).int_value();
349 int ncid;
350
351 check_err(nc_create(filename.c_str(), mode, &ncid));
352
353 return octave_value(ncid);
354 }
355
356 DEFUN_DLD(netcdf_open, args,,
357 "-*- texinfo -*-\n\
358 @deftypefn {Loadable Function} {@var{ncid} =} netcdf_open(@var{filename},@var{mode}) \n\
359 Opens the file named @var{filename} in the mode @var{mode}.\n\
360 @end deftypefn\n\
361 @seealso{netcdf_close}\n")
362 {
363
364 if (args.length() != 2)
365 {
366 print_usage ();
367 return octave_value();
368 }
369
370 std::string filename = args(0).string_value();
371 int mode = netcdf_get_constant(args(1)).int_value();
372 int ncid;
373
374 check_err(nc_open(filename.c_str(), mode, &ncid));
375
376 return octave_value(ncid);
377 }
378
379
380
381 DEFUN_DLD(netcdf_abort, args,,
382 "-*- texinfo -*-\n\
383 @deftypefn {Loadable Function} {} netcdf_abort(@var{ncid}) \n\
384 Aborts all changes since the last time the dataset entered in define mode.\n\
385 @end deftypefn\n\
386 @seealso{netcdf_reDef}\n")
387 {
388
389 if (args.length() != 1)
390 {
391 print_usage ();
392 return octave_value();
393 }
394
395 int ncid = args(0).scalar_value();
396
397 check_err(nc_abort(ncid));
398
399 return octave_value();
400 }
401
402
403 DEFUN_DLD(netcdf_sync, args,,
404 "-*- texinfo -*-\n\
405 @deftypefn {Loadable Function} {} netcdf_sync(@var{ncid}) \n\
406 Writes all changes to the disk and leaves the file open.\n\
407 @end deftypefn\n\
408 @seealso{netcdf_close}\n")
409 {
410
411 if (args.length() != 1)
412 {
413 print_usage ();
414 return octave_value();
415 }
416
417 int ncid = args(0).scalar_value();
418
419 check_err(nc_sync(ncid));
420
421 return octave_value();
422 }
423
424 DEFUN_DLD(netcdf_setFill, args,,
425 "-*- texinfo -*-\n\
426 @deftypefn {Loadable Function} {@var{old_mode} =} netcdf_setFill(@var{ncid},@var{fillmode}) \n\
427 Change the fill mode (@var{fillmode}) of the data set @var{ncid}. The previous value of the fill mode is returned. @var{fillmode} can be either \"fill\" or \"nofill\".\n\
428 @end deftypefn\n\
429 @seealso{netcdf_open}\n")
430 {
431
432 if (args.length() != 2)
433 {
434 print_usage ();
435 return octave_value();
436 }
437
438 int ncid = args(0).scalar_value();
439 int fillmode = netcdf_get_constant(args(1)).int_value();
440 int old_mode;
441
442 check_err (nc_set_fill (ncid, fillmode, &old_mode));
443
444 return octave_value(old_mode);
445 }
446
447
448 //int nc_inq (int ncid, int *ndimsp, int *nvarsp, int *ngattsp,
449 // int *unlimdimidp);
450 DEFUN_DLD(netcdf_inq, args,,
451 "-*- texinfo -*-\n\
452 @deftypefn {Loadable Function} {[@var{ndims},@var{nvars},@var{ngatts},@var{unlimdimid}] =} netcdf_inq(@var{ncid}) \n\
453 Return the number of dimension (@var{ndims}), the number of variables (@var{nvars}), the number of global attributes (@var{ngatts}) and the id of the unlimited dimension (@var{unlimdimid}). \n\
454 If no unlimited dimension is declared -1 is returned. For NetCDF4 files, one should use \n\
455 the function netcdf_inqUnlimDims as multiple unlimite dimension exists. \n\
456 @end deftypefn\n\
457 @seealso{netcdf_inqUnlimDims}\n")
458 {
459 if (args.length() != 1)
460 {
461 print_usage ();
462 return octave_value();
463 }
464
465 int ncid = args(0).scalar_value();
466 int ndims, nvars, ngatts, unlimdimid;
467 octave_value_list retval;
468
469 check_err(nc_inq(ncid,&ndims,&nvars,&ngatts,&unlimdimid));
470
471 retval(0) = octave_value(ndims);
472 retval(1) = octave_value(nvars);
473 retval(2) = octave_value(ngatts);
474 retval(3) = octave_value(unlimdimid);
475 return retval;
476 }
477
478 // int nc_inq_unlimdims(int ncid, int *nunlimdimsp, int *unlimdimidsp);
479 DEFUN_DLD(netcdf_inqUnlimDims, args,,
480 "-*- texinfo -*-\n\
481 @deftypefn {Loadable Function} {@var{unlimdimids} =} netcdf_inqUnlimDims(@var{ncid}) \n\
482 Return the id of all unlimited dimensions of the NetCDF file @var{ncid}.\n\
483 @end deftypefn\n\
484 @seealso{netcdf_inq}\n")
485 {
486 if (args.length() != 1)
487 {
488 print_usage ();
489 return octave_value();
490 }
491
492 int ncid = args(0).scalar_value();
493 int nunlimdims;
494
495 check_err(nc_inq_unlimdims(ncid, &nunlimdims, NULL));
496
497 OCTAVE_LOCAL_BUFFER(int,tmp,nunlimdims);
498 NETCDF_INT_ARRAY unlimdimids = NETCDF_INT_ARRAY(dim_vector(1,nunlimdims));
499 check_err(nc_inq_unlimdims(ncid, &nunlimdims, tmp));
500
501 for (int i=0; i < nunlimdims; i++) {
502 unlimdimids(i) = tmp[i];
503 }
504
505 return octave_value(unlimdimids);
506 }
507
508
509 // int nc_inq_format (int ncid, int *formatp);
510 DEFUN_DLD(netcdf_inqFormat, args,,
511 "-*- texinfo -*-\n\
512 @deftypefn {Loadable Function} {@var{format} =} netcdf_inqFormat(@var{ncid}) \n\
513 Return the NetCDF format of the dataset @var{ncid}.\n\
514 Format might be one of the following \n\
515 \"FORMAT_CLASSIC\", \"FORMAT_64BIT\", \"FORMAT_NETCDF4\" or \"FORMAT_NETCDF4_CLASSIC\" \n\
516 @end deftypefn\n\
517 @seealso{netcdf_inq}\n")
518 {
519
520 if (args.length() != 1)
521 {
522 print_usage ();
523 return octave_value();
524 }
525
526 int ncid = args(0).scalar_value();
527 int format;
528
529 check_err(nc_inq_format(ncid, &format));
530
531 if (format == NC_FORMAT_CLASSIC) {
532 return octave_value("FORMAT_CLASSIC");
533 }
534 if (format == NC_FORMAT_64BIT) {
535 return octave_value("FORMAT_64BIT");
536 }
537 if (format == NC_FORMAT_NETCDF4) {
538 return octave_value("FORMAT_NETCDF4");
539 }
540
541 return octave_value("FORMAT_NETCDF4_CLASSIC");
542 }
543
544 // int nc_def_dim (int ncid, const char *name, size_t len, int *dimidp);
545
546 DEFUN_DLD(netcdf_defDim, args,,
547 "-*- texinfo -*-\n\
548 @deftypefn {Loadable Function} {@var{dimid} =} netcdf_defDim(@var{ncid},@var{name},@var{len}) \n\
549 Define the dimension with the name @var{name} and the length @var{len} in the dataset @var{ncid}. The id of the dimension is returned.\n\
550 @end deftypefn\n\
551 @seealso{netcdf_defVar}\n")
552 {
553
554 if (args.length() != 3)
555 {
556 print_usage ();
557 return octave_value();
558 }
559
560 int ncid = args(0).scalar_value();
561 std::string name = args(1).string_value();
562 size_t len = to_size_t(args(2));
563 int dimid;
564
565 check_err(nc_def_dim (ncid, name.c_str(), len, &dimid));
566
567 return octave_value(dimid);
568 }
569
570
571 // int nc_rename_dim(int ncid, int dimid, const char* name);
572
573
574 DEFUN_DLD(netcdf_renameDim, args,,
575 "-*- texinfo -*-\n\
576 @deftypefn {Loadable Function} {} netcdf_renameDim(@var{ncid},@var{dimid},@var{name}) \n\
577 Renames the dimension with the id @var{dimid} in the data set @var{ncid}. @var{name} is the new name of the dimension.\n\
578 @end deftypefn\n\
579 @seealso{netcdf_defDim}\n")
580 {
581
582 if (args.length() != 3)
583 {
584 print_usage ();
585 return octave_value();
586 }
587
588 int ncid = args(0).scalar_value();
589 int dimid = args(1).scalar_value();
590 std::string name = args(2).string_value();
591
592 check_err(nc_rename_dim (ncid, dimid, name.c_str()));
593
594 return octave_value();
595 }
596
597 // int nc_def_var (int ncid, const char *name, nc_type xtype,
598 // int ndims, const int dimids[], int *varidp);
599
600 DEFUN_DLD(netcdf_defVar, args,,
601 "-*- texinfo -*-\n\
602 @deftypefn {Loadable Function} {@var{varid} = } netcdf_defVar(@var{ncid},@var{name},@var{xtype},@var{dimids}) \n\
603 Defines a variable with the name @var{name} in the dataset @var{ncid}. @var{xtype} can be \"byte\", \"ubyte\", \"short\", \"ushort\", \"int\", \"uint\", \"int64\", \"uint64\", \"float\", \"double\", \"char\" or the corresponding number as returned by netcdf_getConstant. The parameter @var{dimids} define the ids of the dimension. For scalar this parameter is the empty array ([]). The variable id is returned. \n\
604 @end deftypefn\n\
605 @seealso{netcdf_open,netcdf_defDim}\n")
606 {
607
608 if (args.length() != 4)
609 {
610 print_usage ();
611 return octave_value();
612 }
613
614 int ncid = args(0).scalar_value();
615 std::string name = args(1).string_value ();
616 int xtype = netcdf_get_constant(args(2)).int_value();;
617
618 Array<double> tmp;
619
620 if (!args(3).OV_ISEMPTY()) {
621 tmp = args(3).vector_value ();
622 }
623
624 OCTAVE_LOCAL_BUFFER (int, dimids, tmp.numel());
625
626 for (int i = 0; i < tmp.numel(); i++)
627 {
628 dimids[i] = tmp(tmp.numel()-i-1);
629 }
630
631 int varid;
632
633 check_err(nc_def_var (ncid, name.c_str(), xtype, tmp.numel(), dimids, &varid));
634
635 return octave_value(varid);
636 }
637
638
639 // int nc_rename_var(int ncid, int varid, const char* name);
640
641
642 DEFUN_DLD(netcdf_renameVar, args,,
643 "-*- texinfo -*-\n\
644 @deftypefn {Loadable Function} {} netcdf_renameVar(@var{ncid},@var{varid},@var{name}) \n\
645 Renames the variable with the id @var{varid} in the data set @var{ncid}. @var{name} is the new name of the variable.\n\
646 @end deftypefn\n\
647 @seealso{netcdf_defVar}\n")
648 {
649
650 if (args.length() != 3)
651 {
652 print_usage ();
653 return octave_value();
654 }
655
656 int ncid = args(0).scalar_value();
657 int varid = args(1).scalar_value();
658 std::string name = args(2).string_value();
659
660 check_err(nc_rename_var (ncid, varid, name.c_str()));
661
662 return octave_value();
663 }
664
665
666 // int nc_def_var_fill(int ncid, int varid, int no_fill, void *fill_value);
667 DEFUN_DLD(netcdf_defVarFill, args,,
668 "-*- texinfo -*-\n\
669 @deftypefn {Loadable Function} {} netcdf_defVarFill(@var{ncid},@var{varid},@var{no_fill},@var{fillvalue}) \n\
670 Define the fill-value settings of the NetCDF variable @var{varid}.\n\
671 If @var{no_fill} is false, then the values between no-contiguous writes are filled with the value @var{fill_value}. This is disabled by setting @var{no_fill} to true.\n\
672 @end deftypefn\n\
673 @seealso{netcdf_inqVarFill}\n")
674 {
675
676 if (args.length() != 4)
677 {
678 print_usage ();
679 return octave_value();
680 }
681
682 int ncid = args(0).scalar_value();
683 int varid = args(1).scalar_value();
684 int no_fill = args(2).scalar_value(); // boolean
685 octave_value fill_value = args(3);
686 nc_type xtype;
687
688 check_err(nc_inq_vartype (ncid, varid, &xtype));
689
690 switch (xtype)
691 {
692 #define OV_NETCDF_DEF_VAR_FILL(netcdf_type,c_type,method) \
693 case netcdf_type: \
694 { \
695 check_err(nc_def_var_fill(ncid, varid, no_fill, fill_value.method().fortran_vec())); \
696 break; \
697 }
698
699 OV_NETCDF_DEF_VAR_FILL(NC_BYTE, signed char, int8_array_value)
700 OV_NETCDF_DEF_VAR_FILL(NC_UBYTE, unsigned char, uint8_array_value)
701 OV_NETCDF_DEF_VAR_FILL(NC_SHORT, short, int16_array_value)
702 OV_NETCDF_DEF_VAR_FILL(NC_USHORT, unsigned short, uint16_array_value)
703 OV_NETCDF_DEF_VAR_FILL(NC_INT, int, int32_array_value)
704 OV_NETCDF_DEF_VAR_FILL(NC_UINT, unsigned int, uint32_array_value)
705 OV_NETCDF_DEF_VAR_FILL(NC_INT64, long long, int64_array_value)
706 OV_NETCDF_DEF_VAR_FILL(NC_UINT64, unsigned long long, uint64_array_value)
707
708 OV_NETCDF_DEF_VAR_FILL(NC_FLOAT, float, float_array_value)
709 OV_NETCDF_DEF_VAR_FILL(NC_DOUBLE,double,array_value)
710
711 OV_NETCDF_DEF_VAR_FILL(NC_CHAR, char, char_array_value)
712 }
713
714 return octave_value();
715 }
716
717
718
719 // int nc_def_var_fill(int ncid, int varid, int no_fill, void *fill_value);
720 DEFUN_DLD(netcdf_inqVarFill, args,,
721 "-*- texinfo -*-\n\
722 @deftypefn {Loadable Function} {[@var{no_fill},@var{fillvalue}] = } netcdf_inqVarFill(@var{ncid},@var{varid}) \n\
723 Determines the fill-value settings of the NetCDF variable @var{varid}.\n\
724 If @var{no_fill} is false, then the values between no-contiguous writes are filled with the value @var{fill_value}. This is disabled by setting @var{no_fill} to true.\n\
725 @end deftypefn\n\
726 @seealso{netcdf_defVarFill}\n")
727 {
728
729 if (args.length() != 2)
730 {
731 print_usage ();
732 return octave_value();
733 }
734
735 int ncid = args(0).scalar_value();
736 int varid = args(1).scalar_value();
737 int no_fill;
738 nc_type xtype;
739 octave_value_list retval;
740 octave_value data;
741
742 check_err(nc_inq_vartype (ncid, varid, &xtype));
743
744 switch (xtype)
745 {
746 #define OV_NETCDF_INQ_VAR_FILL(netcdf_type,c_type) \
747 case netcdf_type: \
748 { \
749 Array< c_type > fill_value = Array< c_type >(dim_vector(1,1)); \
750 check_err(nc_inq_var_fill(ncid, varid, &no_fill, \
751 fill_value.fortran_vec())); \
752 data = octave_value(fill_value); \
753 break; \
754 }
755
756 OV_NETCDF_INQ_VAR_FILL(NC_BYTE,octave_int8)
757 OV_NETCDF_INQ_VAR_FILL(NC_UBYTE,octave_uint8)
758 OV_NETCDF_INQ_VAR_FILL(NC_SHORT,octave_int16)
759 OV_NETCDF_INQ_VAR_FILL(NC_USHORT,octave_uint16)
760 OV_NETCDF_INQ_VAR_FILL(NC_INT,octave_int32)
761 OV_NETCDF_INQ_VAR_FILL(NC_UINT,octave_uint32)
762 OV_NETCDF_INQ_VAR_FILL(NC_INT64,octave_int64)
763 OV_NETCDF_INQ_VAR_FILL(NC_UINT64,octave_uint64)
764
765 OV_NETCDF_INQ_VAR_FILL(NC_FLOAT,float)
766 OV_NETCDF_INQ_VAR_FILL(NC_DOUBLE,double)
767
768 OV_NETCDF_INQ_VAR_FILL(NC_CHAR,char)
769 }
770
771 //cout << "xtype3 " << xtype << " " << NC_DOUBLE << std::endl;
772 retval(0) = octave_value(no_fill);
773 retval(1) = data;
774 return retval;
775 }
776
777
778
779
780 //nc_def_var_deflate(int ncid, int varid, int shuffle, int deflate,
781 // int deflate_level);
782 DEFUN_DLD(netcdf_defVarDeflate, args,,
783 "-*- texinfo -*-\n\
784 @deftypefn {Loadable Function} {} netcdf_defVarDeflate (@var{ncid},@var{varid},@var{shuffle},@var{deflate},@var{deflate_level}) \n\
785 Define the compression settings NetCDF variable @var{varid}.\n\
786 If @var{deflate} is true, then the variable is compressed. The compression level @var{deflate_level} is an integer between 0 (no compression) and 9 (maximum compression).\n\
787 @end deftypefn\n\
788 @seealso{netcdf_inqVarDeflate}\n")
789 {
790
791 if (args.length() != 5)
792 {
793 print_usage ();
794 return octave_value();
795 }
796
797 int ncid = args(0).scalar_value();
798 int varid = args(1).scalar_value();
799 int shuffle = args(2).scalar_value(); // boolean
800 int deflate = args(3).scalar_value(); // boolean
801 int deflate_level = args(4).scalar_value();
802
803 check_err(nc_def_var_deflate (ncid, varid, shuffle, deflate, deflate_level));
804 return octave_value();
805 }
806
807
808 //nc_inq_var_deflate(int ncid, int varid, int *shufflep,
809 // int *deflatep, int *deflate_levelp);
810 DEFUN_DLD(netcdf_inqVarDeflate, args,,
811 "-*- texinfo -*-\n\
812 @deftypefn {Loadable Function} {[@var{shuffle},@var{deflate},@var{deflate_level}] = } netcdf_inqVarDeflate (@var{ncid},@var{varid}) \n\
813 Determines the compression settings NetCDF variable @var{varid}.\n\
814 If @var{deflate} is true, then the variable is compressed. The compression level @var{deflate_level} is an integer between 0 (no compression) and 9 (maximum compression).\n\
815 @end deftypefn\n\
816 @seealso{netcdf_defVarDeflate}\n")
817 {
818
819 if (args.length() != 2)
820 {
821 print_usage ();
822 return octave_value();
823 }
824
825 int ncid = args(0).scalar_value();
826 int varid = args(1).scalar_value();
827 int shuffle, deflate, deflate_level;
828 octave_value_list retval;
829
830 int format;
831 check_err(nc_inq_format(ncid, &format));
832
833 // nc_inq_var_deflate returns garbage for classic or 64bit files
834 if (format == NC_FORMAT_CLASSIC || format == NC_FORMAT_64BIT) {
835 shuffle = 0;
836 deflate = 0;
837 deflate_level = 0;
838 }
839 else {
840 check_err(nc_inq_var_deflate(ncid, varid,
841 &shuffle,&deflate,&deflate_level));
842 }
843
844 retval(0) = octave_value(shuffle);
845 retval(1) = octave_value(deflate);
846 retval(2) = octave_value(deflate_level);
847
848 return retval;
849 }
850
851 //int nc_def_var_chunking(int ncid, int varid, int storage, size_t *chunksizesp);
852 //chunksizes can be ommited if storage is \"CONTIGUOUS\"
853 DEFUN_DLD(netcdf_defVarChunking, args,,
854 "-*- texinfo -*-\n\
855 @deftypefn {Loadable Function} {} netcdf_defVarChunking (@var{ncid},@var{varid},@var{storage},@var{chunkSizes}) \n\
856 Define the chunking settings of NetCDF variable @var{varid}.\n\
857 If @var{storage} is the string \"chunked\", the variable is stored by chunk of the size @var{chunkSizes}.\n\
858 If @var{storage} is the string \"contiguous\", the variable is stored in a contiguous way.\n\
859 @end deftypefn\n\
860 @seealso{netcdf_inqVarChunking}\n")
861 {
862
863 if (args.length() != 3 && args.length() != 4)
864 {
865 print_usage ();
866 return octave_value();
867 }
868
869 int ncid = args(0).scalar_value();
870 int varid = args(1).scalar_value();
871 std::string storagestr = args(2).string_value();
872 int storage;
873
874 std::transform(storagestr.begin(), storagestr.end(),storagestr.begin(), ::toupper);
875
876 if (storagestr == "CHUNKED") {
877 storage = NC_CHUNKED;
878 }
879 else if (storagestr == "CONTIGUOUS") {
880 storage = NC_CONTIGUOUS;
881 }
882 else {
883 error("unknown storage %s",storagestr.c_str());
884 return octave_value();
885 }
886
887 if (args.length() == 4) {
888 OCTAVE_LOCAL_BUFFER (size_t, chunksizes, args(3).numel());
889 to_size_t_vector(args(3), args(3).numel(),chunksizes);
890
891 check_err(nc_def_var_chunking(ncid, varid, storage, chunksizes));
892 }
893 else {
894 check_err(nc_def_var_chunking(ncid, varid, storage, NULL));
895 }
896
897 return octave_value();
898 }
899
900 //int nc_inq_var_chunking(int ncid, int varid, int *storagep, size_t *chunksizesp);
901 DEFUN_DLD(netcdf_inqVarChunking, args,,
902 "-*- texinfo -*-\n\
903 @deftypefn {Loadable Function} {[@var{storage},@var{chunkSizes}] = } netcdf_inqVarChunking (@var{ncid},@var{varid}) \n\
904 Determines the chunking settings of NetCDF variable @var{varid}.\n\
905 If @var{storage} is the string \"chunked\", the variable is stored by chunk of the size @var{chunkSizes}.\n\
906 If @var{storage} is the string \"contiguous\", the variable is stored in a contiguous way.\n\
907 @end deftypefn\n\
908 @seealso{netcdf_defVarChunking}\n")
909 {
910
911 if (args.length() != 2)
912 {
913 print_usage ();
914 return octave_value();
915 }
916
917 int ncid = args(0).scalar_value();
918 int varid = args(1).scalar_value();
919 int storage;
920 int ndims;
921 octave_value_list retval;
922
923 check_err(nc_inq_varndims (ncid, varid, &ndims));
924 OCTAVE_LOCAL_BUFFER (size_t, chunksizes, ndims);
925
926 check_err(nc_inq_var_chunking(ncid, varid, &storage, chunksizes));
927
928 if (storage == NC_CHUNKED) {
929 retval(0) = octave_value("chunked");
930 // should use uint32NDArray on 32-bit?
931 uint64NDArray chunkSizes = uint64NDArray(dim_vector(1,ndims));
932
933 for (int i = 0; i < ndims; i++)
934 {
935 chunkSizes(ndims-i-1) = chunksizes[i];
936 }
937 retval(1) = octave_value(chunkSizes);
938 }
939 else {
940 retval(0) = octave_value("contiguous");
941 retval(1) = octave_value(Array<double>());
942 }
943
944 return retval;
945 }
946
947 // nc_def_var_fletcher32(int ncid, int varid, int checksum);
948 DEFUN_DLD(netcdf_defVarFletcher32, args,,
949 "-*- texinfo -*-\n\
950 @deftypefn {Loadable Function} {} netcdf_defVarFletcher32(@var{ncid},@var{varid},@var{checksum}) \n\
951 Defines the checksum settings of the variable with the id @var{varid} in the data set @var{ncid}. If @var{checksum} is the string \"FLETCHER32\", then fletcher32 checksums will be turned on for this variable. If @var{checksum} is \"NOCHECKSUM\", then checksums will be disabled. \n\
952 @end deftypefn\n\
953 @seealso{netcdf_defVar,netcdf_inqVarFletcher32}\n")
954 {
955
956 if (args.length() != 3)
957 {
958 print_usage ();
959 return octave_value();
960 }
961
962 int ncid = args(0).scalar_value();
963 int varid = args(1).scalar_value();
964 int checksum = netcdf_get_constant(args(2)).int_value();
965
966 check_err(nc_def_var_fletcher32(ncid, varid, checksum));
967
968 return octave_value();
969 }
970
971
972
973 DEFUN_DLD(netcdf_inqVarFletcher32, args,,
974 "-*- texinfo -*-\n\
975 @deftypefn {Loadable Function} {@var{checksum} =} netcdf_inqVarFletcher32(@var{ncid},@var{varid}) \n\
976 Determines the checksum settings of the variable with the id @var{varid} in the data set @var{ncid}. If fletcher32 checksums is turned on for this variable, then @var{checksum} is the string \"FLETCHER32\". Otherwise it is the string \"NOCHECKSUM\". \n\
977 @end deftypefn\n\
978 @seealso{netcdf_defVar,netcdf_inqVarFletcher32}\n")
979 {
980
981 if (args.length() != 2)
982 {
983 print_usage ();
984 return octave_value();
985 }
986
987 int ncid = args(0).scalar_value();
988 int varid = args(1).scalar_value();
989 int checksum;
990
991 check_err(nc_inq_var_fletcher32(ncid, varid, &checksum));
992
993 if (checksum == NC_FLETCHER32)
994 {
995 return octave_value("FLETCHER32");
996 }
997 else
998 {
999 return octave_value("NOCHECKSUM");
1000 }
1001 }
1002
1003
1004
1005 DEFUN_DLD(netcdf_endDef, args,,
1006 "-*- texinfo -*-\n\
1007 @deftypefn {Loadable Function} {} netcdf_endDef (@var{ncid}) \n\
1008 Leaves define-mode of NetCDF file @var{ncid}.\n\
1009 @end deftypefn\n\
1010 @seealso{netcdf_reDef}\n")
1011 {
1012 if (args.length() != 1)
1013 {
1014 print_usage ();
1015 return octave_value();
1016 }
1017
1018 int ncid = args(0).scalar_value();
1019
1020 check_err(nc_enddef (ncid));
1021
1022 return octave_value();
1023 }
1024
1025 DEFUN_DLD(netcdf_reDef, args,,
1026 "-*- texinfo -*-\n\
1027 @deftypefn {Loadable Function} {} netcdf_reDef (@var{ncid}) \n\
1028 Enter define-mode of NetCDF file @var{ncid}.\n\
1029 @end deftypefn\n\
1030 @seealso{netcdf_endDef}\n")
1031 {
1032 if (args.length() != 1)
1033 {
1034 print_usage ();
1035 return octave_value();
1036 }
1037
1038 int ncid = args(0).scalar_value();
1039
1040 check_err(nc_redef (ncid));
1041
1042 return octave_value();
1043 }
1044
1045 // http://www.unidata.ucar.edu/software/netcdf/docs/netcdf-c/nc_005fput_005fvar_005f-type.html#nc_005fput_005fvar_005f-type
1046
1047 DEFUN_DLD(netcdf_putVar, args,,
1048 "-*- texinfo -*-\n\
1049 @deftypefn {Loadable Function} {} netcdf_putVar (@var{ncid},@var{varid},@var{data}) \n\
1050 @deftypefnx {Loadable Function} {} netcdf_putVar (@var{ncid},@var{varid},@var{start},@var{data}) \n\
1051 @deftypefnx {Loadable Function} {} netcdf_putVar (@var{ncid},@var{varid},@var{start},@var{count},@var{data}) \n\
1052 @deftypefnx {Loadable Function} {} netcdf_putVar (@var{ncid},@var{varid},@var{start},@var{count},@var{stride},@var{data}) \n\
1053 Put data in a NetCDF variable.\n\
1054 The data @var{data} is stored in the variable @var{varid} of the NetCDF file @var{ncid}. \n\
1055 @var{start} is the start index of each dimension (0-based and defaults to a vector of zeros), \n\
1056 @var{count} is the number of elements of to be written along each dimension (default all elements)\n\
1057 and @var{stride} is the sampling interval.\n\
1058 @end deftypefn\n\
1059 @seealso{netcdf_endDef}\n")
1060 {
1061 if (args.length() < 3 || args.length() > 6)
1062 {
1063 print_usage ();
1064 return octave_value();
1065 }
1066
1067 int ncid = args(0).scalar_value();
1068 int varid = args(1).scalar_value ();
1069 octave_value data = args(args.length()-1);
1070
1071 int ndims;
1072 check_err(nc_inq_varndims (ncid, varid, &ndims));
1073
1074 OCTAVE_LOCAL_BUFFER (size_t, start, ndims);
1075 OCTAVE_LOCAL_BUFFER (size_t, count, ndims);
1076 OCTAVE_LOCAL_BUFFER (ptrdiff_t, stride, ndims);
1077
1078 nc_type xtype;
1079
1080
1081 check_err(nc_inq_vartype (ncid, varid, &xtype));
1082 //int sliced_numel = tmp.numel();
1083
1084 start_count_stride(ncid, varid, args, args.length()-1, ndims, start, count, stride);
1085
1086 // check if count matched size(data)
1087
1088 switch (xtype)
1089 {
1090 #define OV_NETCDF_PUT_VAR(netcdf_type,c_type,method) \
1091 case netcdf_type: \
1092 { \
1093 check_err(nc_put_vars (ncid, varid, start, count, stride, (c_type*)data.method().fortran_vec())); \
1094 break; \
1095 }
1096
OV_NETCDF_PUT_VAR(NC_BYTE,signed char,int8_array_value)1097 OV_NETCDF_PUT_VAR(NC_BYTE, signed char, int8_array_value)
1098 OV_NETCDF_PUT_VAR(NC_UBYTE, unsigned char, uint8_array_value)
1099 OV_NETCDF_PUT_VAR(NC_SHORT, short, int16_array_value)
1100 OV_NETCDF_PUT_VAR(NC_USHORT, unsigned short, uint16_array_value)
1101 OV_NETCDF_PUT_VAR(NC_INT, int, int32_array_value)
1102 OV_NETCDF_PUT_VAR(NC_UINT, unsigned int, uint32_array_value)
1103 OV_NETCDF_PUT_VAR(NC_INT64, long long, int64_array_value)
1104 OV_NETCDF_PUT_VAR(NC_UINT64, unsigned long long, uint64_array_value)
1105
1106 OV_NETCDF_PUT_VAR(NC_FLOAT, float, float_array_value)
1107 OV_NETCDF_PUT_VAR(NC_DOUBLE,double,array_value)
1108
1109 OV_NETCDF_PUT_VAR(NC_CHAR, char, char_array_value)
1110 default:
1111 {
1112 error("unknown type %d" ,xtype);
1113 }
1114 }
1115 return octave_value();
1116 }
1117
1118
1119
1120 DEFUN_DLD(netcdf_getVar, args,,
1121 "-*- texinfo -*-\n\
1122 @deftypefn {Loadable Function} {@var{data} =} netcdf_getVar (@var{ncid},@var{varid}) \n\
1123 @deftypefnx {Loadable Function} {@var{data} =} netcdf_getVar (@var{ncid},@var{varid},@var{start}) \n\
1124 @deftypefnx {Loadable Function} {@var{data} =} netcdf_getVar (@var{ncid},@var{varid},@var{start},@var{count}) \n\
1125 @deftypefnx {Loadable Function} {@var{data} =} netcdf_getVar (@var{ncid},@var{varid},@var{start},@var{count},@var{stride}) \n\
1126 Get the data from a NetCDF variable.\n\
1127 The data @var{data} is loaded from the variable @var{varid} of the NetCDF file @var{ncid}. \n\
1128 @var{start} is the start index of each dimension (0-based and defaults to a vector of zeros), \n\
1129 @var{count} is the number of elements of to be written along each dimension (default all elements)\n\
1130 and @var{stride} is the sampling interval.\n\
1131 @end deftypefn\n\
1132 @seealso{netcdf_putVar}\n")
1133 {
1134 if (args.length() < 2 || args.length() > 5)
1135 {
1136 print_usage ();
1137 return octave_value();
1138 }
1139
1140 int ncid = args(0).scalar_value();
1141 int varid = args(1).scalar_value ();
1142 std::list<Range> ranges;
1143 int ndims;
1144 octave_value data;
1145 nc_type xtype;
1146
1147 check_err(nc_inq_vartype (ncid, varid, &xtype));
1148
1149 check_err(nc_inq_varndims (ncid, varid, &ndims));
1150
1151 //std::cout << "ndims " << ndims << std::endl;
1152
1153 OCTAVE_LOCAL_BUFFER (size_t, start, ndims);
1154 OCTAVE_LOCAL_BUFFER (size_t, count, ndims);
1155 OCTAVE_LOCAL_BUFFER (ptrdiff_t, stride, ndims);
1156
1157 int sz = 1;
1158
1159 dim_vector sliced_dim_vector;
1160
1161 if (ndims < 2)
1162 {
1163 sliced_dim_vector.resize(2);
1164 sliced_dim_vector(0) = 1;
1165 sliced_dim_vector(1) = 1;
1166 }
1167 else
1168 {
1169 sliced_dim_vector.resize(ndims);
1170 }
1171
1172 start_count_stride(ncid, varid, args, args.length(), ndims, start, count, stride);
1173 // std::cout << "count " << count[0] << std::endl;
1174 // std::cout << "start " << start[0] << std::endl;
1175 // std::cout << "stide " << stride[0] << std::endl;
1176
1177 // total size sz
1178 for (int i=0; i<ndims; i++) {
1179 sz = sz * count[i];
1180 sliced_dim_vector(i) = count[ndims-i-1];
1181 //sliced_dim_vector(i) = count[i];
1182 }
1183
1184 // std::cout << "sz " << sz << std::endl;
1185 // std::cout << "sliced_dim_vector " << sliced_dim_vector(0) << " x " << sliced_dim_vector(1) << std::endl;
1186
1187 // Array < float > arr = Array < float >(sliced_dim_vector);
1188 // float* time;
1189 // time = (float*)malloc(10 * sizeof(float));
1190 // check_err(nc_get_vars(ncid, varid, start, count, stride, time));
1191 // data = octave_value(arr);
1192 // return data;
1193
1194 switch (xtype)
1195 {
1196 #define OV_NETCDF_GET_VAR_CASE(netcdf_type,c_type) \
1197 case netcdf_type: \
1198 { \
1199 Array < c_type > arr = Array < c_type >(sliced_dim_vector); \
1200 /* necessary for netcdf 4.1.3 */ \
1201 if (sz > 0) { \
1202 check_err(nc_get_vars(ncid, varid, start, count, stride, arr.fortran_vec())); \
1203 } \
1204 data = octave_value(arr); \
1205 break; \
1206 }
1207
OV_NETCDF_GET_VAR_CASE(NC_BYTE,octave_int8)1208 OV_NETCDF_GET_VAR_CASE(NC_BYTE,octave_int8)
1209 OV_NETCDF_GET_VAR_CASE(NC_UBYTE,octave_uint8)
1210 OV_NETCDF_GET_VAR_CASE(NC_SHORT,octave_int16)
1211 OV_NETCDF_GET_VAR_CASE(NC_USHORT,octave_uint16)
1212 OV_NETCDF_GET_VAR_CASE(NC_INT,octave_int32)
1213 OV_NETCDF_GET_VAR_CASE(NC_UINT,octave_uint32)
1214 OV_NETCDF_GET_VAR_CASE(NC_INT64,octave_int64)
1215 OV_NETCDF_GET_VAR_CASE(NC_UINT64,octave_uint64)
1216
1217 OV_NETCDF_GET_VAR_CASE(NC_FLOAT,float)
1218 OV_NETCDF_GET_VAR_CASE(NC_DOUBLE,double)
1219
1220 OV_NETCDF_GET_VAR_CASE(NC_CHAR, char)
1221
1222 default:
1223 {
1224 error("unknown type %d" ,xtype);
1225 }
1226
1227 }
1228
1229 return data;
1230 }
1231
1232 DEFUN_DLD(netcdf_close, args,,
1233 "-*- texinfo -*-\n\
1234 @deftypefn {Loadable Function} {} netcdf_close(@var{ncid}) \n\
1235 Close the NetCDF file with the id @var{ncid}.\n\
1236 @end deftypefn\n\
1237 @seealso{netcdf_open}\n")
1238 {
1239
1240 if (args.length() != 1)
1241 {
1242 print_usage ();
1243 return octave_value();
1244 }
1245
1246 int ncid = args(0).scalar_value();
1247
1248 check_err(nc_close(ncid));
1249 return octave_value ();
1250 }
1251
1252
1253
1254 // int nc_inq_attname(int ncid, int varid, int attnum, char *name);
1255
1256 DEFUN_DLD(netcdf_inqAttName, args,,
1257 "-*- texinfo -*-\n\
1258 @deftypefn {Loadable Function} {@var{name} =} netcdf_inqAttName (@var{ncid},@var{varid},@var{attnum}) \n\
1259 Get the name of a NetCDF attribute.\n\
1260 This function returns the name of the attribute with the id @var{attnum} of the variable \n\
1261 @var{varid} in the NetCDF file @var{ncid}. For global attributes @var{varid} can be \n\
1262 netcdf_getConstant(\"global\").\n\
1263 @seealso{netcdf_inqAttName}\n\
1264 @end deftypefn")
1265 {
1266 if (args.length() != 3) {
1267 print_usage ();
1268 return octave_value();
1269 }
1270
1271 int ncid = args(0).scalar_value();
1272 int varid = args(1).scalar_value();
1273 int attnum = args(2).scalar_value();
1274 char name[NC_MAX_NAME+1];
1275
1276 check_err(nc_inq_attname(ncid, varid, attnum, name));
1277
1278 return octave_value(std::string(name));
1279 }
1280
1281
1282 DEFUN_DLD(netcdf_inqAttID, args,,
1283 "-*- texinfo -*-\n\
1284 @deftypefn {Loadable Function} {@var{attnum} =} netcdf_inqAttID(@var{ncid},@var{varid},@var{attname}) \n\
1285 Return the attribute id @var{attnum} of the attribute named @var{attname} of the variable @var{varid} in the dataset @var{ncid}. \n\
1286 For global attributes @var{varid} can be \n\
1287 netcdf_getConstant(\"global\").\n\
1288 @seealso{netcdf_inqAttName}\n\
1289 @end deftypefn")
1290 {
1291 if (args.length() != 3)
1292 {
1293 print_usage ();
1294 return octave_value ();
1295 }
1296 int ncid = args(0).scalar_value();
1297 int varid = args(1).scalar_value();
1298 std::string attname = args(2).string_value();
1299 int attnum;
1300
1301 check_err (nc_inq_attid (ncid, varid, attname.c_str(), &attnum));
1302
1303 return octave_value(attnum);
1304 }
1305
1306
1307 //int nc_inq_att (int ncid, int varid, const char *name,
1308 // nc_type *xtypep, size_t *lenp);
1309
1310 DEFUN_DLD(netcdf_inqAtt, args,,
1311 "-*- texinfo -*-\n\
1312 @deftypefn {Loadable Function} {[@var{xtype},@var{len}] = } netcdf_inqAtt(@var{ncid},@var{varid},@var{name}) \n\
1313 Get attribute type and length.\n\
1314 @seealso{netcdf_inqAttName}\n\
1315 @end deftypefn")
1316 {
1317 if (args.length() != 3)
1318 {
1319 print_usage ();
1320 return octave_value();
1321 }
1322
1323 int ncid = args(0).scalar_value();
1324 int varid = args(1).scalar_value();
1325 std::string name = args(2).string_value();
1326 int xtype;
1327 size_t len;
1328 octave_value_list retval;
1329
1330 check_err(nc_inq_att(ncid, varid, name.c_str(), &xtype, &len));
1331
1332 retval(0) = octave_value(xtype);
1333 retval(1) = octave_value(len);
1334 return retval;
1335 }
1336
1337
1338 DEFUN_DLD(netcdf_getAtt, args,,
1339 "-*- texinfo -*-\n\
1340 @deftypefn {Loadable Function} {@var{data} =} netcdf_getAtt (@var{ncid},@var{varid},@var{name}) \n\
1341 Get the value of a NetCDF attribute.\n\
1342 This function returns the value of the attribute called @var{name} of the variable \n\
1343 @var{varid} in the NetCDF file @var{ncid}. For global attributes @var{varid} can be \n\
1344 netcdf_getConstant(\"global\").\n\
1345 @seealso{netcdf_putAtt}\n\
1346 @end deftypefn")
1347 {
1348 if (args.length() != 3)
1349 {
1350 print_usage ();
1351 return octave_value();
1352 }
1353
1354 int ncid = args(0).scalar_value();
1355 int varid = args(1).scalar_value();
1356 std::string attname = args(2).string_value();
1357 nc_type xtype;
1358 size_t len;
1359 octave_value data;
1360
1361 check_err(nc_inq_att(ncid, varid, attname.c_str(), &xtype, &len));
1362
1363 #define OV_NETCDF_GET_ATT_CASE(netcdf_type,c_type) \
1364 if (xtype == netcdf_type) \
1365 { \
1366 Array< c_type > arr = Array< c_type >(dim_vector(1,len)); \
1367 check_err(nc_get_att(ncid, varid, attname.c_str(), arr.fortran_vec())); \
1368 data = octave_value(arr); \
1369 }
1370 OV_NETCDF_GET_ATT_CASE(NC_BYTE,octave_int8)
1371 OV_NETCDF_GET_ATT_CASE(NC_UBYTE,octave_uint8)
1372 OV_NETCDF_GET_ATT_CASE(NC_SHORT,octave_int16)
1373 OV_NETCDF_GET_ATT_CASE(NC_USHORT,octave_uint16)
1374 OV_NETCDF_GET_ATT_CASE(NC_INT,octave_int32)
1375 OV_NETCDF_GET_ATT_CASE(NC_UINT,octave_uint32)
1376 OV_NETCDF_GET_ATT_CASE(NC_INT64,octave_int64)
1377 OV_NETCDF_GET_ATT_CASE(NC_UINT64,octave_uint64)
1378
1379 OV_NETCDF_GET_ATT_CASE(NC_FLOAT,float)
1380 OV_NETCDF_GET_ATT_CASE(NC_DOUBLE,double)
1381
1382 OV_NETCDF_GET_ATT_CASE(NC_CHAR, char)
1383
1384
1385 return data;
1386 }
1387
1388
1389 DEFUN_DLD(netcdf_putAtt, args,,
1390 "-*- texinfo -*-\n\
1391 @deftypefn {Loadable Function} {} netcdf_putAtt (@var{ncid},@var{varid},@var{name},@var{data}) \n\
1392 Defines a NetCDF attribute.\n\
1393 This function defines the attribute called @var{name} of the variable \n\
1394 @var{varid} in the NetCDF file @var{ncid}. The value of the attribute will be @var{data}. \n\
1395 For global attributes @var{varid} can be \n\
1396 netcdf_getConstant(\"global\").\n\
1397 @seealso{netcdf_getAtt}\n\
1398 @end deftypefn")
1399 {
1400 if (args.length() != 4)
1401 {
1402 print_usage ();
1403 return octave_value();
1404 }
1405
1406 int ncid = args(0).scalar_value();
1407 int varid = args(1).scalar_value();
1408 std::string attname = args(2).string_value();
1409 octave_value data = args(3);
1410
1411 nc_type xtype;
1412
1413 // get matching netcdf type
1414
1415 if (data.is_string())
1416 xtype = NC_CHAR;
1417 else if (data.is_int8_type())
1418 xtype = NC_BYTE;
1419 else if (data.is_uint8_type())
1420 xtype = NC_UBYTE;
1421 else if (data.is_int16_type())
1422 xtype = NC_SHORT;
1423 else if (data.is_uint16_type())
1424 xtype = NC_USHORT;
1425 else if (data.is_int32_type())
1426 xtype = NC_INT;
1427 else if (data.is_uint32_type())
1428 xtype = NC_UINT;
1429 else if (data.is_int64_type())
1430 xtype = NC_INT64;
1431 else if (data.is_uint64_type())
1432 xtype = NC_UINT64;
1433 else if (data.is_single_type())
1434 xtype = NC_FLOAT;
1435 else
1436 xtype = NC_DOUBLE;
1437
1438 //cout << "xtype " << xtype << endl;
1439 size_t len = data.numel();
1440
1441 switch (xtype)
1442 {
1443 #define OV_NETCDF_PUT_ATT(netcdf_type,c_type,method) \
1444 case netcdf_type: \
1445 { \
1446 check_err(nc_put_att (ncid, varid, attname.c_str(), xtype, len, data.method().fortran_vec())); \
1447 break; \
1448 }
1449
1450 OV_NETCDF_PUT_ATT(NC_BYTE, signed char, int8_array_value)
1451 OV_NETCDF_PUT_ATT(NC_UBYTE, unsigned char, uint8_array_value)
1452 OV_NETCDF_PUT_ATT(NC_SHORT, short, int16_array_value)
1453 OV_NETCDF_PUT_ATT(NC_USHORT, unsigned short, uint16_array_value)
1454 OV_NETCDF_PUT_ATT(NC_INT, int, int32_array_value)
1455 OV_NETCDF_PUT_ATT(NC_UINT, unsigned int, uint32_array_value)
1456 OV_NETCDF_PUT_ATT(NC_INT64, long long, int64_array_value)
1457 OV_NETCDF_PUT_ATT(NC_UINT64, unsigned long long, uint64_array_value)
1458
1459 OV_NETCDF_PUT_ATT(NC_FLOAT, float, float_array_value)
1460 OV_NETCDF_PUT_ATT(NC_DOUBLE,double,array_value)
1461
1462 OV_NETCDF_PUT_ATT(NC_CHAR, char, char_array_value)
1463 }
1464
1465 /* check_err(nc_put_att (int ncid, int varid, const char *name, nc_type xtype,
1466 size_t len, const void *op));*/
1467
1468 return octave_value();
1469
1470 }
1471
1472
1473 DEFUN_DLD(netcdf_copyAtt, args,,
1474 "-*- texinfo -*-\n\
1475 @deftypefn {Loadable Function} {} netcdf_copyAtt (@var{ncid},@var{varid},@var{name},@var{ncid_out},@var{varid_out}) \n\
1476 Copies the attribute named @var{old_name} of the variable @var{varid} in the data set @var{ncid} \n\
1477 to the variable @var{varid_out} in the data set @var{ncid_out}. \n\
1478 To copy a global attribute use netcdf_getConstant(\"global\") for @var{varid} or @var{varid_out}.\n\
1479 @seealso{netcdf_getAtt,netcdf_getConstant}\n\
1480 @end deftypefn")
1481 {
1482
1483 if (args.length() != 5)
1484 {
1485 print_usage ();
1486 return octave_value ();
1487 }
1488
1489 int ncid = args(0).scalar_value();
1490 int varid = args(1).scalar_value();
1491 std::string name = args(2).string_value();
1492 int ncid_out = args(3).scalar_value();
1493 int varid_out = args(4).scalar_value();
1494
1495 check_err (nc_copy_att (ncid, varid, name.c_str(),
1496 ncid_out, varid_out));
1497
1498 return octave_value ();
1499 }
1500
1501
1502 DEFUN_DLD(netcdf_renameAtt, args,,
1503 "-*- texinfo -*-\n\
1504 @deftypefn {Loadable Function} {} netcdf_renameAtt(@var{ncid},@var{varid},@var{old_name},@var{new_name}) \n\
1505 Renames the attribute named @var{old_name} of the variable @var{varid} in the data set @var{ncid}. @var{new_name} is the new name of the attribute.\n\
1506 To rename a global attribute use netcdf_getConstant(\"global\") for @var{varid}.\n\
1507 @seealso{netcdf_copyAtt,netcdf_getConstant}\n\
1508 @end deftypefn")
1509 {
1510
1511 if (args.length() != 4)
1512 {
1513 print_usage ();
1514 return octave_value ();
1515 }
1516
1517 int ncid = args(0).scalar_value();
1518 int varid = args(1).scalar_value();
1519 std::string old_name = args(2).string_value();
1520 std::string new_name = args(3).string_value();
1521
1522 check_err(nc_rename_att (ncid, varid, old_name.c_str(), new_name.c_str()));
1523
1524 return octave_value ();
1525 }
1526
1527
1528 DEFUN_DLD(netcdf_delAtt, args,,
1529 "-*- texinfo -*-\n\
1530 @deftypefn {Loadable Function} {} netcdf_delAtt(@var{ncid},@var{varid},@var{name}) \n\
1531 Deletes the attribute named @var{name} of the variable @var{varid} in the data set @var{ncid}. \n\
1532 To delete a global attribute use netcdf_getConstant(\"global\") for @var{varid}.\n\
1533 @seealso{netcdf_defAtt,netcdf_getConstant}\n\
1534 @end deftypefn")
1535 {
1536
1537 if (args.length() != 3)
1538 {
1539 print_usage ();
1540 return octave_value ();
1541 }
1542
1543 int ncid = args(0).scalar_value();
1544 int varid = args(1).scalar_value();
1545 std::string name = args(2).string_value();
1546
1547 check_err(nc_del_att (ncid, varid, name.c_str()));
1548
1549 return octave_value ();
1550 }
1551
1552
1553 DEFUN_DLD(netcdf_inqVarID, args,,
1554 "-*- texinfo -*-\n\
1555 @deftypefn {Loadable Function} {@var{varid} = } netcdf_inqVarID (@var{ncid},@var{name}) \n\
1556 Return the id of a variable based on its name.\n\
1557 @seealso{netcdf_defVar,netcdf_inqVarIDs}\n\
1558 @end deftypefn")
1559 {
1560
1561 if (args.length() != 2)
1562 {
1563 print_usage ();
1564 return octave_value();
1565 }
1566
1567 int ncid = args(0).scalar_value();
1568 std::string varname = args(1).string_value();
1569 int varid;
1570
1571 check_err(nc_inq_varid(ncid,varname.c_str(), &varid));
1572
1573 return octave_value(varid);
1574 }
1575
1576 DEFUN_DLD(netcdf_inqVarIDs, args,,
1577 "-*- texinfo -*-\n\
1578 @deftypefn {Loadable Function} {@var{varids} = } netcdf_inqVarID (@var{ncid}) \n\
1579 Return all variable ids.\n\
1580 This functions returns all variable ids in a NetCDF file or NetCDF group.\n\
1581 @seealso{netcdf_inqVarID}\n\
1582 @end deftypefn")
1583 {
1584
1585 if (args.length() != 1)
1586 {
1587 print_usage ();
1588 return octave_value ();
1589 }
1590
1591 int ncid = args(0).scalar_value();
1592 int nvars;
1593
1594 check_err(nc_inq_varids(ncid, &nvars, NULL));
1595
1596 OCTAVE_LOCAL_BUFFER(int,tmp,nvars);
1597 NETCDF_INT_ARRAY varids = NETCDF_INT_ARRAY(dim_vector(1,nvars));
1598 check_err(nc_inq_varids(ncid, &nvars, tmp));
1599
1600 for (int i=0; i < nvars; i++) {
1601 varids(i) = tmp[i];
1602 }
1603
1604 return octave_value(varids);
1605 }
1606
1607 DEFUN_DLD(netcdf_inqVar, args,,
1608 "-*- texinfo -*-\n\
1609 @deftypefn {Loadable Function} {[@var{name},@var{nctype},@var{dimids},@var{nattr}] = } netcdf_inqVarID (@var{ncid},@var{varid}) \n\
1610 Inquires information about a NetCDF variable.\n\
1611 This functions returns the @var{name}, the NetCDF type @var{nctype}, an array of dimension ids \n\
1612 @var{dimids} and the number of attributes @var{nattr} of the NetCDF variable. @var{nctype} in an \n\
1613 integer corresponding NetCDF constants.\n\
1614 @seealso{netcdf_inqVarID,netcdf_getConstant}\n\
1615 @end deftypefn")
1616 {
1617
1618 if (args.length() != 2)
1619 {
1620 print_usage ();
1621 return octave_value();
1622 }
1623
1624 int ncid = args(0).scalar_value();
1625 int varid = args(1).scalar_value();
1626 char name[NC_MAX_NAME+1];
1627 int ndims, natts;
1628 nc_type xtype;
1629 octave_value_list retval;
1630
1631 check_err(nc_inq_varndims(ncid, varid, &ndims));
1632 OCTAVE_LOCAL_BUFFER (int, dimids, ndims);
1633
1634 check_err(nc_inq_var(ncid, varid, name, &xtype,
1635 &ndims, dimids, &natts));
1636
1637 retval(0) = octave_value(std::string(name));
1638 retval(1) = octave_value(xtype);
1639
1640 // copy output arguments
1641 Array<double> dimids_ = Array<double>(dim_vector(1,ndims));
1642 for (int i = 0; i < ndims; i++)
1643 {
1644 dimids_(i) = dimids[ndims-i-1];
1645 }
1646
1647 retval(2) = octave_value(dimids_);
1648 retval(3) = octave_value(natts);
1649
1650 return retval;
1651 }
1652
1653
1654
1655 DEFUN_DLD(netcdf_inqDim, args,,
1656 "-*- texinfo -*-\n\
1657 @deftypefn {Loadable Function} {[@var{name},@var{length}] =} netcdf_inqDim(@var{ncid},@var{dimid}) \n\
1658 Returns the name and length of a NetCDF dimension.\n\
1659 @seealso{netcdf_inqDimID}\n\
1660 @end deftypefn")
1661 {
1662
1663 if (args.length() != 2)
1664 {
1665 print_usage ();
1666 return octave_value();
1667 }
1668
1669 int ncid = args(0).scalar_value();
1670 int dimid = args(1).scalar_value();
1671 octave_value_list retval;
1672
1673 char name[NC_MAX_NAME+1];
1674 size_t length;
1675 check_err(nc_inq_dim(ncid, dimid, name, &length));
1676
1677 retval(0) = octave_value(std::string(name));
1678 retval(1) = octave_value(length);
1679
1680 return retval;
1681 }
1682
1683
1684 DEFUN_DLD(netcdf_inqDimID, args,,
1685 "-*- texinfo -*-\n\
1686 @deftypefn {Loadable Function} {@var{dimid} =} netcdf_inqDimID(@var{ncid},@var{dimname}) \n\
1687 Return the id of a NetCDF dimension.\n\
1688 @seealso{netcdf_inqDim}\n\
1689 @end deftypefn")
1690 {
1691
1692 if (args.length() != 2)
1693 {
1694 print_usage ();
1695 return octave_value();
1696 }
1697
1698 int ncid = args(0).scalar_value();
1699 std::string dimname = args(1).string_value();
1700 int id;
1701 octave_value_list retval;
1702
1703 check_err(nc_inq_dimid(ncid, dimname.c_str(), &id));
1704 retval(0) = octave_value(id);
1705
1706 return retval;
1707 }
1708
1709 // int nc_inq_dimids(int ncid, int *ndims, int *dimids, int include_parents);
1710 DEFUN_DLD(netcdf_inqDimIDs, args,,
1711 "-*- texinfo -*-\n\
1712 @deftypefn {Loadable Function} {@var{dimids} =} netcdf_inqDimID(@var{ncid}) \n\
1713 @deftypefnx {Loadable Function} {@var{dimids} =} netcdf_inqDimID(@var{ncid},@var{include_parents}) \n\
1714 Return the dimension ids defined in a NetCDF file.\n\
1715 If @var{include_parents} is 1, the dimension ids of the parent group are also returned.\n\
1716 Per default this is not the case (@var{include_parents} is 0).\n\
1717 @seealso{netcdf_inqDim}\n\
1718 @end deftypefn")
1719 {
1720 if (args.length() != 1 && args.length() != 2)
1721 {
1722 print_usage ();
1723 return octave_value();
1724 }
1725
1726 int ncid = args(0).scalar_value();
1727 int include_parents = 0;
1728 if (args.length() == 2)
1729 {
1730 include_parents = args(0).scalar_value();
1731 }
1732
1733 int ndims;
1734 check_err(nc_inq_ndims(ncid, &ndims));
1735
1736 OCTAVE_LOCAL_BUFFER(int,tmp,ndims);
1737 NETCDF_INT_ARRAY dimids = NETCDF_INT_ARRAY(dim_vector(1,ndims));
1738 check_err(nc_inq_dimids(ncid, &ndims, tmp, include_parents));
1739
1740 for (int i=0; i < ndims; i++) {
1741 dimids(i) = tmp[i];
1742 }
1743
1744 return octave_value(dimids);
1745 }
1746
1747
1748
1749 // groups
1750
1751 //int nc_def_grp(int parent_ncid, const char *name, int *new_ncid);
1752
1753 DEFUN_DLD(netcdf_defGrp, args,,
1754 "-*- texinfo -*-\n\
1755 @deftypefn {Loadable Function} {@var{new_ncid} =} netcdf_defGrp(@var{ncid},@var{name}) \n\
1756 Define a group in a NetCDF file.\n\
1757 @seealso{netcdf_inqGrps}\n\
1758 @end deftypefn")
1759 {
1760
1761 if (args.length() != 2)
1762 {
1763 print_usage ();
1764 return octave_value();
1765 }
1766
1767 int parent_ncid = args(0).scalar_value();
1768 std::string name = args(1).string_value();
1769 int new_ncid;
1770
1771 check_err(nc_def_grp(parent_ncid, name.c_str(), &new_ncid));
1772 return octave_value(new_ncid);
1773 }
1774
1775
1776 // int nc_inq_grps(int ncid, int *numgrps, int *ncids);
1777 DEFUN_DLD(netcdf_inqGrps, args,,
1778 "-*- texinfo -*-\n\
1779 @deftypefn {Loadable Function} {@var{ncids} =} netcdf_inqGrps(@var{ncid}) \n\
1780 Return all groups ids in a NetCDF file.\n\
1781 @seealso{netcdf_inqGrps}\n\
1782 @end deftypefn")
1783 {
1784 if (args.length() != 1)
1785 {
1786 print_usage ();
1787 return octave_value();
1788 }
1789
1790 int ncid = args(0).scalar_value();
1791 int numgrps;
1792
1793 check_err(nc_inq_grps(ncid, &numgrps, NULL));
1794
1795 OCTAVE_LOCAL_BUFFER(int,tmp,numgrps);
1796 NETCDF_INT_ARRAY ncids = NETCDF_INT_ARRAY(dim_vector(1,numgrps));
1797 check_err(nc_inq_grps(ncid, NULL, tmp));
1798
1799 for (int i=0; i < numgrps; i++) {
1800 ncids(i) = tmp[i];
1801 }
1802
1803 return octave_value(ncids);
1804 }
1805
1806 //int nc_inq_grpname(int ncid, char *name);
1807 DEFUN_DLD(netcdf_inqGrpName, args,,
1808 "-*- texinfo -*-\n\
1809 @deftypefn {Loadable Function} {@var{name} =} netcdf_inqGrpName(@var{ncid}) \n\
1810 Return group name in a NetCDF file.\n\
1811 @seealso{netcdf_inqGrps}\n\
1812 @end deftypefn")
1813 {
1814 if (args.length() != 1)
1815 {
1816 print_usage ();
1817 return octave_value();
1818 }
1819
1820 int ncid = args(0).scalar_value();
1821 char name[NC_MAX_NAME+1];
1822
1823 check_err(nc_inq_grpname(ncid, name));
1824
1825 return octave_value(std::string(name));
1826 }
1827
1828 //int nc_inq_grpname_full(int ncid, size_t *lenp, char *full_name);
1829 DEFUN_DLD(netcdf_inqGrpNameFull, args,,
1830 "-*- texinfo -*-\n\
1831 @deftypefn {Loadable Function} {@var{name} =} netcdf_inqGrpNameFull(@var{ncid}) \n\
1832 Return full name of group in NetCDF file.\n\
1833 @seealso{netcdf_inqGrpName}\n\
1834 @end deftypefn")
1835 {
1836 if (args.length() != 1)
1837 {
1838 print_usage ();
1839 return octave_value();
1840 }
1841
1842 int ncid = args(0).scalar_value();
1843 size_t len;
1844
1845 check_err(nc_inq_grpname_len(ncid,&len));
1846
1847 OCTAVE_LOCAL_BUFFER(char, name, len+1);
1848 if (name == 0)
1849 {
1850 error ("netcdf: error allocating buffer for name\n");
1851 return octave_value ();
1852 }
1853 octave_value retval;
1854 check_err(nc_inq_grpname_full(ncid, &len, name));
1855
1856 retval = octave_value(std::string(name));
1857 return retval;
1858 }
1859
1860 // int nc_inq_grp_parent(int ncid, int *parent_ncid);
1861 DEFUN_DLD(netcdf_inqGrpParent, args,,
1862 "-*- texinfo -*-\n\
1863 @deftypefn {Loadable Function} {@var{parent_ncid} =} netcdf_inqGrpParent(@var{ncid}) \n\
1864 Return id of the parent group\n\
1865 @seealso{netcdf_inqGrpName}\n\
1866 @end deftypefn")
1867 {
1868 if (args.length() != 1)
1869 {
1870 print_usage ();
1871 return octave_value();
1872 }
1873
1874 int ncid = args(0).scalar_value();
1875 int parent_ncid;
1876
1877 check_err(nc_inq_grp_parent(ncid, &parent_ncid));
1878 return octave_value(parent_ncid);
1879 }
1880
1881 // int nc_inq_grp_full_ncid(int ncid, char *full_name, int *grp_ncid);
1882 DEFUN_DLD(netcdf_inqGrpFullNcid, args,,
1883 "-*- texinfo -*-\n\
1884 @deftypefn {Loadable Function} {@var{grp_ncid} =} netcdf_inqGrpFullNcid(@var{ncid},@var{name}) \n\
1885 Return the group id based on the full group name.\n\
1886 @seealso{netcdf_inqGrpName}\n\
1887 @end deftypefn")
1888 {
1889 if (args.length() != 2)
1890 {
1891 print_usage ();
1892 return octave_value();
1893 }
1894
1895 int ncid = args(0).scalar_value();
1896 std::string name = args(1).string_value();
1897 int grp_ncid;
1898
1899 int format;
1900 check_err(nc_inq_format(ncid, &format));
1901
1902 if (format == NC_FORMAT_CLASSIC || format == NC_FORMAT_64BIT)
1903 {
1904 if (name == "/")
1905 {
1906 return octave_value(ncid);
1907 }
1908 else
1909 {
1910 error("groups are not supported in this format");
1911 return octave_value();
1912 }
1913 }
1914
1915 // nc_inq_grp_full_ncid makes a segmentation fault if
1916 // file is in non-HDF5 format
1917 check_err(nc_inq_grp_full_ncid(ncid, name.c_str(),&grp_ncid));
1918 return octave_value(grp_ncid);
1919 }
1920
1921
1922
1923 // int nc_inq_ncid(int ncid, const char *name, int *grp_ncid);
1924 DEFUN_DLD(netcdf_inqNcid, args,,
1925 "-*- texinfo -*-\n\
1926 @deftypefn {Loadable Function} {@var{grp_ncid} =} netcdf_inqNcid(@var{ncid},@var{name}) \n\
1927 Return group id based on its name\n\
1928 @seealso{netcdf_inqGrpFullNcid}\n\
1929 @end deftypefn")
1930 {
1931 if (args.length() != 2)
1932 {
1933 print_usage ();
1934 return octave_value();
1935 }
1936
1937 int ncid = args(0).scalar_value();
1938 std::string name = args(1).string_value();
1939 int grp_ncid;
1940
1941 check_err(nc_inq_ncid(ncid, name.c_str(), &grp_ncid));
1942 return octave_value(grp_ncid);
1943 }
1944